[Mitgcm-support] about kpp & flt

[mitgcm-support at dev.mitgcm.org]

Hello,

I have 2 questions about the mixlayer diagnostic in kpp_calc.F :

--------------------------------------
c     Determine mixed layer depth hMix as the shallowest depth at which
c     dB/dz exceeds 5.2e-5 s^-2.
            work1(i,j) = nzmax(i,j,bi,bj)
            DO k = 1, Nr
               IF ( k .LT. nzmax(i,j,bi,bj) .AND.
     &              dbloc(i,j,k) / drC(k+1) .GT. dB_dz )
     &              work1(i,j) = k
            END DO
 
c     Linearly interpolate to find hMix.
            k = work1(i,j)
            IF ( k .EQ. 0 .OR. nzmax(i,j,bi,bj) .EQ. 1 ) THEN
               zRef(i,j) = p0
            ELSEIF ( k .EQ. 1) THEN
               dBdz2 = dbloc(i,j,1) / drC(2)
               zRef(i,j) = drF(1) * dB_dz / dBdz2
            ELSEIF ( k .LT. nzmax(i,j,bi,bj) ) THEN
               dBdz1 = dbloc(i,j,k-1) / drC(k  )
               dBdz2 = dbloc(i,j,k  ) / drC(k+1)
>              zRef(i,j) = rF(k) + drF(k) * (dB_dz - dBdz1) /
>    &                     MAX ( phepsi, dBdz2 - dBdz1 )
            ELSE
               zRef(i,j) = rF(k+1)
            ENDIF                                                                       
--------------------------------------
1) I wonder if in a case where 
the vertical profile of db/dz shows 2 minimum below dB_dz,
 e.g., 1 in the mix layer, stratified just below, then very low 
 stratification (< dB_dz) at some depth in the deep ocean, 
 and finally stratified near the bottom ;
the first loop might capture the 2nd stratification and not the 
1rst one as the comment suggest "the shallowest depth ..."

2) since dBdz1 has to be < dB_dz
  dB_dz - dBdz1 > 0 ; MAX ( phepsi, dBdz2 - dBdz1 ) > 0 
so that drF(k)*(dB_dz - dBdz1) /MAX ( phepsi, dBdz2 - dBdz1 ) 
is also > 0 ; rF(k) < 0 so that |zRef(i,j)| < |rF(k)|
and this does not seems logical to me.

I would expect to find zRef(i,j) = rF(k) - drF(k) * ...

-----------------

and one question about  pkg/flt/flt_functions.F 

      do k=0,Nr
         zlo = rC(k)
         if (k.eq.0)  zlo = rC(k)-delZ(1)*0.5
         zhi = rC(k+1)
         if (k.eq.Nr) zhi = rC(k)+delZ(Nr)*0.5       

Apart from using delZ instead of drF, 3 out-of-bound index 
(in 4 lines !) and rC(k)-delZ(1)*0.5 instead of rC(1)-delZ(1)*0.5,
it seems to me that the signs in front of the 2 delZ are wrong.

I propose to replace this part with the following :

      do k=0,Nr
         if (k.eq.0) then
           zlo = rF(1)
         else
           zlo = rC(k)
         endif
         if (k.eq.Nr) then
           zhi = rF(k+1)
         else
           zhi = rC(k+1)
         endif

Is it OK like this ?

Jean-Michel