File: /nfs/home/0/users/jenkins/mfix.git/model/des/mpi_comm_des_mod.f

1     !------------------------------------------------------------------------
2     ! Module           : desmpi
3     ! Purpose          : Contains wrapper class for mpi communications- send,recv
4     !
5     ! Author           : Pradeep.G
6     !
7     ! Purpose          : Module contains subroutines and variables related to
8     !                    des mpi communication.
9     !
10     ! Comments         : do_nsearch flag should be set to true before calling
11     !                    des_par_exchange; when do_nsearch is true ghost particles of the
12     !                    system will be updated, which will be later used to generate
13     !                    neighbour list.
14     
15     !------------------------------------------------------------------------
16           module mpi_comm_des
17     
18     !-----------------------------------------------
19     ! Modules
20     !-----------------------------------------------
21           use parallel_mpi
22           use mpi_utility
23           use discretelement
24           use desgrid
25           use compar
26           use physprop
27           use sendrecv
28           use des_bc
29           use desmpi_wrapper
30           use sendrecvnode
31           use mfix_pic
32           use des_thermo
33           use run, only: ENERGY_EQ,ANY_SPECIES_EQ
34           use param, only: DIMENSION_N_s
35           use des_rxns
36     
37           use desmpi
38     
39     !-----------------------------------------------
40     
41     ! generic interface definition
42           interface des_gather
43              module procedure des_gather_l,des_gather_i,des_gather_d
44           end interface
45     
46           contains
47     
48     !------------------------------------------------------------------------
49     ! Subroutine       : desmpi_sendrecv_init
50     ! Purpose          : posts asynchronous send and recv and updates the request id
51     !
52     ! Parameters       : pface - face number (1to6)
53     !                    debug - for printing debug statments
54     !------------------------------------------------------------------------
55           subroutine desmpi_sendrecv_init(pface,pdebug)
56     !-----------------------------------------------
57           implicit none
58     !-----------------------------------------------
59     ! dummy variables
60     !-----------------------------------------------
61           integer,intent(in) :: pface
62           integer,intent(in),optional :: pdebug
63     !-----------------------------------------------
64     ! local variables
65     !-----------------------------------------------
66           character(len=80), parameter :: name = 'desmpi_sendrecv_init'
67           integer :: ldebug,ltag,lerr,lrecvface
68     !-----------------------------------------------
69     
70     ! set the debug flag
71           ldebug = 0
72           if (present(pdebug)) then
73             ldebug = pdebug
74           endif
75     
76     !direct copy in case of single processor
77           lrecvface = pface+mod(pface,2)-mod(pface+1,2)
78           if (ineighproc(pface).eq.mype) then
79              drecvbuf(1:isendcnt(pface),lrecvface)= dsendbuf(1:isendcnt(pface),pface)
80           else
81              ltag = message_tag(ineighproc(pface),mype,pface)
82              call des_mpi_irecv(drecvbuf(:,pface),imaxbuf, &
83                                 ineighproc(pface),ltag,irecvreq(pface),lerr)
84              call mpi_check( name //':mpi_irecv ', lerr )
85     
86              ltag = message_tag(mype,ineighproc(pface),lrecvface)
87              call des_mpi_isend(dsendbuf(:,pface),isendcnt(pface), &
88                             ineighproc(pface),ltag,isendreq(pface),lerr)
89              call mpi_check( name //':mpi_isend ', lerr )
90           end if
91           return
92     
93         contains
94     
95           integer function message_tag(lsource,ldest,lrecvface)
96             implicit none
97             integer, intent(in) :: lsource,ldest,lrecvface
98             message_tag = lsource+numpes*ldest+numpes*numpes*lrecvface+100
99           end function message_tag
100     
101         end subroutine desmpi_sendrecv_init
102     
103     !------------------------------------------------------------------------
104     ! Subroutine       : desmpi_sendrecv_wait
105     ! Purpose          : waits for the communication for the specified interface
106     !
107     ! Parameters       : pface - face number (1to6)
108     !                    debug - for printing debug statments
109     !------------------------------------------------------------------------
110           subroutine desmpi_sendrecv_wait(pface,pdebug)
111     !-----------------------------------------------
112           implicit none
113     !-----------------------------------------------
114     ! dummy variables
115     !-----------------------------------------------
116           integer,intent(in) :: pface
117           integer,intent(in),optional :: pdebug
118     !-----------------------------------------------
119     ! local variables
120     !-----------------------------------------------
121           character(len=80), parameter :: name = 'desmpi_sendrecv_wait'
122           integer :: ldebug,lerr
123     !-----------------------------------------------
124     
125     ! set the debug flag
126           ldebug = 0
127           if (present(pdebug)) then
128             ldebug = pdebug
129           endif
130     
131     ! wait for both send and recv request completes
132           if (ineighproc(pface).ne.mype) then
133              call des_mpi_wait(isendreq(pface),lerr)
134              call mpi_check( name //':mpi_wait-send', lerr )
135              call des_mpi_wait(irecvreq(pface),lerr)
136              call mpi_check( name //':mpi_wait-recv', lerr )
137           end if
138           return
139           end subroutine desmpi_sendrecv_wait
140     
141     
142     !------------------------------------------------------------------------
143     ! Subroutine       : desmpi_scatterv
144     ! Purpose          : scatters the particle from PE_IO
145     ! Parameters       : ptype - flag for datatype integer (1) or double precision (2)
146     !                    pdebug - optional flag for debugging
147     !------------------------------------------------------------------------
148           subroutine desmpi_scatterv(ptype,pdebug)
149     !-----------------------------------------------
150           implicit none
151     !-----------------------------------------------
152     ! dummy variables
153     !-----------------------------------------------
154           integer, intent(in) :: ptype
155           integer, intent(in),optional :: pdebug
156     !-----------------------------------------------
157     ! local variables
158     !-----------------------------------------------
159           integer lroot,lidebug,lerr
160           character(len=80), parameter :: name = 'desmpi_scatterv'
161     !-----------------------------------------------
162     
163           lroot = pe_io
164           if (.not. present(pdebug)) then
165              lidebug = 0
166           else
167              lidebug = pdebug
168           endif
169     
170           if (ptype .eq. 1) then
171              call des_MPI_Scatterv(irootbuf,iscattercnts,idispls, &
172                                    iprocbuf,iscr_recvcnt,lroot,lerr)
173           else
174              call des_MPI_Scatterv(drootbuf,iscattercnts,idispls, &
175                                    dprocbuf,iscr_recvcnt,lroot,lerr)
176           end if
177           call MPI_Check( name //':MPI_Scatterv', lerr )
178     
179           return
180           end subroutine desmpi_scatterv
181     
182     
183     !------------------------------------------------------------------------
184     ! Subroutine       : desmpi_gatherv
185     ! Purpose          : gathers the particle from local proc to root proc
186     ! Parameters       : ptype - flag for datatype integer (1) or double precision (2)
187     !                    pdebug - optional flag for debugging
188     !------------------------------------------------------------------------
189           subroutine desmpi_gatherv(ptype,pdebug)
190     !-----------------------------------------------
191           implicit none
192     !-----------------------------------------------
193     ! dummy variables
194     !-----------------------------------------------
195           integer, intent(in) :: ptype
196           integer, intent(in),optional :: pdebug
197     !-----------------------------------------------
198     ! local variables
199     !-----------------------------------------------
200           integer lroot,lidebug,lerr
201           character(len=80), parameter :: name = 'des_gather'
202     !-----------------------------------------------
203     
204           lroot = pe_io
205           if (.not. present(pdebug)) then
206              lidebug = 0
207           else
208              lidebug = pdebug
209           endif
210           if(ptype.eq.1) then
211              call des_MPI_Gatherv(iprocbuf,igath_sendcnt,irootbuf, &
212                                   igathercnts,idispls,lroot,lerr)
213           else
214              call des_MPI_Gatherv(dprocbuf,igath_sendcnt,drootbuf, &
215                                   igathercnts,idispls,lroot,lerr)
216           end if
217           call MPI_Check( name //':MPI_Gatherv', lerr )
218     
219           return
220           end subroutine desmpi_gatherv
221     
222     
223     !------------------------------------------------------------------------
224     ! Subroutine       : des_gather_d
225     ! Purpose          : gathers double precision array from local to root
226     ! Parameters       :
227     !                    parray - array to be writen
228     !------------------------------------------------------------------------
229           subroutine des_gather_d(parray)
230     !-----------------------------------------------
231           implicit none
232     !-----------------------------------------------
233     ! dummy variables
234     !-----------------------------------------------
235           double precision, dimension(:) :: parray
236     !-----------------------------------------------
237     ! local variables
238     !-----------------------------------------------
239           integer :: lcurpar,lparcount,lcount
240     !-----------------------------------------------
241     
242     ! pack the variables in case of
243           lparcount = 1
244           lcount = 0
245           do lcurpar = 1, max_pip
246              if (lparcount.gt.pip) exit
247              if (.not. pea(lcurpar,1)) cycle
248              lparcount = lparcount +1
249              if (pea(lcurpar,4)) cycle
250              lcount = lcount + 1
251              dprocbuf(lcount) = parray(lcurpar)
252           end do
253           call desmpi_gatherv(ptype=2)
254           end subroutine des_gather_d
255     
256     
257     !------------------------------------------------------------------------
258     ! Subroutine       : des_gather_l
259     ! Purpose          : gathers logical array from local to root
260     ! Parameters       :
261     !                    parray - array to be writen
262     !------------------------------------------------------------------------
263           subroutine des_gather_l(parray)
264     !-----------------------------------------------
265           implicit none
266     !-----------------------------------------------
267     ! dummy variables
268     !-----------------------------------------------
269           logical, dimension(:) :: parray
270     !-----------------------------------------------
271     ! local variables
272     !-----------------------------------------------
273           integer :: lcurpar,lparcount,lcount
274     !-----------------------------------------------
275     
276     ! pack the variables in proc buffer
277           lparcount = 1
278           lcount = 0
279           do lcurpar = 1, max_pip
280              if (lparcount.gt.pip) exit
281              if (.not. pea(lcurpar,1)) cycle
282              lparcount = lparcount +1
283              if (pea(lcurpar,4)) cycle
284              lcount = lcount + 1
285              if(parray(lcurpar)) then
286                 iprocbuf(lcount) = 1
287              else
288                 iprocbuf(lcount) = 0
289              end if
290           end do
291           call desmpi_gatherv(ptype=1)
292     
293           end subroutine des_gather_l
294     
295     
296     !------------------------------------------------------------------------
297     ! Subroutine       : des_gather_i
298     ! Purpose          : gathers integer array from local to root
299     ! Parameters       :
300     !                    parray - array to be writen
301     !                    ploc2glb - this flag is used to conver local particle
302     !                    number into global particle number (used for history
303     !                    and neighbour terms)
304     !------------------------------------------------------------------------
305           subroutine des_gather_i(parray,ploc2glb)
306     !-----------------------------------------------
307           implicit none
308     !-----------------------------------------------
309     ! dummy variables
310     !-----------------------------------------------
311           integer, dimension(:) :: parray
312           logical,optional :: ploc2glb
313     !-----------------------------------------------
314     ! local variables
315     !-----------------------------------------------
316           integer :: lcurpar,lparcount,lcount
317           logical :: lloc2glb
318     !-----------------------------------------------
319     
320           if (present(ploc2glb)) then
321              lloc2glb = ploc2glb
322           else
323              lloc2glb = .false.
324           end if
325     ! pack the variables in proc buffer
326           lparcount = 1
327           lcount = 0
328           if (lloc2glb) then
329              do lcurpar = 1, max_pip
330                 if (lparcount.gt.pip) exit
331                 if (.not. pea(lcurpar,1)) cycle
332                 lparcount = lparcount +1
333                 if (pea(lcurpar,4)) cycle
334                 lcount = lcount + 1
335                 if(parray(lcurpar).gt.0) then
336                    iprocbuf(lcount) = iglobal_id(parray(lcurpar))
337                 else
338                    iprocbuf(lcount) = 0
339                 end if
340              end do
341           else
342              do lcurpar = 1, max_pip
343                 if (lparcount.gt.pip) exit
344                 if (.not. pea(lcurpar,1)) cycle
345                 lparcount = lparcount +1
346                 if (pea(lcurpar,4)) cycle
347                 lcount = lcount + 1
348                 iprocbuf(lcount) = parray(lcurpar)
349              end do
350           end if
351           call desmpi_gatherv(ptype=1)
352     
353           end subroutine des_gather_i
354     
355           end module mpi_comm_des
356