New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 12068 for NEMO/branches/2019/UKMO_MERGE_2019/src/OCE/ISF/isfcpl.F90 – NEMO

Ignore:
Timestamp:
2019-12-05T13:18:21+01:00 (4 years ago)
Author:
davestorkey
Message:

2019/UKMO_MERGE_2019 : Merging in changes from ENHANCE-02_ISF_nemo.

Location:
NEMO/branches/2019/UKMO_MERGE_2019/src/OCE/ISF
Files:
1 edited
1 copied

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2019/UKMO_MERGE_2019/src/OCE/ISF/isfcpl.F90

    r11931 r12068  
    4646   !!---------------------------------------------------------------------- 
    4747CONTAINS 
    48    SUBROUTINE isfcpl_init() 
     48   SUBROUTINE isfcpl_init(Kbb, Kmm, Kaa) 
    4949      !!--------------------------------------------------------------------- 
    5050      !! 
     51      !!--------------------------------------------------------------------- 
     52      INTEGER, INTENT(in) :: Kbb, Kmm, Kaa      ! ocean time level indices 
    5153      !!--------------------------------------------------------------------- 
    5254      INTEGER :: id 
     
    7476      ELSE 
    7577         ! extrapolation ssh 
    76          CALL isfcpl_ssh() 
     78         CALL isfcpl_ssh(Kbb, Kmm, Kaa) 
    7779         ! 
    7880         ! extrapolation tracer properties 
    79          CALL isfcpl_tra() 
     81         CALL isfcpl_tra(Kmm) 
    8082         ! 
    8183         ! correction of the horizontal divergence and associated temp. and salt content flux 
    8284         ! Need to : - include in the cpl cons the risfcpl_vol/tsc contribution 
    8385         !           - decide how to manage thickness level change in conservation 
    84          CALL isfcpl_vol() 
     86         CALL isfcpl_vol(Kmm) 
    8587         ! 
    8688         ! apply the 'conservation' method 
    87          IF ( ln_isfcpl_cons ) CALL isfcpl_cons() 
     89         IF ( ln_isfcpl_cons ) CALL isfcpl_cons(Kmm) 
    8890         ! 
    8991      END IF 
    9092      ! 
    9193      ! mask velocity properly (mask used in restart not compatible with new mask) 
    92       un(:,:,:) = un(:,:,:) * umask(:,:,:) 
    93       vn(:,:,:) = vn(:,:,:) * vmask(:,:,:) 
     94      uu(:,:,:,Kmm) = uu(:,:,:,Kmm) * umask(:,:,:) 
     95      vv(:,:,:,Kmm) = vv(:,:,:,Kmm) * vmask(:,:,:) 
    9496      ! 
    9597      ! all before fields set to now values 
    96       tsb  (:,:,:,:) = tsn  (:,:,:,:) 
    97       ub   (:,:,:)   = un   (:,:,:) 
    98       vb   (:,:,:)   = vn   (:,:,:) 
    99       sshb (:,:)     = sshn (:,:) 
    100       e3t_b(:,:,:)   = e3t_n(:,:,:) 
     98      ts  (:,:,:,:,Kbb) = ts  (:,:,:,:,Kmm) 
     99      uu   (:,:,:,Kbb)   = uu   (:,:,:,Kmm) 
     100      vv   (:,:,:,Kbb)   = vv   (:,:,:,Kmm) 
     101      ssh (:,:,Kbb)     = ssh (:,:,Kmm) 
     102      e3t(:,:,:,Kbb)   = e3t(:,:,:,Kmm) 
    101103  
    102104      ! prepare writing restart 
     
    111113   END SUBROUTINE isfcpl_init 
    112114   !  
    113    SUBROUTINE isfcpl_rst_write(kt) 
     115   SUBROUTINE isfcpl_rst_write(kt, Kmm) 
    114116      !!--------------------------------------------------------------------- 
    115117      !!                   ***  ROUTINE iscpl_rst_write  *** 
     
    119121      !!-------------------------- IN  -------------------------------------- 
    120122      INTEGER, INTENT(in) :: kt 
     123      INTEGER, INTENT(in) :: Kmm    ! ocean time level index 
    121124      !!---------------------------------------------------------------------- 
    122125      ! 
     
    124127      CALL iom_rstput( kt, nitrst, numrow, 'tmask'  , tmask , ldxios = lwxios ) 
    125128      CALL iom_rstput( kt, nitrst, numrow, 'ssmask' , ssmask, ldxios = lwxios ) 
    126       CALL iom_rstput( kt, nitrst, numrow, 'e3t_n'  , e3t_n , ldxios = lwxios ) 
    127       CALL iom_rstput( kt, nitrst, numrow, 'e3u_n'  , e3u_n , ldxios = lwxios ) 
    128       CALL iom_rstput( kt, nitrst, numrow, 'e3v_n'  , e3v_n , ldxios = lwxios ) 
    129       CALL iom_rstput( kt, nitrst, numrow, 'gdepw_n', gdepw_n , ldxios = lwxios ) 
     129      CALL iom_rstput( kt, nitrst, numrow, 'e3t_n'  , e3t(:,:,:,Kmm) , ldxios = lwxios ) 
     130      CALL iom_rstput( kt, nitrst, numrow, 'e3u_n'  , e3u(:,:,:,Kmm) , ldxios = lwxios ) 
     131      CALL iom_rstput( kt, nitrst, numrow, 'e3v_n'  , e3v(:,:,:,Kmm) , ldxios = lwxios ) 
     132      CALL iom_rstput( kt, nitrst, numrow, 'gdepw_n', gdepw(:,:,:,Kmm) , ldxios = lwxios ) 
    130133      IF( lwxios ) CALL iom_swap( cxios_context ) 
    131134      ! 
    132135   END SUBROUTINE isfcpl_rst_write 
    133136 
    134    SUBROUTINE isfcpl_ssh() 
     137   SUBROUTINE isfcpl_ssh(Kbb, Kmm, Kaa) 
    135138      !!----------------------------------------------------------------------  
    136139      !!                   ***  ROUTINE iscpl_ssh  *** 
     
    142145      !!---------------------------------------------------------------------- 
    143146      !! 
     147      INTEGER, INTENT(in) :: Kbb, Kmm, Kaa    ! ocean time level indices 
     148      !!---------------------------------------------------------------------- 
    144149      INTEGER :: ji, jj, jd, jk      !! loop index 
    145150      INTEGER :: jip1, jim1, jjp1, jjm1 
     
    154159      ! rude average of the closest neigbourgs (e1e2t not taking into account) 
    155160      ! 
    156       zssh(:,:)     = sshn(:,:) 
     161      zssh(:,:)     = ssh(:,:,Kmm) 
    157162      zssmask0(:,:) = zssmask_b(:,:) 
    158163      ! 
     
    168173               ! 
    169174               IF (zdssmask(ji,jj) == 1._wp .AND. zsummsk /= 0._wp) THEN 
    170                   sshn(ji,jj)=( zssh(jip1,jj)*zssmask0(jip1,jj)     & 
     175                  ssh(ji,jj,Kmm)=( zssh(jip1,jj)*zssmask0(jip1,jj)     & 
    171176                  &           + zssh(jim1,jj)*zssmask0(jim1,jj)     & 
    172177                  &           + zssh(ji,jjp1)*zssmask0(ji,jjp1)     & 
     
    177182         END DO 
    178183         ! 
    179          zssh(:,:) = sshn(:,:) 
     184         zssh(:,:) = ssh(:,:,Kmm) 
    180185         zssmask0(:,:) = zssmask_b(:,:) 
    181186         ! 
     
    184189      END DO 
    185190      ! 
    186       ! update sshn 
    187       sshn(:,:) = zssh(:,:) * ssmask(:,:) 
    188       ! 
    189       sshb(:,:) = sshn(:,:) 
    190       ! 
    191       IF ( ln_isfdebug ) CALL debug('isfcpl_ssh: sshn',sshn(:,:)) 
     191      ! update ssh(:,:,Kmm) 
     192      ssh(:,:,Kmm) = zssh(:,:) * ssmask(:,:) 
     193      ! 
     194      ssh(:,:,Kbb) = ssh(:,:,Kmm) 
     195      ! 
     196      IF ( ln_isfdebug ) CALL debug('isfcpl_ssh: sshn',ssh(:,:,Kmm)) 
    192197      ! 
    193198      ! recompute the vertical scale factor, depth and water thickness 
    194       IF(lwp) write(numout,*) 'isfcpl_ssh : recompute scale factor from sshn (new wet cell)' 
     199      IF(lwp) write(numout,*) 'isfcpl_ssh : recompute scale factor from ssh (new wet cell,Kmm)' 
    195200      IF(lwp) write(numout,*) '~~~~~~~~~~~' 
    196201      DO jk = 1, jpk 
    197          e3t_n(:,:,jk) =  e3t_0(:,:,jk) * ( ht_0(:,:) + sshn(:,:) ) & 
     202         e3t(:,:,jk,Kmm) =  e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kmm) ) & 
    198203             &                          / ( ht_0(:,:) + 1._wp - ssmask(:,:) ) * tmask(:,:,jk)   & 
    199204             &          + e3t_0(:,:,jk)                               * (1._wp -tmask(:,:,jk)) 
    200205      END DO 
    201       e3t_b(:,:,:) = e3t_n(:,:,:) 
    202       CALL dom_vvl_zgr() 
     206      e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 
     207      CALL dom_vvl_zgr(Kbb, Kmm, Kaa) 
    203208      ! 
    204209   END SUBROUTINE isfcpl_ssh 
    205210 
    206    SUBROUTINE isfcpl_tra() 
     211   SUBROUTINE isfcpl_tra(Kmm) 
    207212      !!----------------------------------------------------------------------  
    208213      !!                   ***  ROUTINE iscpl_tra  *** 
     
    213218      !! 
    214219      !!---------------------------------------------------------------------- 
     220      INTEGER, INTENT(in) :: Kmm    ! ocean time level index 
     221      !!---------------------------------------------------------------------- 
    215222      REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztmask_b 
    216223      !REAL(wp), DIMENSION(:,:,:  ), INTENT(in ) :: pdepw_b                         !! depth w before 
     
    237244      !bugged : to be corrected (PM) 
    238245      ! back up original t/s/mask 
    239       !tsb (:,:,:,:) = tsn(:,:,:,:) 
     246      !tsb (:,:,:,:) = ts(:,:,:,:,Kmm) 
    240247      ! 
    241248     ! compute new T/S (interpolation) if vvl only for common wet cell in before and after wmask 
     
    248255! 
    249256!                     !compute weight 
    250 !                     zdzp1 = MAX(0._wp,pdepw_b(ji,jj,jk+1) - gdepw_n(ji,jj,jk+1)) 
    251 !                     zdzm1 = MAX(0._wp,gdepw_n(ji,jj,jk  ) - pdepw_b(ji,jj,jk  )) 
    252 !                     zdz   = e3t_n(ji,jj,jk) - zdzp1 - zdzm1 ! if isf : e3t = gdepw_n(ji,jj,jk+1)- gdepw_n(ji,jj,jk) 
     257!                     zdzp1 = MAX(0._wp,pdepw_b(ji,jj,jk+1) - gdepw(ji,jj,jk+1,Kmm)) 
     258!                     zdzm1 = MAX(0._wp,gdepw(ji,jj,jk  ,Kmm) - pdepw_b(ji,jj,jk  )) 
     259!                     zdz   = e3t(ji,jj,jk,Kmm) - zdzp1 - zdzm1 ! if isf : e3t = gdepw(ji,jj,jk+1,Kmm)- gdepw(ji,jj,jk,Kmm) 
    253260! 
    254261!                     IF (zdz .LT. 0._wp) THEN 
     
    256263!                     END IF 
    257264! 
    258 !                     tsn(ji,jj,jk,jp_tem) = ( zdzp1*tsb(ji,jj,jk+1,jp_tem) & 
    259 !                        &                   + zdz  *tsb(ji,jj,jk  ,jp_tem) & 
    260 !                        &                   + zdzm1*tsb(ji,jj,jk-1,jp_tem) )/e3t_n(ji,jj,jk) 
     265!                     ts(ji,jj,jk,jp_tem,Kmm) = ( zdzp1*ts(ji,jj,jk+1,jp_tem,Kbb) & 
     266!                        &                   + zdz  *ts(ji,jj,jk  ,jp_tem,Kbb) & 
     267!                        &                   + zdzm1*ts(ji,jj,jk-1,jp_tem,Kbb) )/e3t(ji,jj,jk,Kmm) 
    261268! 
    262 !                     tsn(ji,jj,jk,jp_sal) = ( zdzp1*tsb(ji,jj,jk+1,jp_sal) & 
    263 !                        &                   + zdz  *tsb(ji,jj,jk  ,jp_sal) & 
    264 !                        &                   + zdzm1*tsb(ji,jj,jk-1,jp_sal) )/e3t_n(ji,jj,jk) 
     269!                     ts(ji,jj,jk,jp_sal,Kmm) = ( zdzp1*ts(ji,jj,jk+1,jp_sal,Kbb) & 
     270!                        &                   + zdz  *ts(ji,jj,jk  ,jp_sal,Kbb) & 
     271!                        &                   + zdzm1*ts(ji,jj,jk-1,jp_sal,Kbb) )/e3t(ji,jj,jk,Kmm) 
    265272! 
    266273!                  END IF 
     
    270277!      END IF 
    271278 
    272       zts0(:,:,:,:)  = tsn(:,:,:,:) 
     279      zts0(:,:,:,:)  = ts(:,:,:,:,Kmm) 
    273280      ztmask0(:,:,:) = ztmask_b(:,:,:) 
    274281      ztmask1(:,:,:) = ztmask_b(:,:,:) 
     
    294301                     ! 
    295302                     ! horizontal basic extrapolation 
    296                      tsn(ji,jj,jk,1)=( zts0(jip1,jj  ,jk,1) * ztmask0(jip1,jj  ,jk) & 
     303                     ts(ji,jj,jk,1,Kmm)=( zts0(jip1,jj  ,jk,1) * ztmask0(jip1,jj  ,jk) & 
    297304                     &               + zts0(jim1,jj  ,jk,1) * ztmask0(jim1,jj  ,jk) & 
    298305                     &               + zts0(ji  ,jjp1,jk,1) * ztmask0(ji  ,jjp1,jk) & 
    299306                     &               + zts0(ji  ,jjm1,jk,1) * ztmask0(ji  ,jjm1,jk) ) / zsummsk 
    300                      tsn(ji,jj,jk,2)=( zts0(jip1,jj  ,jk,2) * ztmask0(jip1,jj  ,jk) & 
     307                     ts(ji,jj,jk,2,Kmm)=( zts0(jip1,jj  ,jk,2) * ztmask0(jip1,jj  ,jk) & 
    301308                     &               + zts0(jim1,jj  ,jk,2) * ztmask0(jim1,jj  ,jk) & 
    302309                     &               + zts0(ji  ,jjp1,jk,2) * ztmask0(ji  ,jjp1,jk) & 
     
    316323                     zsummsk = ztmask0(ji,jj,jkm1) + ztmask0(ji,jj,jkp1) 
    317324                     IF (zdmask(ji,jj) == 1._wp .AND. zsummsk /= 0._wp ) THEN 
    318                         tsn(ji,jj,jk,1)=( zts0(ji,jj,jkp1,1)*ztmask0(ji,jj,jkp1)     & 
     325                        ts(ji,jj,jk,1,Kmm)=( zts0(ji,jj,jkp1,1)*ztmask0(ji,jj,jkp1)     & 
    319326                        &               + zts0(ji,jj,jkm1,1)*ztmask0(ji,jj,jkm1)) / zsummsk 
    320                         tsn(ji,jj,jk,2)=( zts0(ji,jj,jkp1,2)*ztmask0(ji,jj,jkp1)     & 
     327                        ts(ji,jj,jk,2,Kmm)=( zts0(ji,jj,jkp1,2)*ztmask0(ji,jj,jkp1)     & 
    321328                        &               + zts0(ji,jj,jkm1,2)*ztmask0(ji,jj,jkm1)) / zsummsk 
    322329                        ! 
     
    330337         ! 
    331338         ! update temperature and salinity and mask 
    332          zts0(:,:,:,:)  = tsn(:,:,:,:) 
     339         zts0(:,:,:,:)  = ts(:,:,:,:,Kmm) 
    333340         ztmask0(:,:,:) = ztmask1(:,:,:) 
    334341         ! 
     
    337344      END DO  ! nn_drown 
    338345      ! 
    339       ! mask new tsn field 
    340       tsn(:,:,:,jp_tem) = zts0(:,:,:,jp_tem) * tmask(:,:,:) 
    341       tsn(:,:,:,jp_sal) = zts0(:,:,:,jp_sal) * tmask(:,:,:) 
     346      ! mask new ts(:,:,:,:,Kmm) field 
     347      ts(:,:,:,jp_tem,Kmm) = zts0(:,:,:,jp_tem) * tmask(:,:,:) 
     348      ts(:,:,:,jp_sal,Kmm) = zts0(:,:,:,jp_sal) * tmask(:,:,:) 
    342349      ! 
    343350      ! sanity check 
     
    347354         DO jj = 1,jpj 
    348355            DO ji = 1,jpi 
    349                IF (tmask(ji,jj,jk) == 1._wp .AND. tsn(ji,jj,jk,2) == 0._wp)              & 
     356               IF (tmask(ji,jj,jk) == 1._wp .AND. ts(ji,jj,jk,2,Kmm) == 0._wp)              & 
    350357                  &   CALL ctl_stop('STOP', 'failing to fill all new weet cell,     & 
    351358                  &                          try increase nn_drown or activate XXXX & 
     
    357364   END SUBROUTINE isfcpl_tra 
    358365 
    359    SUBROUTINE isfcpl_vol() 
     366   SUBROUTINE isfcpl_vol(Kmm) 
    360367      !!----------------------------------------------------------------------  
    361368      !!                   ***  ROUTINE iscpl_vol  *** 
     
    370377      !!---------------------------------------------------------------------- 
    371378      !! 
     379      INTEGER, INTENT(in) :: Kmm    ! ocean time level index 
     380      !!---------------------------------------------------------------------- 
    372381      INTEGER :: ji, jj, jk  
    373382      INTEGER :: ikb, ikt 
     
    388397         DO jj = 2, jpjm1 
    389398            DO ji = 2, jpim1 
    390                zqvolb(ji,jj,jk) =  (   e2u(ji,jj) * ze3u_b(ji,jj,jk) * un(ji,jj,jk) - e2u(ji-1,jj  ) * ze3u_b(ji-1,jj  ,jk) * un(ji-1,jj  ,jk)    & 
    391                   &                  + e1v(ji,jj) * ze3v_b(ji,jj,jk) * vn(ji,jj,jk) - e1v(ji  ,jj-1) * ze3v_b(ji  ,jj-1,jk) * vn(ji  ,jj-1,jk)  ) & 
     399               zqvolb(ji,jj,jk) =  (   e2u(ji,jj) * ze3u_b(ji,jj,jk) * uu(ji,jj,jk,Kmm) - e2u(ji-1,jj  ) * ze3u_b(ji-1,jj  ,jk) * uu(ji-1,jj  ,jk,Kmm)    & 
     400                  &                  + e1v(ji,jj) * ze3v_b(ji,jj,jk) * vv(ji,jj,jk,Kmm) - e1v(ji  ,jj-1) * ze3v_b(ji  ,jj-1,jk) * vv(ji  ,jj-1,jk,Kmm)  ) & 
    392401                  &                * ztmask_b(ji,jj,jk) 
    393402            END DO 
     
    397406         ! properly mask velocity  
    398407         ! (velocity are still mask with old mask at this stage) 
    399          un(:,:,jk) = un(:,:,jk) * umask(:,:,jk) 
    400          vn(:,:,jk) = vn(:,:,jk) * vmask(:,:,jk) 
     408         uu(:,:,jk,Kmm) = uu(:,:,jk,Kmm) * umask(:,:,jk) 
     409         vv(:,:,jk,Kmm) = vv(:,:,jk,Kmm) * vmask(:,:,jk) 
    401410         ! compute volume flux divergence after coupling 
    402411         DO jj = 2, jpjm1 
    403412            DO ji = 2, jpim1 
    404                zqvoln(ji,jj,jk) = (   e2u(ji,jj) * e3u_n(ji,jj,jk) * un(ji,jj,jk) - e2u(ji-1,jj  ) * e3u_n(ji-1,jj  ,jk) * un(ji-1,jj  ,jk)    & 
    405                   &                 + e1v(ji,jj) * e3v_n(ji,jj,jk) * vn(ji,jj,jk) - e1v(ji  ,jj-1) * e3v_n(ji  ,jj-1,jk) * vn(ji  ,jj-1,jk)  ) & 
     413               zqvoln(ji,jj,jk) = (   e2u(ji,jj) * e3u(ji,jj,jk,Kmm) * uu(ji,jj,jk,Kmm) - e2u(ji-1,jj  ) * e3u(ji-1,jj  ,jk,Kmm) * uu(ji-1,jj  ,jk,Kmm)    & 
     414                  &                 + e1v(ji,jj) * e3v(ji,jj,jk,Kmm) * vv(ji,jj,jk,Kmm) - e1v(ji  ,jj-1) * e3v(ji  ,jj-1,jk,Kmm) * vv(ji  ,jj-1,jk,Kmm)  ) & 
    406415                  &               * tmask(ji,jj,jk) 
    407416            END DO 
     
    428437      CALL lbc_lnk( 'iscpl', risfcpl_vol, 'T', 1. ) 
    429438      ! 
    430       ! 3.0: set total correction (div, tra, ssh) 
     439      ! 3.0: set total correction (div, tr(:,:,:,:,Krhs), ssh) 
    431440      ! 
    432441      ! 3.1: mask volume flux divergence correction 
    433442      risfcpl_vol(:,:,:) = risfcpl_vol(:,:,:) * tmask(:,:,:) 
    434443      ! 
    435       ! 3.2: get 3d tra increment to apply at the first time step 
    436       ! temperature and salt content flux computed using local tsn  
     444      ! 3.2: get 3d tr(:,:,:,:,Krhs) increment to apply at the first time step 
     445      ! temperature and salt content flux computed using local ts(:,:,:,:,Kmm)  
    437446      ! (very simple advection scheme) 
    438447      ! (>0 out) 
    439       risfcpl_tsc(:,:,:,jp_tem) = -risfcpl_vol(:,:,:) * tsn(:,:,:,jp_tem) 
    440       risfcpl_tsc(:,:,:,jp_sal) = -risfcpl_vol(:,:,:) * tsn(:,:,:,jp_sal) 
     448      risfcpl_tsc(:,:,:,jp_tem) = -risfcpl_vol(:,:,:) * ts(:,:,:,jp_tem,Kmm) 
     449      risfcpl_tsc(:,:,:,jp_sal) = -risfcpl_vol(:,:,:) * ts(:,:,:,jp_sal,Kmm) 
    441450      ! 
    442451      ! 3.3: ssh correction (for dynspg_ts) 
     
    448457   END SUBROUTINE isfcpl_vol 
    449458 
    450    SUBROUTINE isfcpl_cons() 
     459   SUBROUTINE isfcpl_cons(Kmm) 
    451460      !!----------------------------------------------------------------------  
    452461      !!                   ***  ROUTINE iscpl_cons  *** 
     
    463472      TYPE(isfcons), DIMENSION(:),ALLOCATABLE :: zisfpts ! list of point receiving a correction 
    464473      ! 
     474      !!---------------------------------------------------------------------- 
     475      INTEGER, INTENT(in) :: Kmm    ! ocean time level index 
     476      !!---------------------------------------------------------------------- 
    465477      INTEGER  ::   ji   , jj  , jk  , jproc          ! loop index 
    466478      INTEGER  ::   jip1 , jim1, jjp1, jjm1           ! dummy indices 
     
    513525 
    514526               ! volume diff 
    515                zdvol = e3t_n(ji,jj,jk) * tmask(ji,jj,jk) - ze3t_b(ji,jj,jk) * ztmask_b(ji,jj,jk) 
     527               zdvol = e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) - ze3t_b(ji,jj,jk) * ztmask_b(ji,jj,jk) 
    516528 
    517529               ! heat diff 
    518                zdtem = tsn (ji,jj,jk,jp_tem) *  e3t_n(ji,jj,jk) *  tmask  (ji,jj,jk)   & 
     530               zdtem = ts (ji,jj,jk,jp_tem,Kmm) *  e3t(ji,jj,jk,Kmm) *  tmask  (ji,jj,jk)   & 
    519531                     - zt_b(ji,jj,jk)        * ze3t_b(ji,jj,jk) * ztmask_b(ji,jj,jk) 
    520532 
    521533               ! salt diff 
    522                zdsal = tsn(ji,jj,jk,jp_sal) *  e3t_n(ji,jj,jk) *  tmask  (ji,jj,jk)   & 
     534               zdsal = ts(ji,jj,jk,jp_sal,Kmm) *  e3t(ji,jj,jk,Kmm) *  tmask  (ji,jj,jk)   & 
    523535                     - zs_b(ji,jj,jk)       * ze3t_b(ji,jj,jk) * ztmask_b(ji,jj,jk) 
    524536             
Note: See TracChangeset for help on using the changeset viewer.