FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
hecmw_solver_SR_33.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 
6 !C
7 !C***
8 !C*** module hecmw_solver_SR_33
9 !C***
10 !C
12  use hecmw_util
13 
14  private
15 
16  public :: hecmw_solve_send_recv_33
19 
20 #ifndef HECMW_SERIAL
21  type async_buf
22  integer(kind=kint ) :: NEIBPETOT = 0
23  integer(kind=kint ), pointer :: STACK_IMPORT(:)
24  integer(kind=kint ), pointer :: NOD_IMPORT (:)
25  real (kind=kreal), pointer :: ws(:)
26  real (kind=kreal), pointer :: wr(:)
27  real (kind=kreal), pointer :: x(:)
28  integer(kind=kint ), pointer :: req1(:)
29  integer(kind=kint ), pointer :: req2(:)
30  integer(kind=kint ) :: nreq1
31  integer(kind=kint ) :: nreq2
32  end type async_buf
33 
34  integer(kind=kint), parameter :: MAX_NREQ = 8
35  type(async_buf), save :: abuf(MAX_NREQ)
36 #endif
37 
38 contains
39  !C
40  !C*** SOLVER_SEND_RECV
41  !C
43  & ( n, neibpetot, neibpe, stack_import, nod_import, &
44  & stack_export, nod_export, &
45  & ws, wr, x, solver_comm,my_rank)
46 
47  implicit none
48 
49  integer(kind=kint ) , intent(in) :: n
50  integer(kind=kint ) , intent(in) :: neibpetot
51  integer(kind=kint ), pointer :: neibpe (:)
52  integer(kind=kint ), pointer :: stack_import(:)
53  integer(kind=kint ), pointer :: nod_import (:)
54  integer(kind=kint ), pointer :: stack_export(:)
55  integer(kind=kint ), pointer :: nod_export (:)
56  real (kind=kreal), dimension(: ), intent(inout):: ws
57  real (kind=kreal), dimension(: ), intent(inout):: wr
58  real (kind=kreal), dimension(: ), intent(inout):: x
59  integer(kind=kint ) , intent(in) ::solver_comm
60  integer(kind=kint ) , intent(in) :: my_rank
61 
62 #ifndef HECMW_SERIAL
63  integer(kind=kint ), dimension(:,:), allocatable :: sta1
64  integer(kind=kint ), dimension(:,:), allocatable :: sta2
65  integer(kind=kint ), dimension(: ), allocatable :: req1
66  integer(kind=kint ), dimension(: ), allocatable :: req2
67 
68  ! local valiables
69  integer(kind=kint ) :: neib,istart,inum,k,ii,ierr,nreq1,nreq2
70  !C
71  !C-- INIT.
72  allocate (sta1(mpi_status_size,neibpetot))
73  allocate (sta2(mpi_status_size,neibpetot))
74  allocate (req1(neibpetot))
75  allocate (req2(neibpetot))
76 
77  !C
78  !C-- SEND
79  nreq1=0
80  do neib= 1, neibpetot
81  istart= stack_export(neib-1)
82  inum = stack_export(neib ) - istart
83  if (inum==0) cycle
84  nreq1=nreq1+1
85  do k= istart+1, istart+inum
86  ii = 3*nod_export(k)
87  ws(3*k-2)= x(ii-2)
88  ws(3*k-1)= x(ii-1)
89  ws(3*k )= x(ii )
90  enddo
91 
92  call mpi_isend (ws(3*istart+1), 3*inum,mpi_double_precision, &
93  & neibpe(neib), 0, solver_comm, req1(nreq1), ierr)
94  enddo
95 
96  !C
97  !C-- RECEIVE
98  nreq2=0
99  do neib= 1, neibpetot
100  istart= stack_import(neib-1)
101  inum = stack_import(neib ) - istart
102  if (inum==0) cycle
103  nreq2=nreq2+1
104  call mpi_irecv (wr(3*istart+1), 3*inum, mpi_double_precision, &
105  & neibpe(neib), 0, solver_comm, req2(nreq2), ierr)
106  enddo
107 
108  call mpi_waitall (nreq2, req2, sta2, ierr)
109 
110  do neib= 1, neibpetot
111  istart= stack_import(neib-1)
112  inum = stack_import(neib ) - istart
113  do k= istart+1, istart+inum
114  ii = 3*nod_import(k)
115  x(ii-2)= wr(3*k-2)
116  x(ii-1)= wr(3*k-1)
117  x(ii )= wr(3*k )
118  enddo
119  enddo
120 
121  call mpi_waitall (nreq1, req1, sta1, ierr)
122  deallocate (sta1, sta2, req1, req2)
123 #endif
124  end subroutine hecmw_solve_send_recv_33
125 
126  !C
127  !C*** SOLVER_ISEND_IRECV
128  !C
130  & ( n, neibpetot, neibpe, stack_import, nod_import, &
131  & stack_export, nod_export, &
132  & x, solver_comm,my_rank,ireq)
133  implicit none
134  integer(kind=kint ) , intent(in) :: n
135  integer(kind=kint ) , intent(in) :: neibpetot
136  integer(kind=kint ), pointer :: neibpe (:)
137  integer(kind=kint ), pointer :: stack_import(:)
138  integer(kind=kint ), pointer :: nod_import (:)
139  integer(kind=kint ), pointer :: stack_export(:)
140  integer(kind=kint ), pointer :: nod_export (:)
141  real (kind=kreal), target, intent(inout):: x(:)
142  integer(kind=kint ) , intent(in) ::solver_comm
143  integer(kind=kint ) , intent(in) :: my_rank
144  integer(kind=kint ) , intent(out) :: ireq
145 
146 #ifndef HECMW_SERIAL
147  ! local valiables
148  real (kind=kreal), pointer :: ws(:)
149  real (kind=kreal), pointer :: wr(:)
150  integer(kind=kint ), pointer :: req1(:)
151  integer(kind=kint ), pointer :: req2(:)
152  integer(kind=kint ) :: neib,istart,inum,k,ii,ierr,i,nreq1,nreq2
153  !C
154  !C-- INIT.
155  allocate (ws(3*stack_export(neibpetot)))
156  allocate (wr(3*stack_import(neibpetot)))
157  allocate (req1(neibpetot))
158  allocate (req2(neibpetot))
159  !C
160  !C-- SEND
161  nreq1=0
162  do neib= 1, neibpetot
163  istart= stack_export(neib-1)
164  inum = stack_export(neib ) - istart
165  if (inum==0) cycle
166  nreq1=nreq1+1
167  do k= istart+1, istart+inum
168  ii = 3*nod_export(k)
169  ws(3*k-2)= x(ii-2)
170  ws(3*k-1)= x(ii-1)
171  ws(3*k )= x(ii )
172  enddo
173  call mpi_isend (ws(3*istart+1), 3*inum,mpi_double_precision, &
174  & neibpe(neib), 0, solver_comm, req1(nreq1), ierr)
175  enddo
176  !C
177  !C-- RECEIVE
178  nreq2=0
179  do neib= 1, neibpetot
180  istart= stack_import(neib-1)
181  inum = stack_import(neib ) - istart
182  if (inum==0) cycle
183  nreq2=nreq2+1
184  call mpi_irecv (wr(3*istart+1), 3*inum, mpi_double_precision, &
185  & neibpe(neib), 0, solver_comm, req2(nreq2), ierr)
186  enddo
187  !C
188  !C-- Find empty abuf
189  ireq = 0
190  do i = 1, max_nreq
191  if (abuf(i)%NEIBPETOT == 0) then
192  ireq = i
193  exit
194  endif
195  enddo
196  if (ireq == 0) then
197  stop 'Error: hecmw_solve_isend_irecv_33: exceeded maximum num of requests'
198  endif
199  !C
200  !C-- Save in abuf
201  abuf(ireq)%NEIBPETOT = neibpetot
202  abuf(ireq)%STACK_IMPORT=> stack_import
203  abuf(ireq)%NOD_IMPORT => nod_import
204  abuf(ireq)%WS => ws
205  abuf(ireq)%WR => wr
206  abuf(ireq)%X => x
207  abuf(ireq)%req1 => req1
208  abuf(ireq)%req2 => req2
209  abuf(ireq)%nreq1 = nreq1
210  abuf(ireq)%nreq2 = nreq2
211 #endif
212  end subroutine hecmw_solve_isend_irecv_33
213 
214  !C
215  !C*** SOLVER_ISEND_IRECV_WAIT
216  !C
218  implicit none
219  integer(kind=kint ), intent(in) :: ireq
220 
221 #ifndef HECMW_SERIAL
222  ! local valiables
223  integer(kind=kint ) :: neibpetot
224  integer(kind=kint ), pointer :: stack_import(:)
225  integer(kind=kint ), pointer :: nod_import (:)
226  real (kind=kreal), pointer :: ws(:)
227  real (kind=kreal), pointer :: wr(:)
228  real (kind=kreal), pointer :: x(:)
229  integer(kind=kint ), pointer :: req1(:)
230  integer(kind=kint ), pointer :: req2(:)
231  integer(kind=kint ), dimension(:,:), allocatable :: sta1
232  integer(kind=kint ), dimension(:,:), allocatable :: sta2
233  integer(kind=kint ) :: neib,istart,inum,k,ii,ierr,nreq1,nreq2
234  !C-- Check ireq
235  if (ireq < 0 .or. ireq > max_nreq) then
236  stop 'ERROR: hecmw_solve_isend_irecv_33_wait: invalid ireq'
237  endif
238  !C-- Restore from abuf
239  neibpetot = abuf(ireq)%NEIBPETOT
240  stack_import=> abuf(ireq)%STACK_IMPORT
241  nod_import => abuf(ireq)%NOD_IMPORT
242  ws => abuf(ireq)%WS
243  wr => abuf(ireq)%WR
244  x => abuf(ireq)%X
245  req1 => abuf(ireq)%req1
246  req2 => abuf(ireq)%req2
247  nreq1 = abuf(ireq)%nreq1
248  nreq2 = abuf(ireq)%nreq2
249  !C-- Free abuf
250  abuf(ireq)%NEIBPETOT = 0
251  !C-- Allocate
252  allocate (sta1(mpi_status_size,neibpetot))
253  allocate (sta2(mpi_status_size,neibpetot))
254  !C-- Wait irecv
255  call mpi_waitall (nreq2, req2, sta2, ierr)
256  do neib= 1, neibpetot
257  istart= stack_import(neib-1)
258  inum = stack_import(neib ) - istart
259  do k= istart+1, istart+inum
260  ii = 3*nod_import(k)
261  x(ii-2)= wr(3*k-2)
262  x(ii-1)= wr(3*k-1)
263  x(ii )= wr(3*k )
264  enddo
265  enddo
266  !C-- Wait isend
267  call mpi_waitall (nreq1, req1, sta1, ierr)
268  !C-- Deallocate
269  deallocate (sta1, sta2)
270  deallocate (req1, req2)
271  deallocate (ws, wr)
272 #endif
273  end subroutine hecmw_solve_isend_irecv_33_wait
274 
275 end module hecmw_solver_sr_33
subroutine, public hecmw_solve_isend_irecv_33_wait(ireq)
subroutine, public hecmw_solve_send_recv_33(N, NEIBPETOT, NEIBPE, STACK_IMPORT, NOD_IMPORT, STACK_EXPORT, NOD_EXPORT, WS, WR, X, SOLVER_COMM, my_rank)
subroutine, public hecmw_solve_isend_irecv_33(N, NEIBPETOT, NEIBPE, STACK_IMPORT, NOD_IMPORT, STACK_EXPORT, NOD_EXPORT, X, SOLVER_COMM, my_rank, ireq)
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal