Changeset 4007


Ignore:
Timestamp:
2013-08-28T10:10:35+02:00 (8 years ago)
Author:
davestorkey
Message:
  1. Bug fixes for flagu/flagv calculation in bdyini.F90.
  2. Introduce masking of derivatives in radiation velocity calculation in bdylib.F90
  3. Change relaxation term in Orlanski implementation to explicit timestepping in bdylib.F90.
  4. Remove bdyfmask (redundant).
Location:
branches/2013/dev_r3987_METO1_MERCATOR6_OBC_BDY_merge/NEMOGCM/NEMO/OPA_SRC/BDY
Files:
6 edited

Legend:

Unmodified
Added
Removed
  • branches/2013/dev_r3987_METO1_MERCATOR6_OBC_BDY_merge/NEMOGCM/NEMO/OPA_SRC/BDY/bdy_oce.F90

    r3991 r4007  
    109109   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:), TARGET ::   bdyumask   !: Mask defining computational domain at U-points 
    110110   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:), TARGET ::   bdyvmask   !: Mask defining computational domain at V-points 
    111    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:), TARGET ::   bdyfmask   !: Mask defining computational domain at F-points 
    112111 
    113112   REAL(wp)                                    ::   bdysurftot !: Lateral surface of unstructured open boundary 
     
    144143      !!---------------------------------------------------------------------- 
    145144      ! 
    146       ALLOCATE( bdytmask(jpi,jpj) , bdyumask(jpi,jpj), bdyvmask(jpi,jpj), bdyfmask(jpi,jpj),     &   
     145      ALLOCATE( bdytmask(jpi,jpj) , bdyumask(jpi,jpj), bdyvmask(jpi,jpj),     &   
    147146         &      STAT=bdy_oce_alloc ) 
    148147         ! 
  • branches/2013/dev_r3987_METO1_MERCATOR6_OBC_BDY_merge/NEMOGCM/NEMO/OPA_SRC/BDY/bdydyn2d.F90

    r3994 r4007  
    220220      igrd = 2      ! Orlanski bc on u-velocity;  
    221221      !             
    222       CALL bdy_orlanski_2d( idx, igrd, pub2d, pua2d, dta%u2d, umask(:,:,1), ll_npo ) 
     222      CALL bdy_orlanski_2d( idx, igrd, pub2d, pua2d, dta%u2d, ll_npo ) 
    223223 
    224224      igrd = 3      ! Orlanski bc on v-velocity 
    225225      !   
    226       CALL bdy_orlanski_2d( idx, igrd, pvb2d, pva2d, dta%v2d, vmask(:,:,1), ll_npo ) 
     226      CALL bdy_orlanski_2d( idx, igrd, pvb2d, pva2d, dta%v2d, ll_npo ) 
    227227      ! 
    228228      IF( nn_timing == 1 ) CALL timing_stop('bdy_dyn2d_orlanski') 
  • branches/2013/dev_r3987_METO1_MERCATOR6_OBC_BDY_merge/NEMOGCM/NEMO/OPA_SRC/BDY/bdydyn3d.F90

    r3994 r4007  
    239239      igrd = 2      ! Orlanski bc on u-velocity;  
    240240      !             
    241       CALL bdy_orlanski_3d( idx, igrd, ub, ua, dta%u3d, umask, ll_npo ) 
     241      CALL bdy_orlanski_3d( idx, igrd, ub, ua, dta%u3d, ll_npo ) 
    242242 
    243243      igrd = 3      ! Orlanski bc on v-velocity 
    244244      !   
    245       CALL bdy_orlanski_3d( idx, igrd, vb, va, dta%v3d, vmask, ll_npo ) 
     245      CALL bdy_orlanski_3d( idx, igrd, vb, va, dta%v3d, ll_npo ) 
    246246      ! 
    247247      CALL lbc_bdy_lnk( ua, 'U', -1., ib_bdy )    ! Boundary points should be updated 
  • branches/2013/dev_r3987_METO1_MERCATOR6_OBC_BDY_merge/NEMOGCM/NEMO/OPA_SRC/BDY/bdyini.F90

    r3991 r4007  
    2121   !!   bdy_init       : Initialization of unstructured open boundaries 
    2222   !!---------------------------------------------------------------------- 
     23   USE wrk_nemo        ! Memory Allocation 
    2324   USE timing          ! Timing 
    2425   USE oce             ! ocean dynamics and tracers variables 
     
    9293      INTEGER :: com_east_b, com_west_b, com_south_b, com_north_b  ! Flags for boundaries receiving 
    9394      INTEGER :: iw_b(4), ie_b(4), is_b(4), in_b(4)                ! Arrays for neighbours coordinates 
     95      REAL(wp), POINTER, DIMENSION(:,:)       ::   zfmask  ! temporary fmask array excluding coastal boundary condition (shlat) 
    9496 
    9597      !! 
     
    11781180      ENDDO 
    11791181 
    1180       ! bdyfmask required for flagu, flagv calculations below even though F-points  
    1181       ! not defined for BDY grid.  
     1182      ! For the flagu/flagv calculation below we require a version of fmask without 
     1183      ! the land boundary condition (shlat) included: 
     1184      CALL wrk_alloc(jpi,jpj,zfmask)  
    11821185      DO ij = 2, jpjm1 
    11831186         DO ii = 2, jpim1 
    1184             bdyfmask(ii,ij) = bdytmask(ii,ij  ) * bdytmask(ii+1,ij  )   & 
    1185            &                * bdytmask(ii,ij+1) * bdytmask(ii+1,ij+1) 
     1187            zfmask(ii,ij) = tmask(ii,ij  ,1) * tmask(ii+1,ij  ,1)   & 
     1188           &              * tmask(ii,ij+1,1) * tmask(ii+1,ij+1,1) 
    11861189         END DO       
    11871190      END DO 
    11881191 
    11891192      ! Lateral boundary conditions 
     1193      CALL lbc_lnk( zfmask       , 'F', 1. ) 
    11901194      CALL lbc_lnk( fmask        , 'F', 1. )   ;   CALL lbc_lnk( bdytmask(:,:), 'T', 1. ) 
    11911195      CALL lbc_lnk( bdyumask(:,:), 'U', 1. )   ;   CALL lbc_lnk( bdyvmask(:,:), 'V', 1. ) 
     
    12051209            SELECT CASE( igrd ) 
    12061210               CASE( 1 ) 
    1207                   cgrid = 'T'  
    12081211                  pmask => umask(:,:,1) 
    12091212                  i_offset = 0 
    12101213               CASE( 2 )  
    1211                   cgrid = 'U' 
    12121214                  pmask => bdytmask 
    12131215                  i_offset = 1 
    12141216               CASE( 3 )  
    1215                   cgrid = 'V' 
    1216                   pmask => fmask(:,:,1) 
     1217                  pmask => zfmask(:,:) 
    12171218                  i_offset = 0 
    12181219            END SELECT  
     
    12231224               zefl = pmask(nbi+i_offset-1,nbj) 
    12241225               zwfl = pmask(nbi+i_offset,nbj) 
    1225                IF( zefl + zwfl == 2 * i_offset ) THEN 
     1226               ! This error check only works if you are using the bdyXmask arrays 
     1227               IF( i_offset == 1 .and. zefl + zwfl == 2 ) THEN 
    12261228                  icount = icount + 1 
    12271229                  IF(lwp) WRITE(numout,*) 'Problem with igrd = ',igrd,' at (global) nbi, nbj : ',mig(nbi),mjg(nbj) 
    12281230               ELSE 
    1229                   idx_bdy(ib_bdy)%flagu(ib,igrd)=-zefl+zwfl 
     1231                  idx_bdy(ib_bdy)%flagu(ib,igrd) = -zefl + zwfl 
    12301232               ENDIF 
    12311233            END DO 
    12321234            IF( icount /= 0 ) THEN 
    12331235               IF(lwp) WRITE(numout,*) 
    1234                IF(lwp) WRITE(numout,*) ' E R R O R : Some ',cgrid,' grid points,',   & 
     1236               IF(lwp) WRITE(numout,*) ' E R R O R : Some ',cgrid(igrd),' grid points,',   & 
    12351237                  ' are not boundary points (flagu calculation). Check nbi, nbj, indices for boundary set ',ib_bdy 
    12361238               IF(lwp) WRITE(numout,*) ' ========== ' 
     
    12411243 
    12421244         ! Calculate relationship of V direction to the local orientation of the boundary 
    1243          ! flagv = -1 : u component is normal to the dynamical boundary but its direction is outward 
    1244          ! flagv =  0 : u is tangential 
    1245          ! flagv =  1 : u is normal to the boundary and is direction is inward 
     1245         ! flagv = -1 : v component is normal to the dynamical boundary but its direction is outward 
     1246         ! flagv =  0 : v is tangential 
     1247         ! flagv =  1 : v is normal to the boundary and is direction is inward 
    12461248 
    12471249         DO igrd = 1,jpbgrd  
    12481250            SELECT CASE( igrd ) 
    12491251               CASE( 1 ) 
    1250                   cgrid = 'T' 
    12511252                  pmask => vmask(:,:,1) 
    12521253                  j_offset = 0 
    12531254               CASE( 2 ) 
    1254                   cgrid = 'U' 
    1255                   pmask => fmask(:,:,1) 
     1255                  pmask => zfmask(:,:) 
    12561256                  j_offset = 0 
    12571257               CASE( 3 ) 
    1258                   cgrid = 'V' 
    12591258                  pmask => bdytmask 
    12601259                  j_offset = 1 
     
    12661265               znfl = pmask(nbi,nbj+j_offset-1  ) 
    12671266               zsfl = pmask(nbi,nbj+j_offset) 
    1268                IF( znfl + zsfl == 2 * j_offset ) THEN 
     1267               ! This error check only works if you are using the bdyXmask arrays 
     1268               IF( j_offset == 1 .and. znfl + zsfl == 2 ) THEN 
    12691269                  IF(lwp) WRITE(numout,*) 'Problem with igrd = ',igrd,' at (global) nbi, nbj : ',mig(nbi),mjg(nbj) 
    12701270                  icount = icount + 1 
     
    12751275            IF( icount /= 0 ) THEN 
    12761276               IF(lwp) WRITE(numout,*) 
    1277                IF(lwp) WRITE(numout,*) ' E R R O R : Some ',cgrid,' grid points,',   & 
     1277               IF(lwp) WRITE(numout,*) ' E R R O R : Some ',cgrid(igrd),' grid points,',   & 
    12781278                  ' are not boundary points (flagv calculation). Check nbi, nbj, indices for boundary set ',ib_bdy 
    12791279               IF(lwp) WRITE(numout,*) ' ========== ' 
     
    13241324         DEALLOCATE(nbidta, nbjdta, nbrdta) 
    13251325      ENDIF 
     1326 
     1327      CALL wrk_dealloc(jpi,jpj,zfmask)  
    13261328 
    13271329      IF( nn_timing == 1 ) CALL timing_stop('bdy_init') 
  • branches/2013/dev_r3987_METO1_MERCATOR6_OBC_BDY_merge/NEMOGCM/NEMO/OPA_SRC/BDY/bdylib.F90

    r3994 r4007  
    3434CONTAINS 
    3535 
    36    SUBROUTINE bdy_orlanski_2d( idx, igrd, phib, phia, phi_ext, mask, ll_npo ) 
     36   SUBROUTINE bdy_orlanski_2d( idx, igrd, phib, phia, phi_ext, ll_npo ) 
    3737      !!---------------------------------------------------------------------- 
    3838      !!                 ***  SUBROUTINE bdy_orlanski_2d  *** 
     
    5050      REAL(wp), DIMENSION(:,:),   INTENT(inout)  ::   phia     ! model after 2D field (to be updated) 
    5151      REAL(wp), DIMENSION(:),     INTENT(in)     ::   phi_ext  ! external forcing data 
    52       REAL(wp), DIMENSION(:,:),   INTENT(in)     ::   mask     ! land/sea mask 
    5352      LOGICAL,                    INTENT(in)     ::   ll_npo   ! switch for NPO version 
    5453 
     
    5756      INTEGER  ::   iijm1, iijp1, ijjm1, ijjp1             ! 2D addresses 
    5857      INTEGER  ::   iibm1jp1, iibm1jm1, ijbm1jp1, ijbm1jm1 ! 2D addresses 
     58      INTEGER  ::   ii_offset, ij_offset                   ! offsets for mask indices 
    5959      INTEGER  ::   flagu, flagv                           ! short cuts 
    6060      REAL(wp) ::   zdt, zdx, zdy, znor2, zcx, zcy         ! intermediate calculations 
    61       REAL(wp) ::   zout, zwgt, zdy_centred, zsign_ups 
     61      REAL(wp) ::   zout, zwgt, zdy_centred 
     62      REAL(wp) ::   zdy_left, zdy_right, zsign_ups 
     63      REAL(wp), POINTER, DIMENSION(:,:)          :: pmask      ! land/sea mask for field 
     64      REAL(wp), POINTER, DIMENSION(:,:)          :: pmask_xdiv ! land/sea mask for x-derivatives 
     65      REAL(wp), POINTER, DIMENSION(:,:)          :: pmask_ydiv ! land/sea mask for y-derivatives 
    6266      !!---------------------------------------------------------------------- 
    6367 
     
    6872      ! ----------------------------------!  
    6973      
     74      SELECT CASE(igrd) 
     75         CASE(1) 
     76            pmask => tmask(:,:,1) 
     77            pmask_xdiv => umask(:,:,1) 
     78            ii_offset = 0 
     79            pmask_ydiv => vmask(:,:,1) 
     80            ij_offset = 0 
     81         CASE(2) 
     82            pmask => umask(:,:,1) 
     83            pmask_xdiv => tmask(:,:,1) 
     84            ii_offset = 1 
     85            pmask_ydiv => fmask(:,:,1) 
     86            ij_offset = 0 
     87         CASE(3) 
     88            pmask => vmask(:,:,1) 
     89            pmask_xdiv => fmask(:,:,1) 
     90            ii_offset = 0 
     91            pmask_ydiv => tmask(:,:,1) 
     92            ij_offset = 1 
     93         CASE DEFAULT ;   CALL ctl_stop( 'unrecognised value for igrd in bdy_orlanksi_2d' ) 
     94      END SELECT 
    7095      ! 
    7196      DO jb = 1, idx%nblenrim(igrd) 
     
    86111         ijbm1jm1 = ij + flagv - abs(flagu) ; ijbm1jp1 = ij + flagv + abs(flagu)  
    87112         ! 
    88          ! calculate normal (zcx) and tangential (zcy) components of radiation velocities: 
     113         ! Calculate normal (zcx) and tangential (zcy) components of radiation velocities. 
     114         ! Mask derivatives to ensure correct land boundary conditions for each variable. 
     115         ! Centred derivative is calculated as average of "left" and "right" derivatives for  
     116         ! this reason.  
    89117         zdt = phia(iibm1,ijbm1) - phib(iibm1,ijbm1) 
    90          zdx = phia(iibm1,ijbm1) - phia(iibm2,ijbm2) 
    91          zdy_centred = phib(iibm1jp1,ijbm1jp1) - phib(iibm1jm1,ijbm1jm1) 
     118         zdx = ( phia(iibm1,ijbm1) - phia(iibm2,ijbm2) )                          & 
     119        &    * ( abs(iibm1-iibm2) * pmask_xdiv(iibm2+ii_offset,ijbm2          )   & 
     120        &      + abs(ijbm1-ijbm2) * pmask_ydiv(iibm2          ,ijbm2+ij_offset) )  
     121         zdy_left  = phib(iibm1   ,ijbm1   ) - phib(iibm1jm1,ijbm1jm1)                  & 
     122        &    * ( (iibm1-iibm1jm1) * pmask_xdiv(iibm1jm1+ii_offset,ijbm1jm1          )   &  
     123        &      + (ijbm1-ijbm1jm1) * pmask_ydiv(iibm1jm1          ,ijbm1jm1+ij_offset) )  
     124         zdy_right = phib(iibm1jp1,ijbm1jp1) - phib(iibm1   ,ijbm1)                     & 
     125        &    * ( (iibm1jp1-iibm1) * pmask_xdiv(iibm1+ii_offset,ijbm1)                   & 
     126        &      + (ijbm1jp1-ijbm1) * pmask_ydiv(iibm1          ,ijbm1+ij_offset) )  
     127         zdy_centred = 0.5 * ( zdy_left + zdy_right ) 
     128!!$         zdy_centred = phib(iibm1jp1,ijbm1jp1) - phib(iibm1jm1,ijbm1jm1) 
    92129         ! upstream differencing for tangential derivatives 
    93130         zsign_ups = sign( 1., zdt * zdy_centred ) 
    94131         zsign_ups = 0.5*( zsign_ups + abs(zsign_ups) ) 
    95          zdy = zsign_ups        * ( phib(iibm1   ,ijbm1   ) - phib(iibm1jm1,ijbm1jm1) ) & 
    96        &        + (1. - zsign_ups) * ( phib(iibm1jp1,ijbm1jp1) - phib(iibm1   ,ijbm1   ) ) 
     132         zdy = zsign_ups * zdy_left + (1. - zsign_ups) * zdy_right 
    97133         znor2 = zdx * zdx + zdy * zdy 
    98134         znor2 = max(znor2,rsmall) 
     
    106142         ! only apply radiation on outflow points  
    107143         if( ll_npo ) then     !! NPO version !! 
    108             phia(ii,ij) = (1.-zout) * phia(ii,ij)                                                   & 
     144            phia(ii,ij) = (1.-zout) * phib(ii,ij)                                                   & 
    109145           &            + zout      * ( phib(ii,ij) + zcx*phia(iibm1,ijbm1) ) / ( 1. + zcx )  
    110146         else                  !! full oblique radiation !! 
    111147            zsign_ups = sign( 1., zcy ) 
    112148            zsign_ups = 0.5*( zsign_ups + abs(zsign_ups) ) 
    113             phia(ii,ij) = (1.-zout) * phia(ii,ij)                                                   & 
     149            phia(ii,ij) = (1.-zout) * phib(ii,ij)                                                   & 
    114150           &            + zout      * ( phib(ii,ij) + zcx*phia(iibm1,ijbm1)                         & 
    115151           &                    - zsign_ups      * zcy * ( phib(ii   ,ij   ) - phib(iijm1,ijjm1 ) ) & 
    116152           &                    - (1.-zsign_ups) * zcy * ( phib(iijp1,ijjp1) - phib(ii   ,ij    ) ) ) / ( 1. + zcx )  
    117153         end if 
    118          phia(ii,ij) = phia(ii,ij) + zwgt * ( phi_ext(jb) - phia(ii,ij) )  
    119          phia(ii,ij) = phia(ii,ij) * mask(ii,ij) 
     154!!$         phia(ii,ij) = phia(ii,ij) + zwgt * ( phi_ext(jb) - phib(ii,ij) )  
     155         phia(ii,ij) = phia(ii,ij) * pmask(ii,ij) 
    120156      END DO 
    121157      ! 
     
    125161 
    126162 
    127    SUBROUTINE bdy_orlanski_3d( idx, igrd, phib, phia, phi_ext, mask, ll_npo ) 
     163   SUBROUTINE bdy_orlanski_3d( idx, igrd, phib, phia, phi_ext, ll_npo ) 
    128164      !!---------------------------------------------------------------------- 
    129165      !!                 ***  SUBROUTINE bdy_orlanski_3d  *** 
     
    141177      REAL(wp), DIMENSION(:,:,:), INTENT(inout)  ::   phia     ! model after 3D field (to be updated) 
    142178      REAL(wp), DIMENSION(:,:),   INTENT(in)     ::   phi_ext  ! external forcing data 
    143       REAL(wp), DIMENSION(:,:,:), INTENT(in)     ::   mask     ! land/sea mask 
    144179      LOGICAL,                    INTENT(in)     ::   ll_npo   ! switch for NPO version 
    145180 
     
    148183      INTEGER  ::   iijm1, iijp1, ijjm1, ijjp1             ! 2D addresses 
    149184      INTEGER  ::   iibm1jp1, iibm1jm1, ijbm1jp1, ijbm1jm1 ! 2D addresses 
     185      INTEGER  ::   ii_offset, ij_offset                   ! offsets for mask indices 
    150186      INTEGER  ::   flagu, flagv                           ! short cuts 
    151187      REAL(wp) ::   zdt, zdx, zdy, znor2, zcx, zcy         ! intermediate calculations 
    152       REAL(wp) ::   zout, zwgt, zdy_centred, zsign_ups 
     188      REAL(wp) ::   zout, zwgt, zdy_centred 
     189      REAL(wp) ::   zdy_left, zdy_right,  zsign_ups 
     190      REAL(wp), POINTER, DIMENSION(:,:,:)        :: pmask      ! land/sea mask for field 
     191      REAL(wp), POINTER, DIMENSION(:,:,:)        :: pmask_xdiv ! land/sea mask for x-derivatives 
     192      REAL(wp), POINTER, DIMENSION(:,:,:)        :: pmask_ydiv ! land/sea mask for y-derivatives 
    153193      !!---------------------------------------------------------------------- 
    154194 
     
    159199      ! ----------------------------------!  
    160200      
     201      SELECT CASE(igrd) 
     202         CASE(1) 
     203            pmask => tmask(:,:,:) 
     204            pmask_xdiv => umask(:,:,:) 
     205            ii_offset = 0 
     206            pmask_ydiv => vmask(:,:,:) 
     207            ij_offset = 0 
     208         CASE(2) 
     209            pmask => umask(:,:,:) 
     210            pmask_xdiv => tmask(:,:,:) 
     211            ii_offset = 1 
     212            pmask_ydiv => fmask(:,:,:) 
     213            ij_offset = 0 
     214         CASE(3) 
     215            pmask => vmask(:,:,:) 
     216            pmask_xdiv => fmask(:,:,:) 
     217            ii_offset = 0 
     218            pmask_ydiv => tmask(:,:,:) 
     219            ij_offset = 1 
     220         CASE DEFAULT ;   CALL ctl_stop( 'unrecognised value for igrd in bdy_orlanksi_2d' ) 
     221      END SELECT 
     222 
    161223      DO jk = 1, jpk 
    162224         !             
     
    178240            ijbm1jm1 = ij + flagv - abs(flagu) ; ijbm1jp1 = ij + flagv + abs(flagu)  
    179241            ! 
    180             ! calculate normal (zcx) and tangential (zcy) components of radiation velocities: 
     242            ! Calculate normal (zcx) and tangential (zcy) components of radiation velocities. 
     243            ! Mask derivatives to ensure correct land boundary conditions for each variable. 
     244            ! Centred derivative is calculated as average of "left" and "right" derivatives for  
     245            ! this reason.  
    181246            zdt = phia(iibm1,ijbm1,jk) - phib(iibm1,ijbm1,jk) 
    182             zdx = phia(iibm1,ijbm1,jk) - phia(iibm2,ijbm2,jk) 
    183             zdy_centred = phib(iibm1jp1,ijbm1jp1,jk) - phib(iibm1jm1,ijbm1jm1,jk) 
     247            zdx = phia(iibm1,ijbm1,jk) - phia(iibm2,ijbm2,jk)                     & 
     248        &    * ( (iibm1-iibm2) * pmask_xdiv(iibm2+ii_offset,ijbm2          ,jk)   & 
     249        &      + (ijbm1-ijbm2) * pmask_ydiv(iibm2          ,ijbm2+ij_offset,jk) )  
     250            zdy_left  = phib(iibm1   ,ijbm1   ,jk) - phib(iibm1jm1,ijbm1jm1,jk)            & 
     251        &    * ( (iibm1-iibm1jm1) * pmask_xdiv(iibm1jm1+ii_offset,ijbm1jm1          ,jk)   &  
     252        &      + (ijbm1-ijbm1jm1) * pmask_ydiv(iibm1jm1          ,ijbm1jm1+ij_offset,jk) )  
     253            zdy_right = phib(iibm1jp1,ijbm1jp1,jk) - phib(iibm1   ,ijbm1   ,jk)      & 
     254        &    * ( (iibm1jp1-iibm1) * pmask_xdiv(iibm1+ii_offset,ijbm1          ,jk)   & 
     255        &      + (ijbm1jp1-ijbm1) * pmask_ydiv(iibm1          ,ijbm1+ij_offset,jk) )  
     256            zdy_centred = 0.5 * ( zdy_left + zdy_right ) 
     257!!$            zdy_centred = phib(iibm1jp1,ijbm1jp1,jk) - phib(iibm1jm1,ijbm1jm1,jk) 
    184258            ! upstream differencing for tangential derivatives 
    185259            zsign_ups = sign( 1., zdt * zdy_centred ) 
    186260            zsign_ups = 0.5*( zsign_ups + abs(zsign_ups) ) 
    187             zdy = zsign_ups        * ( phib(iibm1   ,ijbm1   ,jk) - phib(iibm1jm1,ijbm1jm1,jk) ) & 
    188            &    + (1. - zsign_ups) * ( phib(iibm1jp1,ijbm1jp1,jk) - phib(iibm1   ,ijbm1   ,jk) ) 
     261            zdy = zsign_ups * zdy_left + (1. - zsign_ups) * zdy_right 
    189262            znor2 = zdx * zdx + zdy * zdy 
    190263            znor2 = max(znor2,rsmall) 
     
    198271            ! only apply radiation on outflow points  
    199272            if( ll_npo ) then     !! NPO version !! 
    200                phia(ii,ij,jk) = (1.-zout) * phia(ii,ij,jk)                                                   & 
     273               phia(ii,ij,jk) = (1.-zout) * phib(ii,ij,jk)                                                   & 
    201274              &               + zout      * ( phib(ii,ij,jk) + zcx*phia(iibm1,ijbm1,jk) ) / ( 1. + zcx )  
    202275            else                  !! full oblique radiation !! 
    203276               zsign_ups = sign( 1., zcy ) 
    204277               zsign_ups = 0.5*( zsign_ups + abs(zsign_ups) ) 
    205                phia(ii,ij,jk) = (1.-zout) * phia(ii,ij,jk)                                                   & 
     278               phia(ii,ij,jk) = (1.-zout) * phib(ii,ij,jk)                                                   & 
    206279              &               + zout      * ( phib(ii,ij,jk) + zcx*phia(iibm1,ijbm1,jk)                      & 
    207280              &                    - zsign_ups      * zcy * ( phib(ii   ,ij   ,jk) - phib(iijm1,ijjm1,jk ) ) & 
    208281              &                    - (1.-zsign_ups) * zcy * ( phib(iijp1,ijjp1,jk) - phib(ii   ,ij   ,jk ) ) ) / ( 1. + zcx )  
    209282            end if 
    210             phia(ii,ij,jk) = phia(ii,ij,jk) + zwgt * ( phi_ext(jb,jk) - phia(ii,ij,jk) )  
    211             phia(ii,ij,jk) = phia(ii,ij,jk) * mask(ii,ij,jk) 
     283!!$            phia(ii,ij,jk) = phia(ii,ij,jk) + zwgt * ( phi_ext(jb,jk) - phib(ii,ij,jk) )  
     284            phia(ii,ij,jk) = phia(ii,ij,jk) * pmask(ii,ij,jk) 
    212285         END DO 
    213286         ! 
  • branches/2013/dev_r3987_METO1_MERCATOR6_OBC_BDY_merge/NEMOGCM/NEMO/OPA_SRC/BDY/bdytra.F90

    r3994 r4007  
    224224      igrd = 1      ! Orlanski bc on temperature;  
    225225      !             
    226       CALL bdy_orlanski_3d( idx, igrd, tsb(:,:,:,jp_tem), tsa(:,:,:,jp_tem), dta%tem, tmask, ll_npo ) 
     226      CALL bdy_orlanski_3d( idx, igrd, tsb(:,:,:,jp_tem), tsa(:,:,:,jp_tem), dta%tem, ll_npo ) 
    227227 
    228228      igrd = 1      ! Orlanski bc on salinity; 
    229229      !   
    230       CALL bdy_orlanski_3d( idx, igrd, tsb(:,:,:,jp_sal), tsa(:,:,:,jp_sal), dta%sal, tmask, ll_npo ) 
     230      CALL bdy_orlanski_3d( idx, igrd, tsb(:,:,:,jp_sal), tsa(:,:,:,jp_sal), dta%sal, ll_npo ) 
    231231      ! 
    232232      IF( nn_timing == 1 ) CALL timing_stop('bdy_tra_orlanski') 
Note: See TracChangeset for help on using the changeset viewer.