Changeset 9250 for branches/2017/dev_merge_2017/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfdrg.F90

Timestamp:

2018-01-17T08:25:53+01:00 (6 years ago)

Author:

gm

Message:

dev_merge_2017 : HPC09_ZDF: drag finalization: move dyn_bfr routine in zdfdrg (renamed zdf_drg_exp) and delete dynbfr.F90 module - result unchanged

File:

• branches/2017/dev_merge_2017/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfdrg.F90 (modified) (4 diffs)

branches/2017/dev_merge_2017/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfdrg.F90

@@ -16 +16 @@
    !!----------------------------------------------------------------------
    !!   zdf_drg       : update bottom friction coefficient (non-linear bottom friction only)
+   !!   zdf_drg_exp   : compute the top & bottom friction in explicit case
    !!   zdf_drg_init  : read in namdrg namelist and control the bottom friction parameters.
    !!       drg_init  :
@@ -23 +24 @@
    USE dom_oce        ! ocean space and time domain variables
    USE zdf_oce        ! ocean vertical physics variables
+   USE trd_oce        ! trends: ocean variables
+   USE trddyn         ! trend manager: dynamics
    !
    USE in_out_manager ! I/O manager
@@ -34 +37 @@
 
    PUBLIC   zdf_drg         ! called by zdf_phy
+   PUBLIC   zdf_drg_exp     ! called by dyn_zdf
    PUBLIC   zdf_drg_init    ! called by zdf_phy_init
 
@@ -140 +144 @@
 
 
+   SUBROUTINE zdf_drg_exp( kt, pub, pvb, pua, pva )
+      !!----------------------------------------------------------------------
+      !!                  ***  ROUTINE zdf_drg_exp  ***
+      !!
+      !! ** Purpose :   compute and add the explicit top and bottom frictions.
+      !!
+      !! ** Method  :   in explicit case, 
+      !!
+      !!              NB: in implicit case the calculation is performed in dynzdf.F90
+      !!
+      !! ** Action  :   (pua,pva)   momentum trend increased by top & bottom friction trend
+      !!---------------------------------------------------------------------
+      INTEGER                         , INTENT(in   ) ::   kt         ! ocean time-step index
+      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pub, pvb   ! the two components of the before velocity
+      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pua, pva   ! the two components of the velocity tendency
+      !! 
+      INTEGER  ::   ji, jj       ! dummy loop indexes
+      INTEGER  ::   ikbu, ikbv   ! local integers
+      REAL(wp) ::   zm1_2dt      ! local scalar
+      REAL(wp) ::   zCdu, zCdv   !   -      -
+      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::   ztrdu, ztrdv
+      !!---------------------------------------------------------------------
+      !
+!!gm bug : time step is only rdt (not 2 rdt if euler start !)
+      zm1_2dt = - 1._wp / ( 2._wp * rdt )
+
+      IF( l_trddyn ) THEN      ! trends: store the input trends
+         ALLOCATE( ztrdu(jpi,jpj,jpk) , ztrdv(jpi,jpj,jpk) )
+         ztrdu(:,:,:) = pua(:,:,:)
+         ztrdv(:,:,:) = pva(:,:,:)
+      ENDIF
+
+      DO jj = 2, jpjm1
+         DO ji = 2, jpim1
+            ikbu = mbku(ji,jj)          ! deepest wet ocean u- & v-levels
+            ikbv = mbkv(ji,jj)
+            !
+            ! Apply stability criteria on absolute value  : abs(bfr/e3) < 1/(2dt) => bfr/e3 > -1/(2dt)
+            zCdu = 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) / e3u_n(ji,jj,ikbu)
+            zCdv = 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) / e3v_n(ji,jj,ikbv)
+            !
+            pua(ji,jj,ikbu) = pua(ji,jj,ikbu) + MAX(  zCdu , zm1_2dt  ) * pub(ji,jj,ikbu)
+            pva(ji,jj,ikbv) = pva(ji,jj,ikbv) + MAX(  zCdv , zm1_2dt  ) * pvb(ji,jj,ikbv)
+         END DO
+      END DO
+      !
+      IF( ln_isfcav ) THEN        ! ocean cavities
+         DO jj = 2, jpjm1
+            DO ji = 2, jpim1
+               ikbu = miku(ji,jj)          ! first wet ocean u- & v-levels
+               ikbv = mikv(ji,jj)
+               !
+               ! Apply stability criteria on absolute value  : abs(bfr/e3) < 1/(2dt) => bfr/e3 > -1/(2dt)
+               zCdu = 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) / e3u_n(ji,jj,ikbu)    ! NB: Cdtop masked
+               zCdv = 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) / e3v_n(ji,jj,ikbv)
+               !
+               pua(ji,jj,ikbu) = pua(ji,jj,ikbu) + MAX(  zCdu , zm1_2dt  ) * pub(ji,jj,ikbu)
+               pva(ji,jj,ikbv) = pva(ji,jj,ikbv) + MAX(  zCdv , zm1_2dt  ) * pvb(ji,jj,ikbv)
+           END DO
+         END DO
+      ENDIF
+      !
+      IF( l_trddyn ) THEN      ! trends: send trends to trddyn for further diagnostics
+         ztrdu(:,:,:) = pua(:,:,:) - ztrdu(:,:,:)
+         ztrdv(:,:,:) = pva(:,:,:) - ztrdv(:,:,:)
+         CALL trd_dyn( ztrdu(:,:,:), ztrdv(:,:,:), jpdyn_bfr, kt )
+         DEALLOCATE( ztrdu, ztrdv )
+      ENDIF
+      !                                          ! print mean trends (used for debugging)
+      IF(ln_ctl)   CALL prt_ctl( tab3d_1=pua, clinfo1=' bfr  - Ua: ', mask1=umask,               &
+         &                       tab3d_2=pva, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' )
+      !
+   END SUBROUTINE zdf_drg_exp
+
+
    SUBROUTINE zdf_drg_init
       !!----------------------------------------------------------------------