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