Ignore:
Timestamp:
2017-04-26T12:03:26+02:00 (3 years ago)
Author:
jchanut
Message:

Add zstar coordinate with AGRIF

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2017/dev_r7963_nemo_v3_6_AGRIF-3_AGRIFVVL/NEMOGCM/NEMO/NST_SRC/agrif_opa_update.F90

    r6204 r7971  
    1313   USE dynspg_oce 
    1414   USE zdf_oce        ! vertical physics: ocean variables  
     15   USE domvvl         ! Need interpolation routines  
    1516 
    1617   IMPLICIT NONE 
    1718   PRIVATE 
    1819 
    19    PUBLIC Agrif_Update_Tra, Agrif_Update_Dyn,Update_Scales 
     20   PUBLIC Agrif_Update_Tra, Agrif_Update_Dyn, Update_Scales, Agrif_Update_vvl 
     21 
    2022# if defined key_zdftke 
    2123   PUBLIC Agrif_Update_Tke 
     
    177179# endif /* key_zdftke */ 
    178180 
     181   RECURSIVE SUBROUTINE Agrif_Update_vvl( ) 
     182      !!--------------------------------------------- 
     183      !!   *** ROUTINE Agrif_Update_vvl *** 
     184      !!--------------------------------------------- 
     185      ! 
     186      IF (Agrif_Root()) RETURN 
     187      ! 
     188#if defined TWO_WAY   
     189      ! 
     190      IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update e3 from grid Number',Agrif_Fixed(), 'Step', Agrif_Nb_Step() 
     191      ! 
     192      Agrif_UseSpecialValueInUpdate = .FALSE. 
     193      Agrif_SpecialValueFineGrid = 0. 
     194      !  
     195# if ! defined DECAL_FEEDBACK 
     196      CALL Agrif_Update_Variable(e3t_id, procname=updatee3t) 
     197# else 
     198      CALL Agrif_Update_Variable(e3t_id, locupdate=(/1,0/), procname=updatee3t) 
     199# endif  
     200      ! 
     201      Agrif_UseSpecialValueInUpdate = .FALSE. 
     202      ! 
     203      CALL Agrif_ChildGrid_To_ParentGrid() 
     204      CALL dom_vvl_update_UVF 
     205      CALL Agrif_ParentGrid_To_ChildGrid() 
     206      ! 
     207      IF ( lk_agrif_doupd ) THEN ! Initialisation: do recursive update: 
     208         CALL Agrif_ChildGrid_To_ParentGrid() 
     209         CALL Agrif_Update_vvl() 
     210         CALL Agrif_ParentGrid_To_ChildGrid() 
     211      ENDIF 
     212      ! 
     213#endif 
     214      ! 
     215   END SUBROUTINE Agrif_Update_vvl 
     216 
     217   SUBROUTINE dom_vvl_update_UVF 
     218      !!--------------------------------------------- 
     219      !!       *** ROUTINE dom_vvl_update_UVF *** 
     220      !!--------------------------------------------- 
     221#  include "domzgr_substitute.h90" 
     222      !! 
     223      INTEGER :: jk 
     224      REAL(wp):: zcoef 
     225      !!--------------------------------------------- 
     226 
     227      IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Finalize e3 on grid Number', & 
     228                  & Agrif_Fixed(), 'Step', Agrif_Nb_Step() 
     229 
     230      CALL dom_vvl_interpol( fse3t_n(:,:,:), fse3u_n(:,:,:),  'U'  ) 
     231      CALL dom_vvl_interpol( fse3t_n(:,:,:), fse3v_n(:,:,:),  'V'  ) 
     232      CALL dom_vvl_interpol( fse3u_n(:,:,:), fse3f_n(:,:,:),  'F'  ) 
     233 
     234      CALL dom_vvl_interpol( fse3u_n(:,:,:), fse3uw_n(:,:,:), 'UW' ) 
     235      CALL dom_vvl_interpol( fse3v_n(:,:,:), fse3vw_n(:,:,:), 'VW' ) 
     236 
     237      ! Update total depths: 
     238      hu(:,:) = 0._wp                        ! Ocean depth at U-points 
     239      hv(:,:) = 0._wp                        ! Ocean depth at V-points 
     240      DO jk = 1, jpkm1 
     241         hu(:,:) = hu(:,:) + fse3u_n(:,:,jk) * umask(:,:,jk) 
     242         hv(:,:) = hv(:,:) + fse3v_n(:,:,jk) * vmask(:,:,jk) 
     243      END DO 
     244      !                                      ! Inverse of the local depth 
     245      hur(:,:) = 1._wp / ( hu(:,:) + 1._wp - umask_i(:,:) ) * umask_i(:,:) 
     246      hvr(:,:) = 1._wp / ( hv(:,:) + 1._wp - vmask_i(:,:) ) * vmask_i(:,:) 
     247 
     248#if ! defined key_dynspg_ts 
     249      IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN 
     250#else 
     251      IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN 
     252#endif 
     253         CALL dom_vvl_interpol( fse3t_b(:,:,:), fse3u_b(:,:,:),  'U'  ) 
     254         CALL dom_vvl_interpol( fse3t_b(:,:,:), fse3v_b(:,:,:),  'V'  ) 
     255 
     256         CALL dom_vvl_interpol( fse3u_b(:,:,:), fse3uw_b(:,:,:), 'UW' ) 
     257         CALL dom_vvl_interpol( fse3v_b(:,:,:), fse3vw_b(:,:,:), 'VW' ) 
     258 
     259         hu_b(:,:) = 0._wp                        ! Ocean depth at U-points 
     260         hv_b(:,:) = 0._wp                        ! Ocean depth at V-points 
     261         DO jk = 1, jpkm1 
     262            hu_b(:,:) = hu_b(:,:) + fse3u_b(:,:,jk) * umask(:,:,jk) 
     263            hv_b(:,:) = hv_b(:,:) + fse3v_b(:,:,jk) * vmask(:,:,jk) 
     264         END DO 
     265         !                                      ! Inverse of the local depth 
     266         hur_b(:,:) = 1._wp / ( hu_b(:,:) + 1._wp - umask_i(:,:) ) * umask_i(:,:) 
     267         hvr_b(:,:) = 1._wp / ( hv_b(:,:) + 1._wp - vmask_i(:,:) ) * vmask_i(:,:) 
     268      ENDIF 
     269 
     270      ! 
     271   END SUBROUTINE dom_vvl_update_UVF 
     272 
    179273   SUBROUTINE updateTS( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before ) 
    180274      !!--------------------------------------------- 
    181275      !!           *** ROUTINE updateT *** 
    182276      !!--------------------------------------------- 
    183 #  include "domzgr_substitute.h90" 
    184277      INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2 
    185278      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres 
     
    188281      INTEGER :: ji,jj,jk,jn 
    189282      !!--------------------------------------------- 
     283      ! 
    190284      ! 
    191285      IF (before) THEN 
     
    194288               DO jj=j1,j2 
    195289                  DO ji=i1,i2 
    196                      tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn) 
     290!> jc tmp 
     291                     tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn)  * fse3t_n(ji,jj,jk) / e3t_0(ji,jj,jk) 
     292!                     tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn)  * fse3t_n(ji,jj,jk) 
     293!< jc tmp 
    197294                  END DO 
    198295               END DO 
     
    200297         END DO 
    201298      ELSE 
     299!> jc tmp 
     300         DO jn = n1,n2 
     301            tabres(i1:i2,j1:j2,k1:k2,jn) =  tabres(i1:i2,j1:j2,k1:k2,jn) * e3t_0(i1:i2,j1:j2,k1:k2) 
     302         ENDDO 
     303!< jc tmp 
     304 
    202305         IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN 
    203306            ! Add asselin part 
     
    209312                           tsb(ji,jj,jk,jn) = tsb(ji,jj,jk,jn) &  
    210313                                 & + atfp * ( tabres(ji,jj,jk,jn) & 
    211                                  &             - tsn(ji,jj,jk,jn) ) * tmask(ji,jj,jk) 
     314                                 &             - tsn(ji,jj,jk,jn)*fse3t_a(ji,jj,jk) ) &  
     315                                 & * tmask(ji,jj,jk) / fse3t_b(ji,jj,jk) 
    212316                        ENDIF 
    213317                     ENDDO 
     
    221325                  DO ji=i1,i2 
    222326                     IF( tabres(ji,jj,jk,jn) .NE. 0. ) THEN  
    223                         tsn(ji,jj,jk,jn) = tabres(ji,jj,jk,jn) * tmask(ji,jj,jk) 
     327#if defined key_vvl 
     328                        tsn(ji,jj,jk,jn) = tabres(ji,jj,jk,jn) / fse3t_n(ji,jj,jk) 
     329#else 
     330                        tsn(ji,jj,jk,jn) = tabres(ji,jj,jk,jn) 
     331#endif 
    224332                     END IF 
    225333                  END DO 
     
    260368            DO jj=j1,j2 
    261369               DO ji=i1,i2 
    262                   tabres(ji,jj,jk) = tabres(ji,jj,jk) / e2u(ji,jj) / fse3u_n(ji,jj,jk) 
     370                  tabres(ji,jj,jk) = tabres(ji,jj,jk) / e2u(ji,jj)  
    263371                  ! 
    264372                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 
    265373                     ub(ji,jj,jk) = ub(ji,jj,jk) &  
    266                            & + atfp * ( tabres(ji,jj,jk) - un(ji,jj,jk) ) * umask(ji,jj,jk) 
     374                           & + atfp * ( tabres(ji,jj,jk)                     &  
     375                           &              - un(ji,jj,jk)*fse3u_a(ji,jj,jk) ) &  
     376                           & * umask(ji,jj,jk) / fse3u_b(ji,jj,jk) 
    267377                  ENDIF 
    268378                  ! 
    269                   un(ji,jj,jk) = tabres(ji,jj,jk) * umask(ji,jj,jk) 
     379                  un(ji,jj,jk) = tabres(ji,jj,jk) * umask(ji,jj,jk) / fse3u_n(ji,jj,jk) 
    270380               END DO 
    271381            END DO 
     
    304414            DO jj=j1,j2 
    305415               DO ji=i1,i2 
    306                   tabres(ji,jj,jk) = tabres(ji,jj,jk) / e1v(ji,jj) / fse3v_n(ji,jj,jk) 
     416                  tabres(ji,jj,jk) = tabres(ji,jj,jk) / e1v(ji,jj) 
    307417                  ! 
    308418                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 
    309419                     vb(ji,jj,jk) = vb(ji,jj,jk) &  
    310                            & + atfp * ( tabres(ji,jj,jk) - vn(ji,jj,jk) ) * vmask(ji,jj,jk) 
     420                           & + atfp * ( tabres(ji,jj,jk)                     &  
     421                           &              - vn(ji,jj,jk)*fse3v_a(ji,jj,jk) ) &  
     422                           & * vmask(ji,jj,jk) / fse3v_b(ji,jj,jk) 
    311423                  ENDIF 
    312424                  ! 
    313                   vn(ji,jj,jk) = tabres(ji,jj,jk) * vmask(ji,jj,jk) 
     425                  vn(ji,jj,jk) = tabres(ji,jj,jk) * vmask(ji,jj,jk) / fse3v_n(ji,jj,jk) 
    314426               END DO 
    315427            END DO 
     
    345457         DO jj=j1,j2 
    346458            DO ji=i1,i2 
    347                tabres(ji,jj) =  tabres(ji,jj) * hur(ji,jj) / e2u(ji,jj)   
     459               tabres(ji,jj) =  tabres(ji,jj) / e2u(ji,jj)   
    348460               !     
    349461               ! Update "now" 3d velocities: 
     
    352464                  spgu(ji,jj) = spgu(ji,jj) + fse3u_n(ji,jj,jk) * un(ji,jj,jk) 
    353465               END DO 
    354                spgu(ji,jj) = spgu(ji,jj) * hur(ji,jj) 
    355466               ! 
    356                zcorr = tabres(ji,jj) - spgu(ji,jj) 
     467               zcorr = (tabres(ji,jj) - spgu(ji,jj)) * hur(ji,jj) 
    357468               DO jk=1,jpkm1               
    358469                  un(ji,jj,jk) = un(ji,jj,jk) + zcorr * umask(ji,jj,jk)            
     
    360471               ! 
    361472               ! Update barotropic velocities: 
    362 #if defined key_dynspg_ts 
    363                IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 
    364                   zcorr = tabres(ji,jj) - un_b(ji,jj) 
     473#if ! defined key_dynspg_ts 
     474               IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN 
     475#else 
     476               IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN 
     477#endif 
     478                  zcorr = (tabres(ji,jj) - un_b(ji,jj) * hu_a(ji,jj)) * hur_b(ji,jj) 
    365479                  ub_b(ji,jj) = ub_b(ji,jj) + atfp * zcorr * umask(ji,jj,1) 
    366                END IF 
    367 #endif                
    368                un_b(ji,jj) = tabres(ji,jj) * umask(ji,jj,1) 
     480               END IF              
     481               un_b(ji,jj) = tabres(ji,jj) * hur(ji,jj) * umask(ji,jj,1) 
    369482               !        
    370483               ! Correct "before" velocities to hold correct bt component: 
     
    373486                  spgu(ji,jj) = spgu(ji,jj) + fse3u_b(ji,jj,jk) * ub(ji,jj,jk) 
    374487               END DO 
    375                spgu(ji,jj) = spgu(ji,jj) * hur_b(ji,jj) 
    376488               ! 
    377                zcorr = ub_b(ji,jj) - spgu(ji,jj) 
     489               zcorr = ub_b(ji,jj) - spgu(ji,jj) * hur_b(ji,jj) 
    378490               DO jk=1,jpkm1               
    379491                  ub(ji,jj,jk) = ub(ji,jj,jk) + zcorr * umask(ji,jj,jk)            
     
    410522         DO jj=j1,j2 
    411523            DO ji=i1,i2 
    412                tabres(ji,jj) =  tabres(ji,jj) * hvr(ji,jj) / e1v(ji,jj)   
     524               tabres(ji,jj) =  tabres(ji,jj) / e1v(ji,jj)   
    413525               !     
    414526               ! Update "now" 3d velocities: 
     
    417529                  spgv(ji,jj) = spgv(ji,jj) + fse3v_n(ji,jj,jk) * vn(ji,jj,jk) 
    418530               END DO 
    419                spgv(ji,jj) = spgv(ji,jj) * hvr(ji,jj) 
    420531               ! 
    421                zcorr = tabres(ji,jj) - spgv(ji,jj) 
     532               zcorr = (tabres(ji,jj) - spgv(ji,jj)) * hvr(ji,jj) 
    422533               DO jk=1,jpkm1               
    423534                  vn(ji,jj,jk) = vn(ji,jj,jk) + zcorr * vmask(ji,jj,jk)            
     
    425536               ! 
    426537               ! Update barotropic velocities: 
    427 #if defined key_dynspg_ts 
    428                IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 
    429                   zcorr = tabres(ji,jj) - vn_b(ji,jj) 
     538#if ! defined key_dynspg_ts 
     539               IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN 
     540#else 
     541               IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN 
     542#endif 
     543                  zcorr = (tabres(ji,jj) - vn_b(ji,jj)*hv_a(ji,jj)) * hvr_b(ji,jj) 
    430544                  vb_b(ji,jj) = vb_b(ji,jj) + atfp * zcorr * vmask(ji,jj,1) 
    431                END IF 
    432 #endif                
    433                vn_b(ji,jj) = tabres(ji,jj) * vmask(ji,jj,1) 
     545               END IF               
     546               vn_b(ji,jj) = tabres(ji,jj) * hvr(ji,jj) * vmask(ji,jj,1) 
    434547               !        
    435548               ! Correct "before" velocities to hold correct bt component: 
     
    438551                  spgv(ji,jj) = spgv(ji,jj) + fse3v_b(ji,jj,jk) * vb(ji,jj,jk) 
    439552               END DO 
    440                spgv(ji,jj) = spgv(ji,jj) * hvr_b(ji,jj) 
    441553               ! 
    442                zcorr = vb_b(ji,jj) - spgv(ji,jj) 
     554               zcorr = vb_b(ji,jj) - spgv(ji,jj) * hvr_b(ji,jj) 
    443555               DO jk=1,jpkm1               
    444556                  vb(ji,jj,jk) = vb(ji,jj,jk) + zcorr * vmask(ji,jj,jk)            
     
    471583      ELSE 
    472584#if ! defined key_dynspg_ts 
    473          IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN 
     585         IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN 
     586#else 
     587         IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN 
     588#endif 
    474589            DO jj=j1,j2 
    475590               DO ji=i1,i2 
     
    479594            END DO 
    480595         ENDIF 
    481 #endif 
     596         ! 
    482597         DO jj=j1,j2 
    483598            DO ji=i1,i2 
     
    655770# endif /* key_zdftke */  
    656771 
     772   SUBROUTINE updatee3t( ptab, i1, i2, j1, j2, k1, k2, before ) 
     773      !!--------------------------------------------- 
     774      !!           *** ROUTINE updatee3t *** 
     775      !!--------------------------------------------- 
     776      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2 
     777      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab 
     778      LOGICAL, INTENT(in) :: before 
     779      INTEGER :: ji,jj,jk 
     780      REAL(wp) :: zcoef 
     781      !!--------------------------------------------- 
     782      ! 
     783      IF (before) THEN 
     784!> jc tmp: 
     785!         ptab(i1:i2,j1:j2,k1:k2) = fse3t_n(i1:i2,j1:j2,k1:k2) 
     786         ptab(i1:i2,j1:j2,k1:k2) = fse3t_n(i1:i2,j1:j2,k1:k2) / e3t_0(i1:i2,j1:j2,k1:k2) 
     787!< jc tmp: 
     788      ELSE 
     789         ! 
     790         ! 1) Updates at before time step: 
     791         ! ------------------------------- 
     792         ! 
     793!> jc tmp: 
     794         ptab(i1:i2,j1:j2,k1:k2) = ptab(i1:i2,j1:j2,k1:k2) * e3t_0(i1:i2,j1:j2,k1:k2) 
     795!< jc tmp: 
     796 
     797#if ! defined key_dynspg_ts 
     798         IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN 
     799#else 
     800         IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN 
     801#endif 
     802            DO jk = 1, jpkm1 
     803               DO jj=j1,j2 
     804                  DO ji=i1,i2 
     805                     fse3t_b(ji,jj,jk) =   fse3t_b(ji,jj,jk) & 
     806                           & + atfp * ( ptab(ji,jj,jk) - fse3t_n(ji,jj,jk) ) 
     807                  END DO 
     808               END DO 
     809            END DO 
     810            ! 
     811            fse3w_b (i1:i2,j1:j2,1) = e3w_0(i1:i2,j1:j2,1) + fse3t_b(i1:i2,j1:j2,1) - e3t_0(i1:i2,j1:j2,1) 
     812            fsdepw_b(i1:i2,j1:j2,1) = 0.0_wp 
     813            fsdept_b(i1:i2,j1:j2,1) = 0.5_wp * fse3w_b(i1:i2,j1:j2,1) 
     814            ! 
     815            DO jk = 2, jpkm1 
     816               DO jj = j1,j2 
     817                  DO ji = i1,i2             
     818                     zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk)) 
     819                     fse3w_b(ji,jj,jk)  = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) *        &  
     820                     &                                        ( fse3t_b(ji,jj,jk-1) - e3t_0(ji,jj,jk-1) )  & 
     821                     &                                    +            0.5_wp * tmask(ji,jj,jk)   *        & 
     822                     &                                        ( fse3t_b(ji,jj,jk  ) - e3t_0(ji,jj,jk  ) ) 
     823                     fsdepw_b(ji,jj,jk) = fsdepw_b(ji,jj,jk-1) + fse3t_b(ji,jj,jk-1) 
     824                     fsdept_b(ji,jj,jk) =      zcoef  * ( fsdepw_b(ji,jj,jk  ) + 0.5 * fse3w_b(ji,jj,jk))  & 
     825                         &                + (1-zcoef) * ( fsdept_b(ji,jj,jk-1) +       fse3w_b(ji,jj,jk))  
     826                  END DO 
     827               END DO 
     828            END DO 
     829         ENDIF         
     830         ! 
     831         ! 2) Updates at now time step: 
     832         ! ---------------------------- 
     833         ! 
     834         ! Update vertical scale factor at T-points: 
     835         fse3t_n(i1:i2,j1:j2,k1:k2) = ptab(i1:i2,j1:j2,k1:k2) 
     836         ! 
     837         ! Update total depth: 
     838         ht(i1:i2,j1:j2) = 0._wp 
     839         DO jk = 1, jpkm1 
     840            ht(i1:i2,j1:j2) = ht(i1:i2,j1:j2) + fse3t_n(i1:i2,j1:j2,jk) * tmask(i1:i2,j1:j2,jk) 
     841         END DO 
     842         ! 
     843         ! Update vertical scale factor at W-points and depths: 
     844         fse3w_n (i1:i2,j1:j2,1) = e3w_0(i1:i2,j1:j2,1) + fse3t_n(i1:i2,j1:j2,1) - e3t_0(i1:i2,j1:j2,1) 
     845         fsdept_n(i1:i2,j1:j2,1) = 0.5_wp * fse3w_n(i1:i2,j1:j2,1) 
     846         fsdepw_n(i1:i2,j1:j2,1) = 0.0_wp 
     847         fsde3w_n(i1:i2,j1:j2,1) = fsdept_n(i1:i2,j1:j2,1) - (ht(i1:i2,j1:j2)-ht_0(i1:i2,j1:j2)) ! Last term in the rhs is ssh 
     848         ! 
     849         DO jk = 2, jpkm1 
     850            DO jj = j1,j2 
     851               DO ji = i1,i2             
     852               zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk)) 
     853               fse3w_n(ji,jj,jk)  = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) * ( fse3t_n(ji,jj,jk-1) - e3t_0(ji,jj,jk-1) )   & 
     854               &                                    +            0.5_wp * tmask(ji,jj,jk)   * ( fse3t_n(ji,jj,jk  ) - e3t_0(ji,jj,jk  ) ) 
     855               fsdepw_n(ji,jj,jk) = fsdepw_n(ji,jj,jk-1) + fse3t_n(ji,jj,jk-1) 
     856               fsdept_n(ji,jj,jk) =      zcoef  * ( fsdepw_n(ji,jj,jk  ) + 0.5 * fse3w_n(ji,jj,jk))  & 
     857                   &                + (1-zcoef) * ( fsdept_n(ji,jj,jk-1) +       fse3w_n(ji,jj,jk))  
     858               fsde3w_n(ji,jj,jk) = fsdept_n(ji,jj,jk) - (ht(ji,jj)-ht_0(ji,jj)) ! Last term in the rhs is ssh 
     859               END DO 
     860            END DO 
     861         END DO 
     862         ! 
     863      ENDIF 
     864      ! 
     865   END SUBROUTINE updatee3t 
     866 
    657867#else 
    658868CONTAINS 
Note: See TracChangeset for help on using the changeset viewer.