Changeset 5038 for branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfmxl.F90

Timestamp:

2015-01-20T15:26:13+01:00 (9 years ago)

Author:

jamesharle

Message:

Merging branch with HEAD of the trunk

File:

• branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfmxl.F90 (modified) (6 diffs)

branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfmxl.F90

@@ -6 +6 @@
    !! History :  1.0  ! 2003-08  (G. Madec)  original code
    !!            3.2  ! 2009-07  (S. Masson, G. Madec)  IOM + merge of DO-loop
+   !!            3.7  ! 2012-03  (G. Madec)  make public the density criteria for trdmxl 
+   !!             -   ! 2014-02  (F. Roquet)  mixed layer depth calculated using N2 instead of rhop 
    !!----------------------------------------------------------------------
    !!   zdf_mxl      : Compute the turbocline and mixed layer depths.
@@ -14 +16 @@
    USE in_out_manager  ! I/O manager
    USE prtctl          ! Print control
+   USE phycst          ! physical constants
    USE iom             ! I/O library
    USE lib_mpp         ! MPP library
@@ -25 +28 @@
    PUBLIC   zdf_mxl       ! called by step.F90
 
-   REAL(wp), PUBLIC ::   rho_c = 0.01_wp    ! density criterion for mixed layer depth
-   REAL(wp), PUBLIC ::   avt_c = 5.e-4_wp   ! Kz criterion for the turbocline depth
-
    INTEGER , PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   nmln    !: number of level in the mixed layer (used by TOP)
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hmld    !: mixing layer depth (turbocline)      [m]
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hmlp    !: mixed layer depth  (rho=rho0+zdcrit) [m]
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hmlpt   !: mixed layer depth at t-points        [m]
+
+   REAL(wp), PUBLIC ::   rho_c = 0.01_wp    !: density criterion for mixed layer depth
+   REAL(wp)         ::   avt_c = 5.e-4_wp   ! Kz criterion for the turbocline depth
 
    !! * Substitutions
@@ -70 +73 @@
       !!      eddy diffusivity coefficient (resulting from the vertical physics
       !!      alone, not the isopycnal part, see trazdf.F) fall below a given
-      !!      value defined locally (avt_c here taken equal to 5 cm/s2)
+      !!      value defined locally (avt_c here taken equal to 5 cm/s2 by default)
       !!
       !! ** Action  :   nmln, hmld, hmlp, hmlpt
       !!----------------------------------------------------------------------
       INTEGER, INTENT(in) ::   kt   ! ocean time-step index
-      !!
-      INTEGER  ::   ji, jj, jk          ! dummy loop indices
-      INTEGER  ::   iikn, iiki          ! temporary integer within a do loop
-      INTEGER, POINTER, DIMENSION(:,:) ::   imld                ! temporary workspace
+      !
+      INTEGER  ::   ji, jj, jk   ! dummy loop indices
+      INTEGER  ::   iikn, iiki, ikt, imkt   ! local integer
+      REAL(wp) ::   zN2_c        ! local scalar
+      INTEGER, POINTER, DIMENSION(:,:) ::   imld   ! 2D workspace
       !!----------------------------------------------------------------------
       !
@@ -94 +98 @@
 
       ! w-level of the mixing and mixed layers
-      nmln(:,:) = mbkt(:,:) + 1        ! Initialization to the number of w ocean point
-      imld(:,:) = mbkt(:,:) + 1
-      DO jk = jpkm1, nlb10, -1         ! from the bottom to nlb10
+      nmln(:,:)  = nlb10               ! Initialization to the number of w ocean point
+      hmlp(:,:)  = 0._wp               ! here hmlp used as a dummy variable, integrating vertically N^2
+      zN2_c = grav * rho_c * r1_rau0   ! convert density criteria into N^2 criteria
+      DO jk = nlb10, jpkm1
+         DO jj = 1, jpj                ! Mixed layer level: w-level 
+            DO ji = 1, jpi
+               ikt = mbkt(ji,jj)
+               hmlp(ji,jj) = hmlp(ji,jj) + MAX( rn2b(ji,jj,jk) , 0._wp ) * fse3w(ji,jj,jk)
+               IF( hmlp(ji,jj) < zN2_c )   nmln(ji,jj) = MIN( jk , ikt ) + 1   ! Mixed layer level
+            END DO
+         END DO
+      END DO
+      !
+      ! w-level of the turbocline
+      imld(:,:) = mbkt(:,:) + 1        ! Initialization to the number of w ocean point
+      DO jk = jpkm1, nlb10, -1         ! from the bottom to nlb10 
          DO jj = 1, jpj
             DO ji = 1, jpi
-               IF( rhop(ji,jj,jk) > rhop(ji,jj,nla10) + rho_c )   nmln(ji,jj) = jk      ! Mixed layer
-               IF( avt (ji,jj,jk) < avt_c                     )   imld(ji,jj) = jk      ! Turbocline 
+               imkt = mikt(ji,jj)
+               IF( avt (ji,jj,jk) < avt_c )   imld(ji,jj) = MAX( imkt, jk )      ! Turbocline 
             END DO
          END DO
@@ -109 +126 @@
             iiki = imld(ji,jj)
             iikn = nmln(ji,jj)
-            hmld (ji,jj) = fsdepw(ji,jj,iiki  ) * tmask(ji,jj,1)    ! Turbocline depth 
-            hmlp (ji,jj) = fsdepw(ji,jj,iikn  ) * tmask(ji,jj,1)    ! Mixed layer depth
-            hmlpt(ji,jj) = fsdept(ji,jj,iikn-1)                     ! depth of the last T-point inside the mixed layer
+            imkt = mikt(ji,jj)
+            hmld (ji,jj) = ( fsdepw(ji,jj,iiki  ) - fsdepw(ji,jj,imkt )            )   * ssmask(ji,jj)    ! Turbocline depth 
+            hmlp (ji,jj) = ( fsdepw(ji,jj,iikn  ) - fsdepw(ji,jj,MAX( imkt,nla10 ) ) ) * ssmask(ji,jj)    ! Mixed layer depth
+            hmlpt(ji,jj) = ( fsdept(ji,jj,iikn-1) - fsdepw(ji,jj,imkt )            )   * ssmask(ji,jj)    ! depth of the last T-point inside the mixed layer
          END DO
       END DO