FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
hecmw_solver_SR_66.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_66
9 !C***
10 !C
12 contains
13  !C
14  !C*** SOLVER_SEND_RECV
15  !C
17  & ( n, neibpetot, neibpe, stack_import, nod_import, &
18  & stack_export, nod_export, &
19  & ws, wr, x, solver_comm,my_rank)
20 
21  use hecmw_util
22  implicit none
23  ! include 'mpif.h'
24  ! include 'hecmw_config_f.h'
25 
26  integer(kind=kint ) , intent(in) :: N
27  integer(kind=kint ) , intent(in) :: NEIBPETOT
28  integer(kind=kint ), pointer :: NEIBPE (:)
29  integer(kind=kint ), pointer :: STACK_IMPORT(:)
30  integer(kind=kint ), pointer :: NOD_IMPORT (:)
31  integer(kind=kint ), pointer :: STACK_EXPORT(:)
32  integer(kind=kint ), pointer :: NOD_EXPORT (:)
33  real (kind=kreal), dimension(: ), intent(inout):: ws
34  real (kind=kreal), dimension(: ), intent(inout):: wr
35  real (kind=kreal), dimension(: ), intent(inout):: x
36  integer(kind=kint ) , intent(in) ::SOLVER_COMM
37  integer(kind=kint ) , intent(in) :: my_rank
38 
39 #ifndef HECMW_SERIAL
40  integer(kind=kint ), dimension(:,:), allocatable :: sta1
41  integer(kind=kint ), dimension(:,:), allocatable :: sta2
42  integer(kind=kint ), dimension(: ), allocatable :: req1
43  integer(kind=kint ), dimension(: ), allocatable :: req2
44 
45  integer(kind=kint ), save :: NFLAG
46  data nflag/0/
47 
48  ! local valiables
49  integer(kind=kint ) :: neib,istart,inum,k,ii,ierr,nreq1,nreq2
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  nreq1=0
60  do neib= 1, neibpetot
61  istart= stack_export(neib-1)
62  inum = stack_export(neib ) - istart
63  if (inum==0) cycle
64  nreq1=nreq1+1
65  do k= istart+1, istart+inum
66  ii = 6*nod_export(k)
67  ws(6*k-5)= x(ii-5)
68  ws(6*k-4)= x(ii-4)
69  ws(6*k-3)= x(ii-3)
70  ws(6*k-2)= x(ii-2)
71  ws(6*k-1)= x(ii-1)
72  ws(6*k )= x(ii )
73  enddo
74 
75  call mpi_isend (ws(6*istart+1), 6*inum,mpi_double_precision, &
76  & neibpe(neib), 0, solver_comm, req1(nreq1), ierr)
77  enddo
78 
79  !C
80  !C-- RECEIVE
81  nreq2=0
82  do neib= 1, neibpetot
83  istart= stack_import(neib-1)
84  inum = stack_import(neib ) - istart
85  if (inum==0) cycle
86  nreq2=nreq2+1
87  call mpi_irecv (wr(6*istart+1), 6*inum, mpi_double_precision, &
88  & neibpe(neib), 0, solver_comm, req2(nreq2), ierr)
89  enddo
90 
91  call mpi_waitall (nreq2, req2, sta2, ierr)
92 
93  do neib= 1, neibpetot
94  istart= stack_import(neib-1)
95  inum = stack_import(neib ) - istart
96  do k= istart+1, istart+inum
97  ii = 6*nod_import(k)
98  x(ii-5)= wr(6*k-5)
99  x(ii-4)= wr(6*k-4)
100  x(ii-3)= wr(6*k-3)
101  x(ii-2)= wr(6*k-2)
102  x(ii-1)= wr(6*k-1)
103  x(ii )= wr(6*k )
104  enddo
105  enddo
106 
107  call mpi_waitall (nreq1, req1, sta1, ierr)
108  deallocate (sta1, sta2, req1, req2)
109 #endif
110  end subroutine hecmw_solve_send_recv_66
111 end module hecmw_solver_sr_66
subroutine hecmw_solve_send_recv_66(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