Changeset 11822 for NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/BDY/bdytra.F90

Timestamp:

2019-10-29T11:41:36+01:00 (4 years ago)

Author:

acc

Message:

Branch 2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps. Sette tested updates to branch to align with trunk changes between 10721 and 11740. Sette tests are passing but results differ from branch before these changes (except for GYRE_PISCES and VORTEX) and branch results already differed from trunk because of algorithmic fixes. Will need more checks to confirm correctness.

File:

• NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/BDY/bdytra.F90 (modified) (2 diffs)

NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/BDY/bdytra.F90

@@ -51 +51 @@
       REAL(wp), DIMENSION(jpi,jpj,jpk,jpts,jpt), INTENT(inout) :: pts       ! tracer fields
       !
-      INTEGER                        :: ib_bdy, jn, igrd   ! Loop indices
-      TYPE(ztrabdy), DIMENSION(jpts) :: zdta               ! Temporary data structure
+      INTEGER                        :: ib_bdy, jn, igrd, ir   ! Loop indeces
+      TYPE(ztrabdy), DIMENSION(jpts) :: zdta                   ! Temporary data structure
+      LOGICAL                        :: llrim0                 ! indicate if rim 0 is treated
+      LOGICAL, DIMENSION(4)          :: llsend1, llrecv1       ! indicate how communications are to be carried out
       !!----------------------------------------------------------------------
       igrd = 1 
-
-      DO ib_bdy=1, nb_bdy
+      llsend1(:) = .false.  ;   llrecv1(:) = .false.
+      DO ir = 1, 0, -1   ! treat rim 1 before rim 0
+         IF( ir == 0 ) THEN   ;   llrim0 = .TRUE.
+         ELSE                 ;   llrim0 = .FALSE.
+         END IF
+         DO ib_bdy=1, nb_bdy
+            !
+            zdta(1)%tra => dta_bdy(ib_bdy)%tem
+            zdta(2)%tra => dta_bdy(ib_bdy)%sal
+            !
+            DO jn = 1, jpts
+               !
+               SELECT CASE( TRIM(cn_tra(ib_bdy)) )
+               CASE('none'        )   ;   CYCLE
+               CASE('frs'         )   ! treat the whole boundary at once
+                  IF( ir == 0 ) CALL bdy_frs ( idx_bdy(ib_bdy),                pts(:,:,:,jn,Kaa), zdta(jn)%tra )
+               CASE('specified'   )   ! treat the whole rim      at once
+                  IF( ir == 0 ) CALL bdy_spe ( idx_bdy(ib_bdy),                pts(:,:,:,jn,Kaa), zdta(jn)%tra )
+               CASE('neumann'     )   ;   CALL bdy_nmn ( idx_bdy(ib_bdy), igrd         , pts(:,:,:,jn,Kaa), llrim0 )   ! tsa masked
+               CASE('orlanski'    )   ;   CALL bdy_orl ( idx_bdy(ib_bdy), pts(:,:,:,jn,Kbb), pts(:,:,:,jn,Kaa), &
+                    & zdta(jn)%tra, llrim0, ll_npo=.false. )
+               CASE('orlanski_npo')   ;   CALL bdy_orl ( idx_bdy(ib_bdy), pts(:,:,:,jn,Kbb), pts(:,:,:,jn,Kaa), &
+                    & zdta(jn)%tra, llrim0, ll_npo=.true.  )
+               CASE('runoff'      )   ;   CALL bdy_rnf ( idx_bdy(ib_bdy),                pts(:,:,:,jn,Kaa), jn, llrim0 )
+               CASE DEFAULT           ;   CALL ctl_stop( 'bdy_tra : unrecognised option for open boundaries for T and S' )
+               END SELECT
+               ! 
+            END DO
+         END DO
          !
-         zdta(1)%tra => dta_bdy(ib_bdy)%tem
-         zdta(2)%tra => dta_bdy(ib_bdy)%sal
+         IF( nn_hls > 1 .AND. ir == 1 ) CYCLE   ! at least 2 halos will be corrected -> no need to correct rim 1 before rim 0
+         IF( nn_hls == 1 ) THEN   ;   llsend1(:) = .false.   ;   llrecv1(:) = .false.   ;   END IF
+         DO ib_bdy=1, nb_bdy
+            SELECT CASE( TRIM(cn_tra(ib_bdy)) )
+            CASE('neumann','runoff')
+               llsend1(:) = llsend1(:) .OR. lsend_bdyint(ib_bdy,1,:,ir)   ! possibly every direction, T points
+               llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(ib_bdy,1,:,ir)   ! possibly every direction, T points
+            CASE('orlanski', 'orlanski_npo')
+               llsend1(:) = llsend1(:) .OR. lsend_bdy(ib_bdy,1,:,ir)   ! possibly every direction, T points
+               llrecv1(:) = llrecv1(:) .OR. lrecv_bdy(ib_bdy,1,:,ir)   ! possibly every direction, T points
+            END SELECT
+         END DO
+         IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN   ! if need to send/recv in at least one direction
+            CALL lbc_lnk( 'bdytra', pts(:,:,:,jn,Kaa), 'T',  1., kfillmode=jpfillnothing ,lsend=llsend1, lrecv=llrecv1 )
+         END IF
          !
-         DO jn = 1, jpts
-            !
-            SELECT CASE( TRIM(cn_tra(ib_bdy)) )
-            CASE('none'        )   ;   CYCLE
-            CASE('frs'         )   ;   CALL bdy_frs ( idx_bdy(ib_bdy),                pts(:,:,:,jn,Kaa), zdta(jn)%tra )
-            CASE('specified'   )   ;   CALL bdy_spe ( idx_bdy(ib_bdy),                pts(:,:,:,jn,Kaa), zdta(jn)%tra )
-            CASE('neumann'     )   ;   CALL bdy_nmn ( idx_bdy(ib_bdy), igrd         , pts(:,:,:,jn,Kaa)               )
-            CASE('orlanski'    )   ;   CALL bdy_orl ( idx_bdy(ib_bdy), pts(:,:,:,jn,Kbb), pts(:,:,:,jn,Kaa), zdta(jn)%tra, ll_npo=.false. )
-            CASE('orlanski_npo')   ;   CALL bdy_orl ( idx_bdy(ib_bdy), pts(:,:,:,jn,Kbb), pts(:,:,:,jn,Kaa), zdta(jn)%tra, ll_npo=.true. )
-            CASE('runoff'      )   ;   CALL bdy_rnf ( idx_bdy(ib_bdy),                pts(:,:,:,jn,Kaa),               jn )
-            CASE DEFAULT           ;   CALL ctl_stop( 'bdy_tra : unrecognised option for open boundaries for T and S' )
-            END SELECT
-            ! Boundary points should be updated
-            CALL lbc_bdy_lnk( 'bdytra', pts(:,:,:,jn,Kaa), 'T', 1., ib_bdy )
-            ! 
-         END DO
-      END DO
+      END DO   ! ir
       !
    END SUBROUTINE bdy_tra
 
 
-   SUBROUTINE bdy_rnf( idx, pt, jpa )
+   SUBROUTINE bdy_rnf( idx, pt, jpa, llrim0 )
       !!----------------------------------------------------------------------
       !!                 ***  SUBROUTINE bdy_rnf  ***
@@ -91 +117 @@
       !! 
       !!----------------------------------------------------------------------
-      TYPE(OBC_INDEX),                     INTENT(in) ::   idx  ! OBC indices
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pt  ! tracer trend
-      INTEGER,                             INTENT(in) ::   jpa  ! TRA index
+      TYPE(OBC_INDEX),                     INTENT(in) ::   idx      ! OBC indices
+      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pt       ! tracer trend
+      INTEGER,                             INTENT(in) ::   jpa      ! TRA index
+      LOGICAL,                             INTENT(in) ::   llrim0   ! indicate if rim 0 is treated
       !
-      REAL(wp) ::   zwgt           ! boundary weight
-      INTEGER  ::   ib, ik, igrd   ! dummy loop indices
-      INTEGER  ::   ii, ij, ip, jp ! 2D addresses
+      INTEGER  ::   ib, ii, ij, igrd   ! dummy loop indices
       !!----------------------------------------------------------------------
       !
       igrd = 1                       ! Everything is at T-points here
-      DO ib = 1, idx%nblenrim(igrd)
-         ii = idx%nbi(ib,igrd)
-         ij = idx%nbj(ib,igrd)
-         DO ik = 1, jpkm1
-            ip = bdytmask(ii+1,ij  ) - bdytmask(ii-1,ij  )
-            jp = bdytmask(ii  ,ij+1) - bdytmask(ii  ,ij-1)
-            if (jpa == jp_tem) pt(ii,ij,ik) = pt(ii+ip,ij+jp,ik) * tmask(ii,ij,ik)
-            if (jpa == jp_sal) pt(ii,ij,ik) =                 0.1 * tmask(ii,ij,ik)
+      IF(      jpa == jp_tem ) THEN
+         CALL bdy_nmn( idx, igrd, pt, llrim0 )
+      ELSE IF( jpa == jp_sal ) THEN
+         IF( .NOT. llrim0 )   RETURN
+         DO ib = 1, idx%nblenrim(igrd)   ! if llrim0 then treat the whole rim
+            ii = idx%nbi(ib,igrd)
+            ij = idx%nbj(ib,igrd)
+            pt(ii,ij,1:jpkm1) = 0.1 * tmask(ii,ij,1:jpkm1)
          END DO
-      END DO
+      END IF
       !
    END SUBROUTINE bdy_rnf