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