FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
fstr_contact_comm.f90
Go to the documentation of this file.
1!-------------------------------------------------------------------------------
2! Copyright (c) 2019 FrontISTR Commons
3! This software is released under the MIT License, see LICENSE.txt
4!-------------------------------------------------------------------------------
5
7
8 use hecmw
9
10 implicit none
11
12 private
13 public :: fstrst_contact_comm
22
24 private
25 integer(kind=kint) :: n_neighbor_pe
26 integer(kind=kint), pointer :: neighbor_pe(:)
27 integer(kind=kint) :: mpi_comm
28 integer(kind=kint), pointer :: ext_index(:)
29 integer(kind=kint), pointer :: ext_item(:)
30 integer(kind=kint), pointer :: int_index(:)
31 integer(kind=kint), pointer :: int_item(:)
32 end type fstrst_contact_comm
33
34 integer(kind=kint), parameter :: op_overwrite = 46810
35
36 integer(kind=kint), parameter :: debug = 0
37
38contains
39
40 subroutine fstr_contact_comm_init(conComm, hecMESH, ndof, n_contact_dof, contact_dofs)
41 implicit none
42 type (fstrst_contact_comm), intent(inout) :: concomm
43 type (hecmwst_local_mesh), intent(in) :: hecmesh
44 integer(kind=kint), intent(in) :: ndof, n_contact_dof
45 integer(kind=kint), intent(in) :: contact_dofs(:)
46 integer(kind=kint), pointer :: ext_index(:), ext_item(:), int_index(:), int_item(:)
47 integer(kind=kint), allocatable :: n_ext_per_dom(:), n_int_per_dom(:), ext_item_remote(:)
48 integer(kind=kint), allocatable :: statuses(:,:), requests(:)
49 integer(kind=kint) :: nn_int, np, ilag, icontact, irow, irank, idom, tag, idof, idx, irow_remote
50 integer(kind=kint) :: n_send, is, ie, len
51 if (hecmesh%n_neighbor_pe == 0) then
52 concomm%n_neighbor_pe = 0
53 return
54 endif
55 nn_int = hecmesh%nn_internal
56 np = hecmesh%n_node
57 ! count external contact dofs
58 allocate(n_ext_per_dom(hecmesh%n_neighbor_pe))
59 n_ext_per_dom(:) = 0
60 do ilag = 1, n_contact_dof
61 icontact = contact_dofs(ilag)
62 if (icontact <= nn_int*ndof) cycle ! skip internal contact dof
63 irow = (icontact+ndof-1) / ndof
64 irank = hecmesh%node_ID(2*irow)
65 call rank_to_idom(hecmesh, irank, idom)
66 n_ext_per_dom(idom) = n_ext_per_dom(idom) + 1
67 enddo
68 ! send external / recv internal contact dofs
69 allocate(statuses(hecmw_status_size, hecmesh%n_neighbor_pe))
70 allocate(requests(hecmesh%n_neighbor_pe))
71 do idom = 1, hecmesh%n_neighbor_pe
72 irank = hecmesh%neighbor_pe(idom)
73 tag = 5001
74 call hecmw_isend_int(n_ext_per_dom(idom), 1, irank, tag, &
75 hecmesh%MPI_COMM, requests(idom))
76 enddo
77 allocate(n_int_per_dom(hecmesh%n_neighbor_pe))
78 do idom = 1, hecmesh%n_neighbor_pe
79 irank = hecmesh%neighbor_pe(idom)
80 tag = 5001
81 call hecmw_recv_int(n_int_per_dom(idom), 1, irank, tag, &
82 hecmesh%MPI_COMM, statuses(:,1))
83 enddo
84 call hecmw_waitall(hecmesh%n_neighbor_pe, requests, statuses)
85 ! make index
86 allocate(ext_index(0:hecmesh%n_neighbor_pe))
87 allocate(int_index(0:hecmesh%n_neighbor_pe))
88 ext_index(0) = 0
89 int_index(0) = 0
90 do idom = 1, hecmesh%n_neighbor_pe
91 ext_index(idom) = ext_index(idom-1) + n_ext_per_dom(idom)
92 int_index(idom) = int_index(idom-1) + n_int_per_dom(idom)
93 enddo
94 ! make ext_item
95 allocate(ext_item(ext_index(hecmesh%n_neighbor_pe)))
96 allocate(ext_item_remote(ext_index(hecmesh%n_neighbor_pe)))
97 n_ext_per_dom(:) = 0
98 do ilag = 1, n_contact_dof
99 icontact = contact_dofs(ilag)
100 if (icontact <= nn_int*ndof) cycle ! skip internal contact dof
101 irow = (icontact+ndof-1) / ndof
102 idof = icontact - ndof*(irow-1)
103 irank = hecmesh%node_ID(2*irow)
104 call rank_to_idom(hecmesh, irank, idom)
105 n_ext_per_dom(idom) = n_ext_per_dom(idom) + 1
106 idx = ext_index(idom-1)+n_ext_per_dom(idom)
107 ext_item(idx) = icontact
108 irow_remote = hecmesh%node_ID(2*irow-1)
109 ext_item_remote(idx) = ndof*(irow_remote-1)+idof
110 enddo
111 deallocate(n_ext_per_dom)
112 deallocate(n_int_per_dom)
113 ! send ext_item_remote and recv int_item
114 n_send = 0
115 do idom = 1, hecmesh%n_neighbor_pe
116 irank = hecmesh%neighbor_pe(idom)
117 is = ext_index(idom-1)+1
118 ie = ext_index(idom)
119 len = ie-is+1
120 if (len == 0) cycle
121 n_send = n_send + 1
122 tag = 5002
123 call hecmw_isend_int(ext_item_remote(is:ie), len, irank, tag, &
124 hecmesh%MPI_COMM, requests(n_send))
125 enddo
126 allocate(int_item(int_index(hecmesh%n_neighbor_pe)))
127 do idom = 1, hecmesh%n_neighbor_pe
128 irank = hecmesh%neighbor_pe(idom)
129 is = int_index(idom-1)+1
130 ie = int_index(idom)
131 len = ie-is+1
132 if (len == 0) cycle
133 tag = 5002
134 call hecmw_recv_int(int_item(is:ie), len, irank, tag, &
135 hecmesh%MPI_COMM, statuses(:,1))
136 enddo
137 call hecmw_waitall(n_send, requests, statuses)
138 deallocate(statuses, requests)
139 if (debug >= 2) then
140 write(0,*) ' DEBUG2: ext_index',ext_index(:)
141 write(0,*) ' DEBUG2: ext_item',ext_item(:)
142 write(0,*) ' DEBUG2: ext_item_remote',ext_item_remote(:)
143 write(0,*) ' DEBUG2: int_index',int_index(:)
144 write(0,*) ' DEBUG2: int_item',int_item(:)
145 endif
146 deallocate(ext_item_remote)
147 !
148 concomm%n_neighbor_pe = hecmesh%n_neighbor_pe
149 allocate(concomm%neighbor_pe(concomm%n_neighbor_pe))
150 concomm%neighbor_pe(:) = hecmesh%neighbor_pe(:)
151 concomm%MPI_COMM = hecmesh%MPI_COMM
152 concomm%ext_index => ext_index
153 concomm%ext_item => ext_item
154 concomm%int_index => int_index
155 concomm%int_item => int_item
156 end subroutine fstr_contact_comm_init
157
158 subroutine fstr_contact_comm_finalize(conComm)
159 implicit none
160 type (fstrst_contact_comm), intent(inout) :: concomm
161 if (concomm%n_neighbor_pe == 0) return
162 if (associated(concomm%neighbor_pe)) deallocate(concomm%neighbor_pe)
163 if (associated(concomm%ext_index)) deallocate(concomm%ext_index)
164 if (associated(concomm%ext_item)) deallocate(concomm%ext_item)
165 if (associated(concomm%int_index)) deallocate(concomm%int_index)
166 if (associated(concomm%int_item)) deallocate(concomm%int_item)
167 concomm%n_neighbor_pe = 0
168 concomm%MPI_COMM = 0
169 end subroutine fstr_contact_comm_finalize
170
171 subroutine fstr_contact_comm_reduce_r(conComm, vec, op)
172 implicit none
173 type (fstrst_contact_comm), intent(in) :: concomm
174 real(kind=kreal), intent(inout) :: vec(:)
175 integer(kind=kint), intent(in) :: op
176 if (concomm%n_neighbor_pe == 0) return
177 call send_recv_contact_info_r(concomm%n_neighbor_pe, concomm%neighbor_pe, concomm%MPI_COMM, &
178 concomm%ext_index, concomm%ext_item, concomm%int_index, concomm%int_item, vec, op)
179 end subroutine fstr_contact_comm_reduce_r
180
181 subroutine fstr_contact_comm_bcast_r(conComm, vec)
182 implicit none
183 type (fstrst_contact_comm), intent(in) :: concomm
184 real(kind=kreal), intent(inout) :: vec(:)
185 integer(kind=kint) :: op
186 if (concomm%n_neighbor_pe == 0) return
187 op = op_overwrite
188 call send_recv_contact_info_r(concomm%n_neighbor_pe, concomm%neighbor_pe, concomm%MPI_COMM, &
189 concomm%int_index, concomm%int_item, concomm%ext_index, concomm%ext_item, vec, op)
190 end subroutine fstr_contact_comm_bcast_r
191
192 subroutine fstr_contact_comm_reduce_i(conComm, vec, op)
193 implicit none
194 type (fstrst_contact_comm), intent(in) :: concomm
195 integer(kind=kint), intent(inout) :: vec(:)
196 integer(kind=kint), intent(in) :: op
197 if (concomm%n_neighbor_pe == 0) return
198 call send_recv_contact_info_i(concomm%n_neighbor_pe, concomm%neighbor_pe, concomm%MPI_COMM, &
199 concomm%ext_index, concomm%ext_item, concomm%int_index, concomm%int_item, vec, op)
200 end subroutine fstr_contact_comm_reduce_i
201
202 subroutine fstr_contact_comm_bcast_i(conComm, vec)
203 implicit none
204 type (fstrst_contact_comm), intent(in) :: concomm
205 integer(kind=kint), intent(inout) :: vec(:)
206 integer(kind=kint) :: op
207 if (concomm%n_neighbor_pe == 0) return
208 op = op_overwrite
209 call send_recv_contact_info_i(concomm%n_neighbor_pe, concomm%neighbor_pe, concomm%MPI_COMM, &
210 concomm%int_index, concomm%int_item, concomm%ext_index, concomm%ext_item, vec, op)
211 end subroutine fstr_contact_comm_bcast_i
212
213 subroutine fstr_contact_comm_allreduce_r(conComm, vec, op)
214 implicit none
215 type (fstrst_contact_comm), intent(in) :: concomm
216 real(kind=kreal), intent(inout) :: vec(:)
217 integer(kind=kint), intent(in) :: op
218 call fstr_contact_comm_reduce_r(concomm, vec, op)
219 call fstr_contact_comm_bcast_r(concomm, vec)
220 end subroutine fstr_contact_comm_allreduce_r
221
222 subroutine fstr_contact_comm_allreduce_i(conComm, vec, op)
223 implicit none
224 type (fstrst_contact_comm), intent(in) :: concomm
225 integer(kind=kint), intent(inout) :: vec(:)
226 integer(kind=kint), intent(in) :: op
227 call fstr_contact_comm_reduce_i(concomm, vec, op)
228 call fstr_contact_comm_bcast_i(concomm, vec)
229 end subroutine fstr_contact_comm_allreduce_i
230
231 !
232 ! private subroutines
233 !
234
235 subroutine rank_to_idom(hecMESH, rank, idom)
236 implicit none
237 type (hecmwst_local_mesh), intent(in) :: hecmesh
238 integer(kind=kint), intent(in) :: rank
239 integer(kind=kint), intent(out) :: idom
240 integer(kind=kint) :: i
241 do i = 1, hecmesh%n_neighbor_pe
242 if (hecmesh%neighbor_pe(i) == rank) then
243 idom = i
244 return
245 endif
246 enddo
247 stop 'ERROR: exp_rank not found in neighbor_pe'
248 end subroutine rank_to_idom
249
250 subroutine send_recv_contact_info_r(n_neighbor_pe, neighbor_pe, MPI_COMM, &
251 send_index, send_item, recv_index, recv_item, vec, op)
252 implicit none
253 integer(kind=kint), intent(in) :: n_neighbor_pe
254 integer(kind=kint), intent(in) :: neighbor_pe(:)
255 integer(kind=kint), intent(in) :: mpi_comm
256 integer(kind=kint), pointer, intent(in) :: send_index(:), send_item(:), recv_index(:), recv_item(:)
257 real(kind=kreal), intent(inout) :: vec(:)
258 integer(kind=kint), intent(in) :: op
259 real(kind=kreal), allocatable :: send_buf(:), recv_buf(:)
260 integer(kind=kint) :: i, n_send, idom, irank, is, ie, len, tag
261 integer(kind=kint), allocatable :: requests(:), statuses(:,:)
262 if (n_neighbor_pe == 0) return
263 allocate(requests(n_neighbor_pe))
264 allocate(statuses(hecmw_status_size, n_neighbor_pe))
265 allocate(send_buf(send_index(n_neighbor_pe)))
266 allocate(recv_buf(recv_index(n_neighbor_pe)))
267 do i = 1, send_index(n_neighbor_pe)
268 send_buf(i) = vec(send_item(i))
269 enddo
270 n_send = 0
271 do idom = 1, n_neighbor_pe
272 irank = neighbor_pe(idom)
273 is = send_index(idom-1)+1
274 ie = send_index(idom)
275 len = ie-is+1
276 if (len == 0) cycle
277 n_send = n_send + 1
278 tag = 5011
279 call hecmw_isend_r(send_buf(is:ie), len, irank, tag, &
280 mpi_comm, requests(n_send))
281 enddo
282 do idom = 1, n_neighbor_pe
283 irank = neighbor_pe(idom)
284 is = recv_index(idom-1)+1
285 ie = recv_index(idom)
286 len = ie-is+1
287 if (len == 0) cycle
288 tag = 5011
289 call hecmw_recv_r(recv_buf(is:ie), len, irank, tag, &
290 mpi_comm, statuses(:,1))
291 enddo
292 call hecmw_waitall(n_send, requests, statuses)
293 if (op == hecmw_sum) then
294 do i = 1, recv_index(n_neighbor_pe)
295 vec(recv_item(i)) = vec(recv_item(i)) + recv_buf(i)
296 enddo
297 elseif (op == hecmw_prod) then
298 do i = 1, recv_index(n_neighbor_pe)
299 vec(recv_item(i)) = vec(recv_item(i)) * recv_buf(i)
300 enddo
301 elseif (op == hecmw_max) then
302 do i = 1, recv_index(n_neighbor_pe)
303 vec(recv_item(i)) = max(vec(recv_item(i)), recv_buf(i))
304 enddo
305 elseif (op == hecmw_min) then
306 do i = 1, recv_index(n_neighbor_pe)
307 vec(recv_item(i)) = min(vec(recv_item(i)), recv_buf(i))
308 enddo
309 else ! overwrite
310 do i = 1, recv_index(n_neighbor_pe)
311 vec(recv_item(i)) = recv_buf(i)
312 enddo
313 endif
314 deallocate(requests)
315 deallocate(statuses)
316 if (debug >= 2) then
317 write(0,*) ' DEBUG2: send_buf',send_buf(:)
318 write(0,*) ' DEBUG2: recv_buf',recv_buf(:)
319 endif
320 deallocate(send_buf)
321 deallocate(recv_buf)
322 end subroutine send_recv_contact_info_r
323
324 subroutine send_recv_contact_info_i(n_neighbor_pe, neighbor_pe, MPI_COMM, &
325 send_index, send_item, recv_index, recv_item, vec, op)
326 implicit none
327 integer(kind=kint), intent(in) :: n_neighbor_pe
328 integer(kind=kint), intent(in) :: neighbor_pe(:)
329 integer(kind=kint), intent(in) :: mpi_comm
330 integer(kind=kint), pointer, intent(in) :: send_index(:), send_item(:), recv_index(:), recv_item(:)
331 integer(kind=kint), intent(inout) :: vec(:)
332 integer(kind=kint), intent(in) :: op
333 integer(kind=kint), allocatable :: send_buf(:), recv_buf(:)
334 integer(kind=kint) :: i, n_send, idom, irank, is, ie, len, tag
335 integer(kind=kint), allocatable :: requests(:), statuses(:,:)
336 if (n_neighbor_pe == 0) return
337 allocate(requests(n_neighbor_pe))
338 allocate(statuses(hecmw_status_size, n_neighbor_pe))
339 allocate(send_buf(send_index(n_neighbor_pe)))
340 allocate(recv_buf(recv_index(n_neighbor_pe)))
341 do i = 1, send_index(n_neighbor_pe)
342 send_buf(i) = vec(send_item(i))
343 enddo
344 n_send = 0
345 do idom = 1, n_neighbor_pe
346 irank = neighbor_pe(idom)
347 is = send_index(idom-1)+1
348 ie = send_index(idom)
349 len = ie-is+1
350 if (len == 0) cycle
351 n_send = n_send + 1
352 tag = 5011
353 call hecmw_isend_int(send_buf(is:ie), len, irank, tag, &
354 mpi_comm, requests(n_send))
355 enddo
356 do idom = 1, n_neighbor_pe
357 irank = neighbor_pe(idom)
358 is = recv_index(idom-1)+1
359 ie = recv_index(idom)
360 len = ie-is+1
361 if (len == 0) cycle
362 tag = 5011
363 call hecmw_recv_int(recv_buf(is:ie), len, irank, tag, &
364 mpi_comm, statuses(:,1))
365 enddo
366 call hecmw_waitall(n_send, requests, statuses)
367 if (op == hecmw_sum) then
368 do i = 1, recv_index(n_neighbor_pe)
369 vec(recv_item(i)) = vec(recv_item(i)) + recv_buf(i)
370 enddo
371 elseif (op == hecmw_prod) then
372 do i = 1, recv_index(n_neighbor_pe)
373 vec(recv_item(i)) = vec(recv_item(i)) * recv_buf(i)
374 enddo
375 elseif (op == hecmw_max) then
376 do i = 1, recv_index(n_neighbor_pe)
377 vec(recv_item(i)) = max(vec(recv_item(i)), recv_buf(i))
378 enddo
379 elseif (op == hecmw_min) then
380 do i = 1, recv_index(n_neighbor_pe)
381 vec(recv_item(i)) = min(vec(recv_item(i)), recv_buf(i))
382 enddo
383 else ! overwrite
384 do i = 1, recv_index(n_neighbor_pe)
385 vec(recv_item(i)) = recv_buf(i)
386 enddo
387 endif
388 deallocate(requests)
389 deallocate(statuses)
390 if (debug >= 2) then
391 write(0,*) ' DEBUG2: send_buf',send_buf(:)
392 write(0,*) ' DEBUG2: recv_buf',recv_buf(:)
393 endif
394 deallocate(send_buf)
395 deallocate(recv_buf)
396 end subroutine send_recv_contact_info_i
397
398end module m_fstr_contact_comm
Definition: hecmw.f90:6
subroutine, public fstr_contact_comm_allreduce_r(concomm, vec, op)
subroutine, public fstr_contact_comm_reduce_r(concomm, vec, op)
integer(kind=kint), parameter op_overwrite
subroutine, public fstr_contact_comm_reduce_i(concomm, vec, op)
subroutine, public fstr_contact_comm_init(concomm, hecmesh, ndof, n_contact_dof, contact_dofs)
subroutine, public fstr_contact_comm_bcast_i(concomm, vec)
subroutine, public fstr_contact_comm_allreduce_i(concomm, vec, op)
subroutine, public fstr_contact_comm_finalize(concomm)
subroutine, public fstr_contact_comm_bcast_r(concomm, vec)