File: /nfs/home/0/users/jenkins/mfix.git/model/des/mpi_node_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     ! Contains following subroutines:
16     !    des_par_exchange
17     !    desmpi_sendrecv_wait, desmpi_gatherv, desmpi_scatterv
18     !    desmpi_check_sendrecvbuf,
19     !    desmpi_pack_ghostpar, desmpi_unpack_ghostpar, desmpi_cleanup,
20     !    desmpi_pack_parcross, desmpi_unpack_parcross,
21     !    des_addnodevalues, des_addnodevalues2,
22     !    des_gather_d,l,i
23     
24     !    des_restart_neigh, redim_par, des_dbgmpi
25     ! Contains following functions:
26     !    locate_par, exten_locate_par
27     !------------------------------------------------------------------------
28           module mpi_node_des
29     
30     !-----------------------------------------------
31     ! Modules
32     !-----------------------------------------------
33           use parallel_mpi
34           use mpi_utility
35           use discretelement
36           use desgrid
37           use compar
38           use physprop
39           use sendrecv
40           use des_bc
41           use desmpi_wrapper
42           use sendrecvnode
43           use mfix_pic
44           use des_thermo
45           use run, only: ENERGY_EQ,ANY_SPECIES_EQ
46           use param, only: DIMENSION_N_s
47           use des_rxns
48           use desmpi
49     
50           use mpi_comm_des, only: desmpi_sendrecv_init
51           use mpi_comm_des, only: desmpi_sendrecv_wait
52     
53           contains
54     
55     
56     !------------------------------------------------------------------------
57     ! Subroutine       : des_addnodevalues_mean_fields
58     ! Purpose          : This routine is specially used for computing mean
59     !                    fields by backward interpolation.
60     !
61     ! Parameters       : None
62     !------------------------------------------------------------------------
63           subroutine des_addnodevalues_mean_fields()
64     
65           use functions
66     !-----------------------------------------------
67           implicit none
68     !-----------------------------------------------
69     ! local variables
70     !-----------------------------------------------
71           integer :: li, lj, lk
72           integer :: lm,ijk,lface,lijkmin,lijkmax
73           integer :: linode,ljnode,lknode,lijknode
74     !-----------------------------------------------
75     
76     ! fill the temporary buffer
77           DO LM = 1,DES_MMAX
78              CALL DES_EXCHANGENODE(DES_ROPS_NODE(:,LM),PADD=.TRUE.)
79              DO LI =1,DIMN
80                 CALL DES_EXCHANGENODE(DES_VEL_NODE(:,LI,LM),PADD=.TRUE.)
81              END DO
82           END DO
83     
84     ! adjust for periodic boundaries with no domain decomposition
85           if (des_periodic_walls_x .and. nodesi.eq.1) then
86              do lk = kstart2,kend2
87              do lj = jstart2,jend2
88                 lijkmin = funijk(1,lj,lk)
89                 lijkmax = funijk(imax1,lj,lk)
90                 des_rops_node(lijkmin,:)  = des_rops_node(lijkmin,:)+des_rops_node(lijkmax,:)
91                 des_vel_node(lijkmin,:,:) = des_vel_node(lijkmin,:,:)+des_vel_node(lijkmax,:,:)
92                 des_rops_node(lijkmax,:)  = des_rops_node(lijkmin,:)
93                 des_vel_node(lijkmax,:,:) = des_vel_node(lijkmin,:,:)
94              end do
95              end do
96           end if
97           if (des_periodic_walls_y .and. nodesj.eq.1) then
98              do lk = kstart2,kend2
99              do li = istart2,iend2
100                 lijkmin = funijk(li,1,lk)
101                 lijkmax = funijk(li,jmax1,lk)
102                 des_rops_node(lijkmin,:)  = des_rops_node(lijkmin,:)+des_rops_node(lijkmax,:)
103                 des_vel_node(lijkmin,:,:) = des_vel_node(lijkmin,:,:)+des_vel_node(lijkmax,:,:)
104                 des_rops_node(lijkmax,:)  = des_rops_node(lijkmin,:)
105                 des_vel_node(lijkmax,:,:) = des_vel_node(lijkmin,:,:)
106              end do
107              end do
108           end if
109           if (des_periodic_walls_z .and. nodesk.eq.1 .and. do_K) then
110              do li = istart2,iend2
111              do lj = jstart2,jend2
112                 lijkmin = funijk(li,lj,1)
113                 lijkmax = funijk(li,lj,kmax1)
114                 des_rops_node(lijkmin,:)  = des_rops_node(lijkmin,:)+des_rops_node(lijkmax,:)
115                 des_vel_node(lijkmin,:,:) = des_vel_node(lijkmin,:,:)+des_vel_node(lijkmax,:,:)
116                 des_rops_node(lijkmax,:)  = des_rops_node(lijkmin,:)
117                 des_vel_node(lijkmax,:,:) = des_vel_node(lijkmin,:,:)
118              end do
119              end do
120           end if
121     
122           return
123     
124           end subroutine des_addnodevalues_mean_fields
125     
126     
127     
128     !------------------------------------------------------------------------
129     ! Subroutine       : des_addnodevalues
130     ! Purpose          : This routine is specially used for des_drag_gs
131     !                    The backward interpolation in des_drag_gs computes
132     !                    the grid node values of drag_am and drag_bm
133     !                    node values are from istart2 to iend1;
134     !                    hence a separate module is created to exchange
135     !                    node values
136     !
137     ! Parameters       : None
138     !------------------------------------------------------------------------
139           subroutine des_addnodevalues()
140     
141           use functions
142     !-----------------------------------------------
143           implicit none
144     !-----------------------------------------------
145     ! local variables
146     !-----------------------------------------------
147           integer :: li, lj, lk
148           integer :: lm,ijk,lface,lijkmin,lijkmax
149           integer :: linode,ljnode,lknode,lijknode
150     !-----------------------------------------------
151     
152     ! fill the temporary buffer
153           call des_exchangenode(drag_am, padd=.true.)
154           do li =1,dimn
155              call des_exchangenode(drag_bm(:,li), padd=.true.)
156           end do
157     
158     ! adjust for periodic boundaries with no domain decomposition
159           if (des_periodic_walls_x .and. nodesi.eq.1) then
160              do lk = kstart2,kend2
161              do lj = jstart2,jend2
162                 lijkmin = funijk(1,lj,lk)
163                 lijkmax = funijk(imax1,lj,lk)
164                 drag_am(lijkmin) = drag_am(lijkmin)+drag_am(lijkmax)
165                 drag_bm(lijkmin,:) = drag_bm(lijkmin,:)+drag_bm(lijkmax,:)
166                 drag_am(lijkmax) = drag_am(lijkmin)
167                 drag_bm(lijkmax,:) = drag_bm(lijkmin,:)
168              end do
169              end do
170           end if
171           if (des_periodic_walls_y .and. nodesj.eq.1) then
172              do lk = kstart2,kend2
173              do li = istart2,iend2
174                 lijkmin = funijk(li,1,lk)
175                 lijkmax = funijk(li,jmax1,lk)
176                 drag_am(lijkmin) = drag_am(lijkmin)+drag_am(lijkmax)
177                 drag_bm(lijkmin,:) = drag_bm(lijkmin,:)+drag_bm(lijkmax,:)
178                 drag_am(lijkmax) = drag_am(lijkmin)
179                 drag_bm(lijkmax,:) = drag_bm(lijkmin,:)
180              end do
181              end do
182           end if
183           if (des_periodic_walls_z .and. nodesk.eq.1 .and. do_K) then
184              do li = istart2,iend2
185              do lj = jstart2,jend2
186                 lijkmin = funijk(li,lj,1)
187                 lijkmax = funijk(li,lj,kmax1)
188                 drag_am(lijkmin) = drag_am(lijkmin)+drag_am(lijkmax)
189                 drag_bm(lijkmin,:) = drag_bm(lijkmin,:)+drag_bm(lijkmax,:)
190                 drag_am(lijkmax) = drag_am(lijkmin)
191                 drag_bm(lijkmax,:) = drag_bm(lijkmin,:)
192              end do
193              end do
194           end if
195     
196           return
197     
198           end subroutine des_addnodevalues
199     
200     
201     !------------------------------------------------------------------------
202     ! Subroutine       : des_addnodevalues2
203     ! Purpose          : This routine is specially used for calc_des_rop_s
204     !                    The backward interpolation in calc_des_rop_s computes
205     !                    the grid node values of des_rops_node
206     !                    node values are from istart2 to iend1;
207     !                    hence a separate module is created to exchange
208     !                    node values
209     !
210     ! Parameters       : None
211     !------------------------------------------------------------------------
212           subroutine des_addnodevalues2()
213     
214           use functions
215     !-----------------------------------------------
216           implicit none
217     !-----------------------------------------------
218     ! local variables
219     !-----------------------------------------------
220           integer :: li, lj, lk
221           integer :: lm,ijk,lface,lijkmin,lijkmax
222           integer :: linode,ljnode,lknode,lijknode
223     !-----------------------------------------------
224     
225     ! fill the temporary buffer
226           do lm = 1,DES_MMAX
227              call des_exchangenode(des_rops_node(:,lm),padd=.true.)
228           end do
229     
230     ! adjust for periodic boundaries with no domain decomposition
231           if (des_periodic_walls_x .and. nodesi.eq.1) then
232              do lk = kstart2,kend2
233              do lj = jstart2,jend2
234                 lijkmin = funijk(1,lj,lk)
235                 lijkmax = funijk(imax1,lj,lk)
236                 des_rops_node(lijkmin,:) = des_rops_node(lijkmin,:)+des_rops_node(lijkmax,:)
237                 des_rops_node(lijkmax,:) = des_rops_node(lijkmin,:)
238              end do
239              end do
240           end if
241           if (des_periodic_walls_y .and. nodesj.eq.1) then
242              do lk = kstart2,kend2
243              do li = istart2,iend2
244                 lijkmin = funijk(li,1,lk)
245                 lijkmax = funijk(li,jmax1,lk)
246                 des_rops_node(lijkmin,:) = des_rops_node(lijkmin,:)+des_rops_node(lijkmax,:)
247                 des_rops_node(lijkmax,:) = des_rops_node(lijkmin,:)
248              end do
249              end do
250           end if
251           if (des_periodic_walls_z .and. nodesk.eq.1 .and. do_K) then
252              do li = istart2,iend2
253              do lj = jstart2,jend2
254                 lijkmin = funijk(li,lj,1)
255                 lijkmax = funijk(li,lj,kmax1)
256                 des_rops_node(lijkmin,:) = des_rops_node(lijkmin,:)+des_rops_node(lijkmax,:)
257                 des_rops_node(lijkmax,:) = des_rops_node(lijkmin,:)
258              end do
259              end do
260           end if
261     
262           return
263     
264           end subroutine des_addnodevalues2
265     
266     
267     
268           end module mpi_node_des
269