FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_adapt_int_sr.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!-------------------------------------------------------------------------------
6
8contains
9 !C
10 !C***
11 !C*** hecmw_adapt_INT_SEND_RECV
12 !C***
13 !C
15 & ( n, neibpetot,neibpe,stack_import, nod_import, &
16 & stack_export, nod_export, &
17 & ws, wr, x, solver_comm,my_rank, nb, m)
18
19 use hecmw_util
20 implicit real*8 (a-h,o-z)
21
22 integer(kind=kint) , intent(in) :: N, m
23 integer(kind=kint) , intent(in) :: NEIBPETOT
24 integer(kind=kint), pointer :: NEIBPE (:)
25 integer(kind=kint), pointer :: STACK_IMPORT(:)
26 integer(kind=kint), pointer :: NOD_IMPORT (:)
27 integer(kind=kint), pointer :: STACK_EXPORT(:)
28 integer(kind=kint), pointer :: NOD_EXPORT (:)
29 integer(kind=kint), dimension(NB*m), intent(inout):: WS
30 integer(kind=kint), dimension(NB*m), intent(inout):: WR
31 integer(kind=kint), dimension(NB*N), intent(inout):: X
32 integer(kind=kint) , intent(in) ::SOLVER_COMM
33 integer(kind=kint) , intent(in) :: my_rank
34
35 integer(kind=kint ), dimension(:,:), save, allocatable :: sta1
36 integer(kind=kint ), dimension(:,:), save, allocatable :: sta2
37 integer(kind=kint ), dimension(: ), save, allocatable :: req1
38 integer(kind=kint ), dimension(: ), save, allocatable :: req2
39 integer(kind=kint ), save :: NFLAG
40 data nflag/0/
41
42 !C
43 !C-- INIT.
44 if (nflag.eq.0) then
45 allocate (sta1(mpi_status_size,neibpetot))
46 allocate (sta2(mpi_status_size,neibpetot))
47 allocate (req1(neibpetot))
48 allocate (req2(neibpetot))
49 nflag= 1
50 endif
51
52 !C
53 !C-- SEND
54 do neib= 1, neibpetot
55 istart= stack_export(neib-1)
56 inum = stack_export(neib ) - istart
57
58 do k= istart+1, istart+inum
59 ii= nb*nod_export(k) - nb
60 ik= nb*k - nb
61 do j= 1, nb
62 ws(ik+j)= x(ii+j)
63 enddo
64 enddo
65 call mpi_isend (ws(nb*istart+1), nb*inum, mpi_integer, &
66 & neibpe(neib), 0, solver_comm, req1(neib), ierr)
67 enddo
68
69 !C
70 !C-- RECEIVE
71 do neib= 1, neibpetot
72 istart= stack_import(neib-1)
73 inum = stack_import(neib ) - istart
74 call mpi_irecv (wr(nb*istart+1), nb*inum, mpi_integer, &
75 & neibpe(neib), 0, solver_comm, req2(neib), ierr)
76 enddo
77
78 call mpi_waitall (neibpetot, req2, sta2, ierr)
79
80 do neib= 1, neibpetot
81 istart= stack_import(neib-1)
82 inum = stack_import(neib ) - istart
83 do k= istart+1, istart+inum
84 ii= nb*nod_import(k) - nb
85 ik= nb*k - nb
86 do j= 1, nb
87 x(ii+j)= wr(ik+j)
88 enddo
89 enddo
90 enddo
91
92 call mpi_waitall (neibpetot, req1, sta1, ierr)
93
94 end subroutine hecmw_adapt_int_send_recv
95end module hecmw_adapt_int_sr
96
97
98
Adaptive Mesh Refinement.
subroutine hecmw_adapt_int_send_recv(n, neibpetot, neibpe, stack_import, nod_import, stack_export, nod_export, ws, wr, x, solver_comm, my_rank, nb, m)
I/O and Utility.
Definition: hecmw_util_f.F90:7