Changeset 10253 for branches/UKMO/dev_r5518_AMM15_package/NEMOGCM/NEMO/OPA_SRC/BDY/bdylib.F90

Timestamp:

2018-10-29T15:55:40+01:00 (5 years ago)

Author:

kingr

Message:

Merged AMM15_v3_6_STABLE_package_collate@10237

File:

• branches/UKMO/dev_r5518_AMM15_package/NEMOGCM/NEMO/OPA_SRC/BDY/bdylib.F90 (modified) (3 diffs)

branches/UKMO/dev_r5518_AMM15_package/NEMOGCM/NEMO/OPA_SRC/BDY/bdylib.F90

@@ -26 +26 @@
    PUBLIC   bdy_orlanski_2d     ! routine called where?
    PUBLIC   bdy_orlanski_3d     ! routine called where?
+   PUBLIC   bdy_nmn     ! routine called where? 
 
    !!----------------------------------------------------------------------
@@ -354 +355 @@
    END SUBROUTINE bdy_orlanski_3d
 
+   SUBROUTINE bdy_nmn( idx, igrd, phia ) 
+      !!---------------------------------------------------------------------- 
+      !!                 ***  SUBROUTINE bdy_nmn  *** 
+      !!                     
+      !! ** Purpose : Duplicate the value at open boundaries, zero gradient. 
+      !!  
+      !!---------------------------------------------------------------------- 
+      INTEGER,                    INTENT(in)     ::   igrd     ! grid index 
+      REAL(wp), DIMENSION(:,:,:), INTENT(inout)  ::   phia     ! model after 3D field (to be updated) 
+      TYPE(OBC_INDEX), INTENT(in) ::   idx  ! OBC indices 
+      !!  
+      REAL(wp) ::   zcoef, zcoef1, zcoef2 
+      REAL(wp), POINTER, DIMENSION(:,:,:)        :: pmask      ! land/sea mask for field 
+      REAL(wp), POINTER, DIMENSION(:,:)        :: bdypmask      ! land/sea mask for field 
+      INTEGER  ::   ib, ik   ! dummy loop indices 
+      INTEGER  ::   ii, ij, ip, jp   ! 2D addresses 
+      !!---------------------------------------------------------------------- 
+      ! 
+      IF( nn_timing == 1 ) CALL timing_start('bdy_nmn') 
+      ! 
+      SELECT CASE(igrd) 
+         CASE(1) 
+            pmask => tmask(:,:,:) 
+            bdypmask => bdytmask(:,:) 
+         CASE(2) 
+            pmask => umask(:,:,:) 
+            bdypmask => bdyumask(:,:) 
+         CASE(3) 
+            pmask => vmask(:,:,:) 
+            bdypmask => bdyvmask(:,:) 
+         CASE DEFAULT ;   CALL ctl_stop( 'unrecognised value for igrd in bdy_nmn' ) 
+      END SELECT 
+      DO ib = 1, idx%nblenrim(igrd) 
+         ii = idx%nbi(ib,igrd) 
+         ij = idx%nbj(ib,igrd) 
+         DO ik = 1, jpkm1 
+            ! search the sense of the gradient 
+            zcoef1 = bdypmask(ii-1,ij  )*pmask(ii-1,ij,ik) +  bdypmask(ii+1,ij  )*pmask(ii+1,ij,ik) 
+            zcoef2 = bdypmask(ii  ,ij-1)*pmask(ii,ij-1,ik) +  bdypmask(ii  ,ij+1)*pmask(ii,ij+1,ik) 
+            IF ( nint(zcoef1+zcoef2) == 0) THEN 
+               ! corner **** we probably only want to set the tangentail component for the dynamics here 
+               zcoef = pmask(ii-1,ij,ik) + pmask(ii+1,ij,ik) +  pmask(ii,ij-1,ik) +  pmask(ii,ij+1,ik) 
+               IF (zcoef > .5_wp) THEN ! Only set none isolated points. 
+                 phia(ii,ij,ik) = phia(ii-1,ij  ,ik) * pmask(ii-1,ij  ,ik) + & 
+                   &              phia(ii+1,ij  ,ik) * pmask(ii+1,ij  ,ik) + & 
+                   &              phia(ii  ,ij-1,ik) * pmask(ii  ,ij-1,ik) + & 
+                   &              phia(ii  ,ij+1,ik) * pmask(ii  ,ij+1,ik) 
+                 phia(ii,ij,ik) = ( phia(ii,ij,ik) / zcoef ) * pmask(ii,ij,ik) 
+               ELSE 
+                 phia(ii,ij,ik) = phia(ii,ij  ,ik) * pmask(ii,ij  ,ik) 
+               ENDIF 
+            ELSEIF ( nint(zcoef1+zcoef2) == 2) THEN 
+               ! oblique corner **** we probably only want to set the normal component for the dynamics here 
+               zcoef = pmask(ii-1,ij,ik)*bdypmask(ii-1,ij  ) + pmask(ii+1,ij,ik)*bdypmask(ii+1,ij  ) + & 
+                   &   pmask(ii,ij-1,ik)*bdypmask(ii,ij -1 ) +  pmask(ii,ij+1,ik)*bdypmask(ii,ij+1  ) 
+               phia(ii,ij,ik) = phia(ii-1,ij  ,ik) * pmask(ii-1,ij  ,ik)*bdypmask(ii-1,ij  ) + & 
+                   &            phia(ii+1,ij  ,ik) * pmask(ii+1,ij  ,ik)*bdypmask(ii+1,ij  )  + & 
+                   &            phia(ii  ,ij-1,ik) * pmask(ii  ,ij-1,ik)*bdypmask(ii,ij -1 ) + & 
+                   &            phia(ii  ,ij+1,ik) * pmask(ii  ,ij+1,ik)*bdypmask(ii,ij+1  ) 
+   
+               phia(ii,ij,ik) = ( phia(ii,ij,ik) / MAX(1._wp, zcoef) ) * pmask(ii,ij,ik) 
+            ELSE 
+               ip = nint(bdypmask(ii+1,ij  )*pmask(ii+1,ij,ik) - bdypmask(ii-1,ij  )*pmask(ii-1,ij,ik)) 
+               jp = nint(bdypmask(ii  ,ij+1)*pmask(ii,ij+1,ik) - bdypmask(ii  ,ij-1)*pmask(ii,ij-1,ik)) 
+               phia(ii,ij,ik) = phia(ii+ip,ij+jp,ik) * pmask(ii+ip,ij+jp,ik) 
+            ENDIF 
+         END DO 
+      END DO 
+      ! 
+      IF( nn_timing == 1 ) CALL timing_stop('bdy_nmn') 
+      ! 
+   END SUBROUTINE bdy_nmn 
+
 
 #else
@@ -366 +440 @@
       WRITE(*,*) 'bdy_orlanski_3d: You should not have seen this print! error?', kt
    END SUBROUTINE bdy_orlanski_3d
+   SUBROUTINE bdy_nmn( idx, igrd, phia )      ! Empty routine 
+      WRITE(*,*) 'bdy_nmn: You should not have seen this print! error?', kt 
+   END SUBROUTINE bdy_nmn 
 #endif