FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_adapt_real_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*** REAL_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 real (kind=kreal), dimension(NB*m), intent(inout):: ws
30 real (kind=kreal), dimension(NB*m), intent(inout):: wr
31 real (kind=kreal), 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
40 integer(kind=kint ), save :: NFLAG
41 data nflag/0/
42
43 !C
44 !C-- INIT.
45 if (nflag.eq.0) then
46 allocate (sta1(mpi_status_size,neibpetot))
47 allocate (sta2(mpi_status_size,neibpetot))
48 allocate (req1(neibpetot))
49 allocate (req2(neibpetot))
50 nflag= 1
51 endif
52
53 !C
54 !C-- SEND
55 do neib= 1, neibpetot
56 istart= stack_export(neib-1)
57 inum = stack_export(neib ) - istart
58
59 do k= istart+1, istart+inum
60 ii= nb*nod_export(k) - nb
61 ik= nb*k - nb
62 do j= 1, nb
63 ws(ik+j)= x(ii+j)
64 enddo
65 enddo
66 call mpi_isend (ws(nb*istart+1), nb*inum,mpi_double_precision, &
67 & neibpe(neib), 0, solver_comm, req1(neib), ierr)
68 enddo
69
70 !C
71 !C-- RECEIVE
72 do neib= 1, neibpetot
73 istart= stack_import(neib-1)
74 inum = stack_import(neib ) - istart
75 call mpi_irecv (wr(nb*istart+1), nb*inum, mpi_double_precision, &
76 & neibpe(neib), 0, solver_comm, req2(neib), ierr)
77 enddo
78
79 call mpi_waitall (neibpetot, req2, sta2, ierr)
80
81 do neib= 1, neibpetot
82 istart= stack_import(neib-1)
83 inum = stack_import(neib ) - istart
84 do k= istart+1, istart+inum
85 ii= nb*nod_import(k) - nb
86 ik= nb*k - nb
87 do j= 1, nb
88 x(ii+j)= wr(ik+j)
89 enddo
90 enddo
91 enddo
92
93 call mpi_waitall (neibpetot, req1, sta1, ierr)
94
95 end subroutine hecmw_adapt_real_send_recv
96end module hecmw_adapt_real_sr
97
98
99
Adaptive Mesh Refinement.
subroutine hecmw_adapt_real_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
integer(kind=4), parameter kreal