1 !C
2 !C*** SOLVER_SEND_RECV
3 !C
4 subroutine HECMOLVE_SEND_RECV_33 &
5 & ( NW_S, NEIBPETOT, NEIBPE, STACK_IMPORT, NOD_IMPORT, &
6 & STACK_EXPORT, NOD_EXPORT, &
7 & WS, WR, X, SOLVER_COMM,my_rank)
8
9 implicit none
10
11 integer(kind=kint ) , intent(in) :: N
12 integer(kind=kint ) , intent(in) :: NEIBPETOT
13 integer(kind=kint ), pointer :: NEIBPE (:)
14 integer(kind=kint ), pointer :: STACK_IMPORT(:)
15 integer(kind=kint ), pointer :: NOD_IMPORT (:)
16 integer(kind=kint ), pointer :: STACK_EXPORT(:)
17 integer(kind=kint ), pointer :: NOD_EXPORT (:)
18 real (kind=kreal), dimension(: ), intent(inout):: WS
19 real (kind=kreal), dimension(: ), intent(inout):: WR
20 real (kind=kreal), dimension(: ), intent(inout):: X
21 integer(kind=kint ) , intent(in) ::SOLVER_COMM
22 integer(kind=kint ) , intent(in) :: my_rank
23
24 #ifndef HECMW_SERIAL
25 integer(kind=kint ), dimension(:,:), allocatable :: sta1
26 integer(kind=kint ), dimension(:,:), allocatable :: sta2
27 integer(kind=kint ), dimension(: ), allocatable :: req1
28 integer(kind=kint ), dimension(: ), allocatable :: req2
29
30 ! local valiables
31 integer(kind=kint ) :: neib,istart,inum,k,ii,ierr
32 !C
33 !C-- INIT.
34 allocate (sta1(MPI_STATUS_SIZE,NEIBPETOT))
35 allocate (sta2(MPI_STATUS_SIZE,NEIBPETOT))
36 allocate (req1(NEIBPETOT))
37 allocate (req2(NEIBPETOT))
38
39 !C
40 !C-- SEND
41 do neib= 1, NEIBPETOT
42 istart= STACK_EXPORT(neib-1)
43 inum = STACK_EXPORT(neib ) - istart
44 do k= istart+1, istart+inum
45 ii = 3*NOD_EXPORT(k)
46 WS(3*k-2)= X(ii-2)
47 WS(3*k-1)= X(ii-1)
48 WS(3*k )= X(ii )
49 enddo
50
51 call MPI_ISEND (WS(3*istart+1), 3*inum,MPI_DOUBLE_PRECISION, &
52 & NEIBPE(neib), 0, SOLVER_COMM, req1(neib), ierr)
53 enddo
54
55 !C
56 !C-- RECEIVE
57 do neib= 1, NEIBPETOT
58 istart= STACK_IMPORT(neib-1)
59 inum = STACK_IMPORT(neib ) - istart
60 call MPI_IRECV (WR(3*istart+1), 3*inum, MPI_DOUBLE_PRECISION, &
61 & NEIBPE(neib), 0, SOLVER_COMM, req2(neib), ierr)
62 enddo
63
64 call MPI_WAITALL (NEIBPETOT, req2, sta2, ierr)
65
66 do neib= 1, NEIBPETOT
67 istart= STACK_IMPORT(neib-1)
68 inum = STACK_IMPORT(neib ) - istart
69 do k= istart+1, istart+inum
70 ii = 3*NOD_IMPORT(k)
71 X(ii-2)= WR(3*k-2)
72 X(ii-1)= WR(3*k-1)
73 X(ii )= WR(3*k )
74 enddo
75 enddo
76
77 call MPI_WAITALL (NEIBPETOT, req1, sta1, ierr)
78 deallocate (sta1, sta2, req1, req2)
79 #endif
80 end subroutine hecmw_solve_send_recv_33