File: N:\mfix\model\des\desmpi_wrapper_mod.f
1
2
3
4
5
6
7 #include "version.inc"
8
9 module desmpi_wrapper
10
11 use parallel_mpi
12 use mpi_utility
13 use compar
14
15 interface des_mpi_irecv
16 module procedure des_mpi_irecv_db
17 end interface
18
19 interface des_mpi_isend
20 module procedure des_mpi_isend_db
21 end interface
22
23 interface des_mpi_scatterv
24 module procedure des_mpi_scatterv_i
25 module procedure des_mpi_scatterv_db
26 end interface
27
28 interface des_mpi_gatherv
29 module procedure des_mpi_gatherv_i
30 module procedure des_mpi_gatherv_db
31 end interface
32
33 contains
34
35
36
37
38
39
40 subroutine des_mpi_barrier
41 implicit none
42
43 character(len=80), parameter :: name = 'des_mpi_barrier'
44 #ifdef MPI
45 integer lerr
46 call mpi_barrier(mpi_comm_world, lerr)
47 call mpi_check( name //':mpi_barrier ', lerr )
48 #endif
49 return
50 end subroutine
51
52
53
54
55
56 subroutine des_mpi_irecv_db(precvbuf,precvcnt,ptoproc,ptag,preq,perr)
57 implicit none
58
59 double precision, dimension(:) :: precvbuf
60 integer :: precvcnt,ptoproc,ptag,preq,perr
61 #ifdef MPI
62 call mpi_irecv(precvbuf,precvcnt,mpi_double_precision,ptoproc,ptag,mpi_comm_world,preq,perr)
63 #endif
64 return
65 end subroutine
66
67
68
69
70
71 subroutine des_mpi_isend_db(psendbuf,psendcnt,ptoproc,ptag,preq,perr)
72 implicit none
73
74 double precision, dimension(:) :: psendbuf
75 integer :: psendcnt,ptoproc,ptag,preq,perr
76 #ifdef MPI
77 call mpi_isend(psendbuf,psendcnt,mpi_double_precision,ptoproc,ptag,mpi_comm_world,preq,perr)
78 #endif
79 return
80 end subroutine
81
82
83
84
85
86 subroutine des_mpi_wait(preq,perr)
87 implicit none
88
89 integer :: preq,perr
90 #ifdef MPI
91
92 integer :: lmpi_status(mpi_status_size)
93 call mpi_wait(preq,lmpi_status,perr)
94 #endif
95 return
96 end subroutine
97
98
99
100
101
102 subroutine des_mpi_scatterv_db(prootbuf,pscattercnts,pdispls, &
103 pprocbuf,precvcnt,proot,perr )
104 implicit none
105
106 double precision, dimension(:):: prootbuf,pprocbuf
107 integer, dimension (:) :: pdispls,pscattercnts
108 integer :: precvcnt,proot,perr
109
110 #ifdef MPI
111 call mpi_scatterv(prootbuf,pscattercnts,pdispls,mpi_double_precision, &
112 pprocbuf,precvcnt,mpi_double_precision,proot,mpi_comm_world,perr )
113 #else
114 pprocbuf = prootbuf
115 #endif
116 return
117 end subroutine
118
119
120
121
122 subroutine des_mpi_scatterv_i(prootbuf,pscattercnts,pdispls, &
123 pprocbuf, precvcnt,proot,perr )
124 implicit none
125
126 integer, dimension(:) :: prootbuf,pprocbuf
127 integer, dimension (:) :: pdispls,pscattercnts
128 integer :: precvcnt,proot,perr
129
130 #ifdef MPI
131 call mpi_scatterv(prootbuf,pscattercnts,pdispls,mpi_integer, &
132 pprocbuf,precvcnt,mpi_integer,proot,mpi_comm_world,perr )
133 #else
134 pprocbuf = prootbuf
135 #endif
136 return
137 end subroutine
138
139
140
141
142 subroutine des_mpi_gatherv_db(psendbuf,psendcnt,precvbuf, &
143 precvcnts, pdispls,proot,perr )
144 implicit none
145
146 double precision, dimension(:) :: psendbuf,precvbuf
147 integer, dimension (:) :: pdispls,precvcnts
148 integer :: psendcnt,proot,perr
149
150 #ifdef MPI
151 call mpi_gatherv(psendbuf,psendcnt,mpi_double_precision,precvbuf,precvcnts, &
152 pdispls,mpi_double_precision,proot,mpi_comm_world,perr )
153 #else
154 precvbuf = psendbuf
155 #endif
156 return
157 end subroutine
158
159
160
161
162 subroutine des_mpi_gatherv_i(psendbuf,psendcnt,precvbuf, &
163 precvcnts, pdispls,proot,perr )
164 implicit none
165
166 integer, dimension(:) :: psendbuf,precvbuf
167 integer, dimension (:) :: pdispls,precvcnts
168 integer :: psendcnt,proot,perr
169
170 #ifdef MPI
171 call mpi_gatherv(psendbuf,psendcnt,mpi_integer,precvbuf,precvcnts, &
172 pdispls,mpi_integer,proot,mpi_comm_world,perr )
173 #else
174 precvbuf = psendbuf
175 #endif
176 return
177 end subroutine
178
179
180
181
182
183 subroutine des_mpi_stop(myid)
184
185 USE funits
186 implicit none
187
188 integer, optional, intent(in) :: myid
189
190 #ifdef MPI
191
192 INTEGER :: mylid, ERRORCODE
193
194 if (.not. present(myid)) then
195 mylid = myPE
196 else
197 mylid = myid
198 endif
199
200 write(*,100) mylid
201 write(UNIT_LOG,100) mylid
202
203 100 format(/,'*****************',&
204 '********************************************',/, &
205 '(PE ',I2,') : A fatal error occurred in des routines',/,9X, &
206 '*.LOG file may contain other error messages ',/,'*****************', &
207 '********************************************',/)
208
209
210 call MPI_ABORT(MPI_COMM_WORLD, ERRORCODE, mpierr)
211 write(*,"('(PE ',I2,') : MPI_ABORT return = ',I2)") &
212 mylid,mpierr
213
214 call MPI_Finalize(mpierr)
215
216 ERROR_STOP 'MPI terminated from des_mpi_stop'
217 #else
218 ERROR_STOP 'terminated from des_mpi_stop'
219 #endif
220 end subroutine
221
222 end module
223