[2] | 1 | ! |
---|
| 2 | ! $Id: MessPass3D.F,v 1.2 2003/12/17 13:55:58 pmarches Exp $ |
---|
| 3 | ! |
---|
| 4 | #include "cppdefs.h" |
---|
| 5 | #if defined MPI && defined SOLVE3D |
---|
| 6 | |
---|
| 7 | subroutine MessPass3D_tile (Istr,Iend,Jstr,Jend, A, nmax) |
---|
| 8 | ! |
---|
| 9 | ! This subroutine is designed for ROMS-MPI code. It exchanges domain |
---|
| 10 | ! boundary information, including 2 ghost-cells in each direction. |
---|
| 11 | ! Ping Wang 9/20/99. |
---|
| 12 | ! |
---|
| 13 | ! implicit none |
---|
| 14 | # include "param.h" |
---|
| 15 | # include "scalars.h" |
---|
| 16 | include 'mpif.h' |
---|
| 17 | |
---|
| 18 | integer nmax |
---|
| 19 | real A(GLOBAL_2D_ARRAY,nmax) |
---|
| 20 | CSDISTRIBUTE_RESHAPE A(BLOCK_PATTERN) BLOCK_CLAUSE |
---|
| 21 | integer Istr,Iend,Jstr,Jend, i,j,k, isize,jsize,ksize, |
---|
| 22 | & req(8), status(MPI_STATUS_SIZE,8), ierr |
---|
| 23 | |
---|
| 24 | integer sub_X,size_X, sub_E,size_E, size_Z |
---|
| 25 | #if ! defined AGRIF |
---|
| 26 | parameter (size_Z=4*(N+1), |
---|
| 27 | & sub_X=(Lm+NSUB_X-1)/NSUB_X, size_X=2*(N+1)*(4+sub_X), |
---|
| 28 | & sub_E=(Mm+NSUB_E-1)/NSUB_E, size_E=2*(N+1)*(4+sub_E)) |
---|
| 29 | |
---|
| 30 | real buf_snd4(size_Z), ibuf_sndN(size_X), buf_snd2(size_Z), |
---|
| 31 | & buf_rev4(size_Z), ibuf_revN(size_X), buf_rev2(size_Z), |
---|
| 32 | & jbuf_sndW(size_E), jbuf_sndE(size_E), |
---|
| 33 | & jbuf_revW(size_E), jbuf_revE(size_E), |
---|
| 34 | & buf_snd1(size_Z), ibuf_sndS(size_X), buf_snd3(size_Z), |
---|
| 35 | & buf_rev1(size_Z), ibuf_revS(size_X), buf_rev3(size_Z) |
---|
| 36 | #else |
---|
| 37 | real, dimension(:), allocatable :: |
---|
| 38 | & buf_snd4, ibuf_sndN, buf_snd2, |
---|
| 39 | & buf_rev4, ibuf_revN, buf_rev2, |
---|
| 40 | & jbuf_sndW, jbuf_sndE, |
---|
| 41 | & jbuf_revW, jbuf_revE, |
---|
| 42 | & buf_snd1, ibuf_sndS, buf_snd3, |
---|
| 43 | & buf_rev1, ibuf_revS, buf_rev3 |
---|
| 44 | #endif |
---|
| 45 | |
---|
| 46 | c** |
---|
| 47 | c common /buffers_3D/ |
---|
| 48 | c & buf_snd4, ibuf_sndN, buf_snd2, |
---|
| 49 | c & buf_rev4, ibuf_revN, buf_rev2, |
---|
| 50 | c |
---|
| 51 | c & jbuf_sndW, jbuf_sndE, |
---|
| 52 | c & jbuf_revW, jbuf_revE, |
---|
| 53 | c |
---|
| 54 | c & buf_snd1, ibuf_sndS, buf_snd3, |
---|
| 55 | c & buf_rev1, ibuf_revS, buf_rev3 |
---|
| 56 | c** |
---|
| 57 | ! |
---|
| 58 | #include "compute_message_bounds.h" |
---|
| 59 | #if defined AGRIF |
---|
| 60 | size_Z=4*(N+1) |
---|
| 61 | sub_X=(Lm+NSUB_X-1)/NSUB_X |
---|
| 62 | size_X=2*(N+1)*(4+sub_X) |
---|
| 63 | sub_E=(Mm+NSUB_E-1)/NSUB_E |
---|
| 64 | size_E=2*(N+1)*(4+sub_E) |
---|
| 65 | |
---|
| 66 | Allocate(buf_snd4(size_Z), ibuf_sndN(size_X), buf_snd2(size_Z), |
---|
| 67 | & buf_rev4(size_Z), ibuf_revN(size_X), buf_rev2(size_Z), |
---|
| 68 | & jbuf_sndW(size_E), jbuf_sndE(size_E), |
---|
| 69 | & jbuf_revW(size_E), jbuf_revE(size_E), |
---|
| 70 | & buf_snd1(size_Z), ibuf_sndS(size_X), buf_snd3(size_Z), |
---|
| 71 | & buf_rev1(size_Z), ibuf_revS(size_X), buf_rev3(size_Z)) |
---|
| 72 | #endif |
---|
| 73 | ! |
---|
| 74 | ksize=4*nmax ! message sizes for |
---|
| 75 | isize=2*ishft*nmax ! corner messages and sides |
---|
| 76 | jsize=2*jshft*nmax ! in XI and ETA directions |
---|
| 77 | ! |
---|
| 78 | ! Prepare to receive and send: sides.... |
---|
| 79 | ! |
---|
| 80 | |
---|
| 81 | if (WEST_INTER) then |
---|
| 82 | do k=1,nmax |
---|
| 83 | do j=jmin,jmax |
---|
| 84 | jbuf_sndW(k+nmax*(j-jmin ))=A(1,j,k) |
---|
| 85 | jbuf_sndW(k+nmax*(j-jmin+jshft))=A(2,j,k) |
---|
| 86 | enddo |
---|
| 87 | enddo |
---|
| 88 | call MPI_Irecv (jbuf_revW, jsize, MPI_DOUBLE_PRECISION, |
---|
| 89 | & p_W, 2, MPI_COMM_WORLD, req(1), ierr) |
---|
| 90 | call MPI_Send (jbuf_sndW, jsize, MPI_DOUBLE_PRECISION, |
---|
| 91 | & p_W, 1, MPI_COMM_WORLD, ierr) |
---|
| 92 | endif |
---|
| 93 | |
---|
| 94 | if (EAST_INTER) then |
---|
| 95 | do k=1,nmax |
---|
| 96 | do j=jmin,jmax |
---|
| 97 | jbuf_sndE(k+nmax*(j-jmin ))=A(Lmmpi-1,j,k) |
---|
| 98 | jbuf_sndE(k+nmax*(j-jmin+jshft))=A(Lmmpi ,j,k) |
---|
| 99 | enddo |
---|
| 100 | enddo |
---|
| 101 | call MPI_Irecv (jbuf_revE, jsize, MPI_DOUBLE_PRECISION, |
---|
| 102 | & p_E, 1, MPI_COMM_WORLD, req(2), ierr) |
---|
| 103 | call MPI_Send (jbuf_sndE, jsize, MPI_DOUBLE_PRECISION, |
---|
| 104 | & p_E, 2, MPI_COMM_WORLD, ierr) |
---|
| 105 | endif |
---|
| 106 | |
---|
| 107 | if (SOUTH_INTER) then |
---|
| 108 | ibuf_snds = 0. |
---|
| 109 | |
---|
| 110 | do k=1,nmax |
---|
| 111 | do i=imin,imax |
---|
| 112 | ibuf_sndS(k+nmax*(i-imin ))=A(i,1,k) |
---|
| 113 | ibuf_sndS(k+nmax*(i-imin+ishft))=A(i,2,k) |
---|
| 114 | enddo |
---|
| 115 | enddo |
---|
| 116 | call MPI_Irecv (ibuf_revS, isize, MPI_DOUBLE_PRECISION, |
---|
| 117 | & p_S, 4, MPI_COMM_WORLD, req(3), ierr) |
---|
| 118 | call MPI_Send (ibuf_sndS, isize, MPI_DOUBLE_PRECISION, |
---|
| 119 | & p_S, 3, MPI_COMM_WORLD, ierr) |
---|
| 120 | endif |
---|
| 121 | |
---|
| 122 | if (NORTH_INTER) then |
---|
| 123 | ibuf_sndn = 0. |
---|
| 124 | do k=1,nmax |
---|
| 125 | do i=imin,imax |
---|
| 126 | ibuf_sndN(k+nmax*(i-imin ))=A(i,Mmmpi-1,k) |
---|
| 127 | ibuf_sndN(k+nmax*(i-imin+ishft))=A(i,Mmmpi ,k) |
---|
| 128 | enddo |
---|
| 129 | enddo |
---|
| 130 | |
---|
| 131 | call MPI_Irecv (ibuf_revN, isize, MPI_DOUBLE_PRECISION, |
---|
| 132 | & p_N, 3, MPI_COMM_WORLD, req(4), ierr) |
---|
| 133 | call MPI_Send (ibuf_sndN, isize, MPI_DOUBLE_PRECISION, |
---|
| 134 | & p_N, 4, MPI_COMM_WORLD, ierr) |
---|
| 135 | endif |
---|
| 136 | |
---|
| 137 | ! |
---|
| 138 | ! ...corners: |
---|
| 139 | ! |
---|
| 140 | if (SOUTH_INTER .and. WEST_INTER) then |
---|
| 141 | do k=1,nmax |
---|
| 142 | buf_snd1(k )=A(1,1,k) |
---|
| 143 | buf_snd1(k+nmax )=A(2,1,k) |
---|
| 144 | buf_snd1(k+2*nmax)=A(1,2,k) |
---|
| 145 | buf_snd1(k+3*nmax)=A(2,2,k) |
---|
| 146 | enddo |
---|
| 147 | call MPI_Irecv (buf_rev1, ksize, MPI_DOUBLE_PRECISION, p_SW, |
---|
| 148 | & 6, MPI_COMM_WORLD, req(5),ierr) |
---|
| 149 | call MPI_Send (buf_snd1, ksize, MPI_DOUBLE_PRECISION, p_SW, |
---|
| 150 | & 5, MPI_COMM_WORLD, ierr) |
---|
| 151 | endif |
---|
| 152 | |
---|
| 153 | if (NORTH_INTER .and. EAST_INTER) then |
---|
| 154 | do k=1,nmax |
---|
| 155 | buf_snd2(k )=A(Lmmpi-1,Mmmpi-1,k) |
---|
| 156 | buf_snd2(k+nmax )=A(Lmmpi ,Mmmpi-1,k) |
---|
| 157 | buf_snd2(k+2*nmax)=A(Lmmpi-1,Mmmpi ,k) |
---|
| 158 | buf_snd2(k+3*nmax)=A(Lmmpi ,Mmmpi ,k) |
---|
| 159 | enddo |
---|
| 160 | call MPI_Irecv (buf_rev2, ksize, MPI_DOUBLE_PRECISION, p_NE, |
---|
| 161 | & 5, MPI_COMM_WORLD, req(6),ierr) |
---|
| 162 | call MPI_Send (buf_snd2, ksize, MPI_DOUBLE_PRECISION, p_NE, |
---|
| 163 | & 6, MPI_COMM_WORLD, ierr) |
---|
| 164 | endif |
---|
| 165 | |
---|
| 166 | if (SOUTH_INTER .and. EAST_INTER) then |
---|
| 167 | do k=1,nmax |
---|
| 168 | buf_snd3(k )=A(Lmmpi-1,1,k) |
---|
| 169 | buf_snd3(k+nmax )=A(Lmmpi ,1,k) |
---|
| 170 | buf_snd3(k+2*nmax)=A(Lmmpi-1,2,k) |
---|
| 171 | buf_snd3(k+3*nmax)=A(Lmmpi ,2,k) |
---|
| 172 | enddo |
---|
| 173 | call MPI_Irecv (buf_rev3, ksize, MPI_DOUBLE_PRECISION, p_SE, |
---|
| 174 | & 8, MPI_COMM_WORLD, req(7),ierr) |
---|
| 175 | call MPI_Send (buf_snd3, ksize, MPI_DOUBLE_PRECISION, p_SE, |
---|
| 176 | & 7, MPI_COMM_WORLD, ierr) |
---|
| 177 | endif |
---|
| 178 | |
---|
| 179 | if (NORTH_INTER .and. WEST_INTER) then |
---|
| 180 | do k=1,nmax |
---|
| 181 | buf_snd4(k )=A(1,Mmmpi-1,k) |
---|
| 182 | buf_snd4(k+nmax )=A(2,Mmmpi-1,k) |
---|
| 183 | buf_snd4(k+2*nmax)=A(1,Mmmpi ,k) |
---|
| 184 | buf_snd4(k+3*nmax)=A(2,Mmmpi ,k) |
---|
| 185 | enddo |
---|
| 186 | call MPI_Irecv (buf_rev4, ksize, MPI_DOUBLE_PRECISION, p_NW, |
---|
| 187 | & 7, MPI_COMM_WORLD, req(8),ierr) |
---|
| 188 | call MPI_Send (buf_snd4, ksize, MPI_DOUBLE_PRECISION, p_NW, |
---|
| 189 | & 8, MPI_COMM_WORLD, ierr) |
---|
| 190 | endif |
---|
| 191 | ! |
---|
| 192 | ! Wait for completion of receive and fill ghost points: sides... |
---|
| 193 | ! |
---|
| 194 | if (WEST_INTER) then |
---|
| 195 | call MPI_Wait (req(1),status(1,1),ierr) |
---|
| 196 | |
---|
| 197 | do k=1,nmax |
---|
| 198 | do j=jmin,jmax |
---|
| 199 | A(-1,j,k)=jbuf_revW(k+nmax*(j-jmin )) |
---|
| 200 | A( 0,j,k)=jbuf_revW(k+nmax*(j-jmin+jshft)) |
---|
| 201 | enddo |
---|
| 202 | enddo |
---|
| 203 | endif |
---|
| 204 | |
---|
| 205 | if (EAST_INTER) then |
---|
| 206 | call MPI_Wait (req(2),status(1,2),ierr) |
---|
| 207 | |
---|
| 208 | do k=1,nmax |
---|
| 209 | do j=jmin,jmax |
---|
| 210 | A(Lmmpi+1,j,k)=jbuf_revE(k+nmax*(j-jmin )) |
---|
| 211 | A(Lmmpi+2,j,k)=jbuf_revE(k+nmax*(j-jmin+jshft)) |
---|
| 212 | enddo |
---|
| 213 | enddo |
---|
| 214 | endif |
---|
| 215 | |
---|
| 216 | if (SOUTH_INTER) then |
---|
| 217 | call MPI_Wait (req(3),status(1,3),ierr) |
---|
| 218 | |
---|
| 219 | do k=1,nmax |
---|
| 220 | do i=imin,imax |
---|
| 221 | A(i,-1,k)=ibuf_revS(k+nmax*(i-imin )) |
---|
| 222 | A(i, 0,k)=ibuf_revS(k+nmax*(i-imin+ishft)) |
---|
| 223 | enddo |
---|
| 224 | enddo |
---|
| 225 | endif |
---|
| 226 | |
---|
| 227 | if (NORTH_INTER) then |
---|
| 228 | call MPI_Wait (req(4),status(1,4),ierr) |
---|
| 229 | |
---|
| 230 | do k=1,nmax |
---|
| 231 | do i=imin,imax |
---|
| 232 | A(i,Mmmpi+1,k)=ibuf_revN(k+nmax*(i-imin )) |
---|
| 233 | A(i,Mmmpi+2,k)=ibuf_revN(k+nmax*(i-imin+ishft)) |
---|
| 234 | enddo |
---|
| 235 | enddo |
---|
| 236 | endif |
---|
| 237 | ! |
---|
| 238 | ! ...corners: |
---|
| 239 | ! |
---|
| 240 | if (SOUTH_INTER .and. WEST_INTER) then |
---|
| 241 | call MPI_Wait (req(5),status(1,5),ierr) |
---|
| 242 | |
---|
| 243 | do k=1,nmax |
---|
| 244 | A(-1,-1,k)=buf_rev1(k ) |
---|
| 245 | A( 0,-1,k)=buf_rev1(k+nmax ) |
---|
| 246 | A(-1, 0,k)=buf_rev1(k+2*nmax) |
---|
| 247 | A( 0, 0,k)=buf_rev1(k+3*nmax) |
---|
| 248 | enddo |
---|
| 249 | endif |
---|
| 250 | |
---|
| 251 | if (NORTH_INTER.and.EAST_INTER) then |
---|
| 252 | call MPI_Wait (req(6),status(1,6),ierr) |
---|
| 253 | |
---|
| 254 | do k=1,nmax |
---|
| 255 | A(Lmmpi+1,Mmmpi+1,k)=buf_rev2(k ) |
---|
| 256 | A(Lmmpi+2,Mmmpi+1,k)=buf_rev2(k+nmax ) |
---|
| 257 | A(Lmmpi+1,Mmmpi+2,k)=buf_rev2(k+2*nmax) |
---|
| 258 | A(Lmmpi+2,Mmmpi+2,k)=buf_rev2(k+3*nmax) |
---|
| 259 | enddo |
---|
| 260 | endif |
---|
| 261 | |
---|
| 262 | if (SOUTH_INTER .and. EAST_INTER) then |
---|
| 263 | call MPI_Wait (req(7),status(1,7),ierr) |
---|
| 264 | |
---|
| 265 | do k=1,nmax |
---|
| 266 | A(Lmmpi+1,-1,k)=buf_rev3(k ) |
---|
| 267 | A(Lmmpi+2,-1,k)=buf_rev3(k+nmax ) |
---|
| 268 | A(Lmmpi+1, 0,k)=buf_rev3(k+2*nmax) |
---|
| 269 | A(Lmmpi+2, 0,k)=buf_rev3(k+3*nmax) |
---|
| 270 | enddo |
---|
| 271 | endif |
---|
| 272 | |
---|
| 273 | if (NORTH_INTER .and. WEST_INTER) then |
---|
| 274 | call MPI_Wait (req(8),status(1,8),ierr) |
---|
| 275 | |
---|
| 276 | do k=1,nmax |
---|
| 277 | A(-1,Mmmpi+1,k)=buf_rev4(k ) |
---|
| 278 | A( 0,Mmmpi+1,k)=buf_rev4(k+nmax ) |
---|
| 279 | A(-1,Mmmpi+2,k)=buf_rev4(k+2*nmax) |
---|
| 280 | A( 0,Mmmpi+2,k)=buf_rev4(k+3*nmax) |
---|
| 281 | enddo |
---|
| 282 | endif |
---|
| 283 | |
---|
| 284 | #if defined AGRIF |
---|
| 285 | |
---|
| 286 | DeAllocate(buf_snd4, ibuf_sndN, buf_snd2, |
---|
| 287 | & buf_rev4, ibuf_revN, buf_rev2, |
---|
| 288 | & jbuf_sndW, jbuf_sndE, |
---|
| 289 | & jbuf_revW, jbuf_revE, |
---|
| 290 | & buf_snd1, ibuf_sndS, buf_snd3, |
---|
| 291 | & buf_rev1, ibuf_revS, buf_rev3) |
---|
| 292 | #endif |
---|
| 293 | |
---|
| 294 | return |
---|
| 295 | end |
---|
| 296 | #else |
---|
| 297 | subroutine MessPass3D_empty |
---|
| 298 | return |
---|
| 299 | end |
---|
| 300 | #endif |
---|
| 301 | |
---|