Ignore:
Timestamp:
Mar 22, 2007 9:54:05 AM (15 years ago)
Author:
raasch
Message:

preliminary update for changes concerning non-cyclic boundary conditions

File:
1 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/SOURCE/exchange_horiz.f90

    r4 r75  
    1  SUBROUTINE exchange_horiz( ar, xrp, ynp )
     1 SUBROUTINE exchange_horiz( ar )
    22
    33!------------------------------------------------------------------------------!
    44! Actual revisions:
    55! -----------------
    6 !
     6! Special cases for additional gridpoints along x or y in case of non-cyclic
     7! boundary conditions are not regarded any more
    78!
    89! Former revisions:
     
    3233    IMPLICIT NONE
    3334
    34     INTEGER ::  xrp, ynp
    35 
    3635#if defined( __parallel )
    37     INTEGER                               ::  typexz
    3836    INTEGER, DIMENSION(4)                 ::  req
    3937    INTEGER, DIMENSION(MPI_STATUS_SIZE,4) ::  wait_stat
    4038#endif
    4139
    42     REAL ::  ar(nzb:nzt+1,nys-1:nyn+ynp+1,nxl-1:nxr+xrp+1)
     40    REAL ::  ar(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1)
    4341
    4442
     
    5957
    6058    ELSE
     59
    6160       req = 0
    6261!
     
    7675               ar(nzb,nys-1,nxl-1), ngp_yz(grid_level), MPI_REAL, pleft,  1, &
    7776                          comm2d, req(4), ierr )
    78        call MPI_Waitall (4,req,wait_stat,ierr)
     77       CALL MPI_WAITALL( 4, req, wait_stat, ierr )
     78
    7979    ENDIF
    8080
     
    9090
    9191    ELSE
    92 !
    93 !--    Set the MPI data type, which depends on the size of the array
    94 !--    (the v array has an additional gridpoint along y in case of non-cyclic
    95 !--    boundary conditions)
    96        IF ( ynp == 0 )  THEN
    97           typexz = type_xz(grid_level)
    98        ELSE
    99           typexz = type_xz_p
    100        ENDIF
    10192
    10293       req = 0
    10394!
    10495!--    Send front boundary, receive rear one
    105        CALL MPI_ISEND( ar(nzb,nys,nxl-1),   1, typexz, psouth, 0, comm2d, &
    106                        req(1), ierr )
    107        CALL MPI_IRECV( ar(nzb,nyn+1,nxl-1), 1, typexz, pnorth, 0, comm2d, &
    108                        req(2), ierr )
     96       CALL MPI_ISEND( ar(nzb,nys,nxl-1),   1, type_xz(grid_level), psouth, 0, &
     97                       comm2d, req(1), ierr )
     98       CALL MPI_IRECV( ar(nzb,nyn+1,nxl-1), 1, type_xz(grid_level), pnorth, 0, &
     99                       comm2d, req(2), ierr )
    109100!
    110101!--    Send rear boundary, receive front one
    111        CALL MPI_ISEND( ar(nzb,nyn,nxl-1),   1, typexz, pnorth, 1, comm2d, &
    112                        req(3), ierr )
    113        CALL MPI_IRECV( ar(nzb,nys-1,nxl-1), 1, typexz, psouth, 1, comm2d, &
    114                        req(4), ierr )
    115        call MPI_Waitall (4,req,wait_stat,ierr)
     102       CALL MPI_ISEND( ar(nzb,nyn,nxl-1),   1, type_xz(grid_level), pnorth, 1, &
     103                       comm2d, req(3), ierr )
     104       CALL MPI_IRECV( ar(nzb,nys-1,nxl-1), 1, type_xz(grid_level), psouth, 1, &
     105                       comm2d, req(4), ierr )
     106       call MPI_WAITALL( 4, req, wait_stat, ierr )
     107
    116108    ENDIF
    117109
Note: See TracChangeset for help on using the changeset viewer.