Observatoire de Paris - Code Fortran  version1.0
mod_all_object.f90
Aller à la documentation de ce fichier.
2 
3  implicit none
4  private
5 
9  public :: blocklist, size
10  public :: linecomm
11 
12  interface scalafx_cpl2g
13  module procedure cpl2g_real, cpl2g_dreal, cpl2g_complex, cpl2g_dcomplex
14  module procedure cpl2g_int
15  end interface scalafx_cpl2g
16 
17  interface scalafx_cpg2l
18  module procedure cpg2l_real, cpg2l_dreal, cpg2l_complex, cpg2l_dcomplex
19  module procedure cpg2l_int
20  end interface scalafx_cpg2l
21 
22  interface scalafx_addl2g
23  module procedure addl2g_real, addl2g_dreal, addl2g_complex, addl2g_dcomplex
24  module procedure addl2g_int
25  end interface scalafx_addl2g
26 
27  interface scalafx_addg2l
28  module procedure addg2l_real, addg2l_dreal, addg2l_complex, addg2l_dcomplex
29  module procedure addg2l_int
30  end interface scalafx_addg2l
31 
33  module procedure writearray_master_int
34  module procedure writearray_master_real, writearray_master_dreal
35  module procedure writearray_master_complex, writearray_master_dcomplex
36  end interface writearray_master
37 
38  interface writearray_slave
39  module procedure writearray_slave_int
40  module procedure writearray_slave_real, writearray_slave_dreal
41  module procedure writearray_slave_complex, writearray_slave_dcomplex
42  end interface writearray_slave
43 
44  interface readarray_master
45  module procedure readarray_master_int
46  module procedure readarray_master_real, readarray_master_dreal
47  module procedure readarray_master_complex, readarray_master_dcomplex
48  end interface readarray_master
49 
50  interface readarray_slave
51  module procedure readarray_slave_int
52  module procedure readarray_slave_real, readarray_slave_dreal
53  module procedure readarray_slave_complex, readarray_slave_dcomplex
54  end interface readarray_slave
55 
56 
77 type :: blocklist
78  private
79  integer :: nn
80  integer :: nb
81  integer :: nproc
82  integer :: myproc
83  integer :: srcproc
84  integer :: nblock
85  public
86  integer :: init_val
87 contains
89  procedure :: init => blocklist_init
90 
92  procedure :: getsize => blocklist_getsize
93 
95  procedure :: getblock => blocklist_getblock
96 end type blocklist
97 
98 interface size
99  module procedure blocklist_getsize
100 end interface size
101 
102 contains
103 
104 
105 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
106 !!! Blocklist
107 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
108 
115  subroutine blocklist_init(self, mygrid, desc, rowcol)
116  class(blocklist), intent(inout) :: self
117  type(blacsgrid), intent(in) :: mygrid
118  integer, intent(in) :: desc(DLEN_)
119  character, intent(in) :: rowcol
120 
121  integer :: nblockall, nextrablock, mydist,mysize
122 
123  if (rowcol == "c" .or. rowcol == "C") then
124  self%nn = desc(n_)
125  self%nb = desc(nb_)
126  self%nproc = mygrid%ncol
127  self%myproc = mygrid%mycol
128  self%srcproc = desc(csrc_)
129  else
130  self%nn = desc(m_)
131  self%nb = desc(mb_)
132  self%nproc = mygrid%nrow
133  self%myproc = mygrid%myrow
134  self%srcproc = desc(rsrc_)
135  end if
136  nblockall = self%nn / self%nb
137  self%nblock = nblockall / self%nproc
138  nextrablock = mod(nblockall, self%nproc)
139  mydist = mod(self%nproc + self%myproc - self%srcproc, self%nproc)
140  if (mydist < nextrablock) then
141  self%nblock = self%nblock + 1
142  elseif (mydist == nextrablock .and. mod(self%nn, self%nb) > 0) then
143  self%nblock = self%nblock +1
144  end if
145  mysize= blocklist_getsize(self)
146  end subroutine blocklist_init
147 
151  function blocklist_getsize(self) result(res)
152  class(blocklist), intent(in) :: self
153  integer :: res
154 
155  res = self%nblock
156 
157  end function blocklist_getsize
158 
159 
166  subroutine blocklist_getblock(self, iblock, iglob, iloc, bsize)
167  class(blocklist), intent(in) :: self
168  integer, intent(in) :: iblock
169  integer, intent(out) :: iglob, iloc, bsize
170 
171  integer :: mydist
172 
173  if (iblock >= 1 .and. iblock <= self%nblock) then
174  mydist = mod(self%nproc + self%myproc - self%srcproc, self%nproc)
175  iglob = ((iblock - 1) * self%nproc + mydist) * self%nb + 1
176  iloc = (iblock - 1) * self%nb + mod(iglob - 1, self%nb) + 1
177  bsize = min(self%nb, self%nn - iglob + 1)
178  else
179  iglob = 0
180  iloc = 0
181  bsize = 0
182  end if
183 
184  end subroutine blocklist_getblock
185 
186 
187 
188 
189 
190 
191 end module mod_scalapackfx_tools
subroutine blocklist_getblock(self, iblock, iglob, iloc, bsize)
Returns the indices (local and global) of a local block.
subroutine blocklist_init(self, mygrid, desc, rowcol)
integer function blocklist_getsize(self)