Changeset 7971


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

Add zstar coordinate with AGRIF

Location:
branches/2017/dev_r7963_nemo_v3_6_AGRIF-3_AGRIFVVL/NEMOGCM/NEMO
Files:
5 edited

Legend:

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

    r6204 r7971  
    137137         DO jk=1,jpkm1 
    138138            DO jj=1,jpj 
    139                spgu(2,jj)=spgu(2,jj)+fse3u(2,jj,jk)*ua(2,jj,jk) 
     139               spgu(2,jj)=spgu(2,jj)+fse3u_a(2,jj,jk)*ua(2,jj,jk) 
    140140            END DO 
    141141         END DO 
     
    143143         DO jj=1,jpj 
    144144            IF (umask(2,jj,1).NE.0.) THEN 
    145                spgu(2,jj)=spgu(2,jj)/hu(2,jj) 
     145               spgu(2,jj)=spgu(2,jj)*hur_a(2,jj) 
    146146            ENDIF 
    147147         END DO 
     
    161161         DO jk=1,jpkm1 
    162162            DO jj=1,jpj 
    163                spgu1(2,jj)=spgu1(2,jj)+fse3u(2,jj,jk)*ua(2,jj,jk) 
     163               spgu1(2,jj)=spgu1(2,jj)+fse3u_a(2,jj,jk)*ua(2,jj,jk) 
    164164            END DO 
    165165         END DO 
     
    167167         DO jj=1,jpj 
    168168            IF (umask(2,jj,1).NE.0.) THEN 
    169                spgu1(2,jj)=spgu1(2,jj)/hu(2,jj) 
     169               spgu1(2,jj)=spgu1(2,jj)*hur_a(2,jj) 
    170170            ENDIF 
    171171         END DO 
     
    207207         DO jk=1,jpkm1 
    208208            DO jj=1,jpj 
    209                spgu(nlci-2,jj)=spgu(nlci-2,jj)+fse3u(nlci-2,jj,jk)*ua(nlci-2,jj,jk) 
     209               spgu(nlci-2,jj)=spgu(nlci-2,jj)+fse3u_a(nlci-2,jj,jk)*ua(nlci-2,jj,jk) 
    210210            ENDDO 
    211211         ENDDO 
    212212         DO jj=1,jpj 
    213213            IF (umask(nlci-2,jj,1).NE.0.) THEN 
    214                spgu(nlci-2,jj)=spgu(nlci-2,jj)/hu(nlci-2,jj) 
     214               spgu(nlci-2,jj)=spgu(nlci-2,jj)*hur_a(nlci-2,jj) 
    215215            ENDIF 
    216216         END DO 
     
    229229         DO jk=1,jpkm1 
    230230            DO jj=1,jpj 
    231                spgu1(nlci-2,jj)=spgu1(nlci-2,jj)+fse3u(nlci-2,jj,jk)*ua(nlci-2,jj,jk)*umask(nlci-2,jj,jk) 
     231               spgu1(nlci-2,jj)=spgu1(nlci-2,jj)+fse3u_a(nlci-2,jj,jk)*ua(nlci-2,jj,jk)*umask(nlci-2,jj,jk) 
    232232            END DO 
    233233         END DO 
    234234         DO jj=1,jpj 
    235235            IF (umask(nlci-2,jj,1).NE.0.) THEN 
    236                spgu1(nlci-2,jj)=spgu1(nlci-2,jj)/hu(nlci-2,jj) 
     236               spgu1(nlci-2,jj)=spgu1(nlci-2,jj)*hur_a(nlci-2,jj) 
    237237            ENDIF 
    238238         END DO 
     
    278278         DO jk=1,jpkm1 
    279279            DO ji=1,jpi 
    280                spgv(ji,2)=spgv(ji,2)+fse3v(ji,2,jk)*va(ji,2,jk) 
     280               spgv(ji,2)=spgv(ji,2)+fse3v_a(ji,2,jk)*va(ji,2,jk) 
    281281            END DO 
    282282         END DO 
     
    284284         DO ji=1,jpi 
    285285            IF (vmask(ji,2,1).NE.0.) THEN 
    286                spgv(ji,2)=spgv(ji,2)/hv(ji,2) 
     286               spgv(ji,2)=spgv(ji,2)*hvr_a(ji,2) 
    287287            ENDIF 
    288288         END DO 
     
    302302         DO jk=1,jpkm1 
    303303            DO ji=1,jpi 
    304                spgv1(ji,2)=spgv1(ji,2)+fse3v(ji,2,jk)*va(ji,2,jk)*vmask(ji,2,jk) 
     304               spgv1(ji,2)=spgv1(ji,2)+fse3v_a(ji,2,jk)*va(ji,2,jk)*vmask(ji,2,jk) 
    305305            END DO 
    306306         END DO 
     
    308308         DO ji=1,jpi 
    309309            IF (vmask(ji,2,1).NE.0.) THEN 
    310                spgv1(ji,2)=spgv1(ji,2)/hv(ji,2) 
     310               spgv1(ji,2)=spgv1(ji,2)*hvr_a(ji,2) 
    311311            ENDIF 
    312312         END DO 
     
    353353         DO jk=1,jpkm1 
    354354            DO ji=1,jpi 
    355                spgv(ji,nlcj-2)=spgv(ji,nlcj-2)+fse3v(ji,nlcj-2,jk)*va(ji,nlcj-2,jk) 
     355               spgv(ji,nlcj-2)=spgv(ji,nlcj-2)+fse3v_a(ji,nlcj-2,jk)*va(ji,nlcj-2,jk) 
    356356            END DO 
    357357         END DO 
     
    359359         DO ji=1,jpi 
    360360            IF (vmask(ji,nlcj-2,1).NE.0.) THEN 
    361                spgv(ji,nlcj-2)=spgv(ji,nlcj-2)/hv(ji,nlcj-2) 
     361               spgv(ji,nlcj-2)=spgv(ji,nlcj-2)*hvr_a(ji,nlcj-2) 
    362362            ENDIF 
    363363         END DO 
     
    378378         DO jk=1,jpkm1 
    379379            DO ji=1,jpi 
    380                spgv1(ji,nlcj-2)=spgv1(ji,nlcj-2)+fse3v(ji,nlcj-2,jk)*va(ji,nlcj-2,jk) 
     380               spgv1(ji,nlcj-2)=spgv1(ji,nlcj-2)+fse3v_a(ji,nlcj-2,jk)*va(ji,nlcj-2,jk) 
    381381            END DO 
    382382         END DO 
     
    384384         DO ji=1,jpi 
    385385            IF (vmask(ji,nlcj-2,1).NE.0.) THEN 
    386                spgv1(ji,nlcj-2)=spgv1(ji,nlcj-2)/hv(ji,nlcj-2) 
     386               spgv1(ji,nlcj-2)=spgv1(ji,nlcj-2)*hvr_a(ji,nlcj-2) 
    387387            ENDIF 
    388388         END DO 
     
    503503         zt = REAL(Agrif_NbStepint(),wp) / zrhot 
    504504      ENDIF 
    505  
    506       ! Linear interpolation of sea level 
    507       Agrif_SpecialValue    = 0.e0 
    508       Agrif_UseSpecialValue = .TRUE. 
    509       CALL Agrif_Bc_variable(sshn_id,calledweight=zt, procname=interpsshn ) 
    510       Agrif_UseSpecialValue = .FALSE. 
    511505 
    512506      ! Interpolate barotropic fluxes 
     
    539533   SUBROUTINE Agrif_ssh( kt ) 
    540534      !!---------------------------------------------------------------------- 
    541       !!                  ***  ROUTINE Agrif_DYN  *** 
     535      !!                  ***  ROUTINE Agrif_ssh  *** 
    542536      !!----------------------------------------------------------------------   
    543537      INTEGER, INTENT(in) ::   kt 
    544538      !! 
     539      INTEGER :: ji, jj 
    545540      !!----------------------------------------------------------------------   
    546541 
    547542      IF( Agrif_Root() )   RETURN 
    548543 
     544      ! Linear interpolation of sea level 
     545      Agrif_SpecialValue    = 0.e0 
     546      Agrif_UseSpecialValue = .TRUE. 
     547      CALL Agrif_Bc_variable(sshn_id, procname=interpsshn ) 
     548      Agrif_UseSpecialValue = .FALSE. 
     549 
    549550      IF((nbondi == -1).OR.(nbondi == 2)) THEN 
    550          ssha(2,:)=ssha(3,:) 
    551          sshn(2,:)=sshn(3,:) 
     551         DO jj=1,jpj 
     552            ssha(2,jj) = hbdy_w(jj) 
     553         END DO 
    552554      ENDIF 
    553555 
    554556      IF((nbondi == 1).OR.(nbondi == 2)) THEN 
    555          ssha(nlci-1,:)=ssha(nlci-2,:) 
    556          sshn(nlci-1,:)=sshn(nlci-2,:) 
     557         DO jj=1,jpj 
     558            ssha(nlci-1,jj) = hbdy_e(jj) 
     559         END DO 
    557560      ENDIF 
    558561 
    559562      IF((nbondj == -1).OR.(nbondj == 2)) THEN 
    560          ssha(:,2)=ssha(:,3) 
    561          sshn(:,2)=sshn(:,3) 
     563         DO ji=1,jpi 
     564            ssha(ji,2) = hbdy_s(ji) 
     565         END DO 
    562566      ENDIF 
    563567 
    564568      IF((nbondj == 1).OR.(nbondj == 2)) THEN 
    565          ssha(:,nlcj-1)=ssha(:,nlcj-2) 
    566          sshn(:,nlcj-1)=sshn(:,nlcj-2) 
     569         DO ji=1,jpi 
     570            ssha(ji,nlcj-1) = hbdy_n(ji) 
     571         END DO 
    567572      ENDIF 
    568573 
     
    812817               DO ji=i1,i2 
    813818                  ptab(ji,jj,jk) = e2u(ji,jj) * un(ji,jj,jk) 
    814                   ptab(ji,jj,jk) = ptab(ji,jj,jk) * fse3u(ji,jj,jk) 
     819                  ptab(ji,jj,jk) = ptab(ji,jj,jk) * fse3u_n(ji,jj,jk) 
    815820               END DO 
    816821            END DO 
     
    821826            DO jj=j1,j2 
    822827               ua(i1:i2,jj,jk) = (ptab(i1:i2,jj,jk)/(zrhoy*e2u(i1:i2,jj))) 
    823                ua(i1:i2,jj,jk) = ua(i1:i2,jj,jk) / fse3u(i1:i2,jj,jk) 
     828               ua(i1:i2,jj,jk) = ua(i1:i2,jj,jk) / fse3u_a(i1:i2,jj,jk) 
    824829            END DO 
    825830         END DO 
     
    880885               DO ji=i1,i2 
    881886                  ptab(ji,jj,jk) = e1v(ji,jj) * vn(ji,jj,jk) 
    882                   ptab(ji,jj,jk) = ptab(ji,jj,jk) * fse3v(ji,jj,jk) 
     887                  ptab(ji,jj,jk) = ptab(ji,jj,jk) * fse3v_n(ji,jj,jk) 
    883888               END DO 
    884889            END DO 
     
    889894            DO jj=j1,j2 
    890895               va(i1:i2,jj,jk) = (ptab(i1:i2,jj,jk)/(zrhox*e1v(i1:i2,jj))) 
    891                va(i1:i2,jj,jk) = va(i1:i2,jj,jk) / fse3v(i1:i2,jj,jk) 
     896               va(i1:i2,jj,jk) = va(i1:i2,jj,jk) / fse3v_a(i1:i2,jj,jk) 
    892897            END DO 
    893898         END DO 
     
    11101115         ! Polynomial interpolation coefficients: 
    11111116         zat = zrhot * (  zt1**2._wp * (-2._wp*zt1 + 3._wp)        & 
    1112                &      - zt0**2._wp * (-2._wp*zt0 + 3._wp)        )  
     1117                 &      - zt0**2._wp * (-2._wp*zt0 + 3._wp)        )  
    11131118         !  
    11141119         IF(western_side ) ubdy_w(j1:j2) = zat * ptab(i1,j1:j2)   
     
    11511156         ! Polynomial interpolation coefficients: 
    11521157         zat = zrhot * (  zt1**2._wp * (-2._wp*zt1 + 3._wp)        & 
    1153                &      - zt0**2._wp * (-2._wp*zt0 + 3._wp)        )  
     1158                 &      - zt0**2._wp * (-2._wp*zt0 + 3._wp)        )  
    11541159         ! 
    11551160         IF(western_side ) vbdy_w(j1:j2) = zat * ptab(i1,j1:j2)   
  • 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 
  • branches/2017/dev_r7963_nemo_v3_6_AGRIF-3_AGRIFVVL/NEMOGCM/NEMO/NST_SRC/agrif_user.F90

    r6204 r7971  
    309309!     or the absolute maximum nesting level...TBC                         
    310310   IF ( Agrif_Level().EQ.Agrif_MaxLevel() ) THEN  
    311       ! NB: Do tracers first, dynamics after because nbcline incremented in dynamics 
     311      ! NB: Order matters below:  
     312      CALL Agrif_Update_vvl()  
    312313      CALL Agrif_Update_tra() 
    313314      CALL Agrif_Update_dyn() 
     
    445446   CALL Agrif_Set_Updatetype(avm_id, update = AGRIF_Update_Average) 
    446447# endif 
     448 
     449   CALL Agrif_Set_Updatetype(e3t_id, update = AGRIF_Update_Average) 
    447450 
    448451! High order updates 
  • branches/2017/dev_r7963_nemo_v3_6_AGRIF-3_AGRIFVVL/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90

    r6348 r7971  
    10061006 
    10071007#if defined key_agrif 
    1008       IF (.NOT.Agrif_Root()) CALL ctl_stop( 'AGRIF not implemented with non-linear free surface (key_vvl)' ) 
     1008      IF ((.NOT.Agrif_Root()).AND.(.NOT.ln_vvl_zstar)) CALL ctl_stop( 'AGRIF implemented with zstar coordinate only (key_vvl)' ) 
    10091009#endif 
    10101010 
  • branches/2017/dev_r7963_nemo_v3_6_AGRIF-3_AGRIFVVL/NEMOGCM/NEMO/OPA_SRC/step.F90

    r7494 r7971  
    350350 
    351351      IF ( Agrif_NbStepint().EQ.0 ) THEN 
     352         IF( lk_vvl )          CALL Agrif_Update_vvl()      ! Update vert. grid  
    352353                               CALL Agrif_Update_Tra()      ! Update active tracers 
    353354                               CALL Agrif_Update_Dyn()      ! Update momentum 
Note: See TracChangeset for help on using the changeset viewer.