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 8993 for branches/2017/dev_METO_MERCATOR_2017/NEMOGCM/NEMO/NST_SRC – NEMO

Ignore:
Timestamp:
2017-12-12T16:42:29+01:00 (6 years ago)
Author:
timgraham
Message:

Merged Mercator branch in

Location:
branches/2017/dev_METO_MERCATOR_2017/NEMOGCM/NEMO/NST_SRC
Files:
9 edited

Legend:

Unmodified
Added
Removed
  • branches/2017/dev_METO_MERCATOR_2017/NEMOGCM/NEMO/NST_SRC/agrif_lim2_update.F90

    r5656 r8993  
    5959      Agrif_SpecialValueFineGrid = 0. 
    6060# if defined TWO_WAY 
    61       IF( MOD(nbcline,nbclineupdate) == 0) THEN 
    62          CALL Agrif_Update_Variable( adv_ice_id , procname = update_adv_ice  ) 
    63          CALL Agrif_Update_Variable( u_ice_id   , procname = update_u_ice    ) 
    64          CALL Agrif_Update_Variable( v_ice_id   , procname = update_v_ice    ) 
    65       ELSE 
    66          CALL Agrif_Update_Variable( adv_ice_id , locupdate=(/0,2/), procname = update_adv_ice  ) 
    67          CALL Agrif_Update_Variable( u_ice_id   , locupdate=(/0,1/), procname = update_u_ice    ) 
    68          CALL Agrif_Update_Variable( v_ice_id   , locupdate=(/0,1/), procname = update_v_ice    ) 
    69       ENDIF 
     61      CALL Agrif_Update_Variable( adv_ice_id , procname = update_adv_ice  ) 
     62      CALL Agrif_Update_Variable( u_ice_id   , procname = update_u_ice    ) 
     63      CALL Agrif_Update_Variable( v_ice_id   , procname = update_v_ice    ) 
     64!      CALL Agrif_Update_Variable( adv_ice_id , locupdate=(/0,2/), procname = update_adv_ice  ) 
     65!      CALL Agrif_Update_Variable( u_ice_id   , locupdate=(/0,1/), procname = update_u_ice    ) 
     66!      CALL Agrif_Update_Variable( v_ice_id   , locupdate=(/0,1/), procname = update_v_ice    ) 
    7067# endif 
    7168      ! 
  • branches/2017/dev_METO_MERCATOR_2017/NEMOGCM/NEMO/NST_SRC/agrif_lim3_update.F90

    r7761 r8993  
    5959      Agrif_SpecialValueFineGrid = -9999. 
    6060# if defined TWO_WAY 
    61       IF( MOD(nbcline,nbclineupdate) == 0) THEN ! update the whole basin at each nbclineupdate (=nn_cln_update) baroclinic parent time steps 
    62                                                 ! nbcline is incremented (+1) at the end of each parent time step from 0 (1st time step) 
    63          CALL Agrif_Update_Variable( tra_ice_id , procname = update_tra_ice  ) 
    64          CALL Agrif_Update_Variable( u_ice_id   , procname = update_u_ice    ) 
    65          CALL Agrif_Update_Variable( v_ice_id   , procname = update_v_ice    ) 
    66       ELSE                                      ! update only the boundaries defined par locupdate 
    67          CALL Agrif_Update_Variable( tra_ice_id , locupdate=(/0,2/), procname = update_tra_ice  ) 
    68          CALL Agrif_Update_Variable( u_ice_id   , locupdate=(/0,1/), procname = update_u_ice    ) 
    69          CALL Agrif_Update_Variable( v_ice_id   , locupdate=(/0,1/), procname = update_v_ice    ) 
    70       ENDIF 
     61      CALL Agrif_Update_Variable( tra_ice_id , procname = update_tra_ice  ) 
     62      CALL Agrif_Update_Variable( u_ice_id   , procname = update_u_ice    ) 
     63      CALL Agrif_Update_Variable( v_ice_id   , procname = update_v_ice    ) 
     64 
     65!      CALL Agrif_Update_Variable( tra_ice_id , locupdate=(/0,2/), procname = update_tra_ice  ) 
     66!      CALL Agrif_Update_Variable( u_ice_id   , locupdate=(/0,1/), procname = update_u_ice    ) 
     67!      CALL Agrif_Update_Variable( v_ice_id   , locupdate=(/0,1/), procname = update_v_ice    ) 
    7168# endif 
    7269      Agrif_UseSpecialValueInUpdate = .FALSE. 
  • branches/2017/dev_METO_MERCATOR_2017/NEMOGCM/NEMO/NST_SRC/agrif_oce.F90

    r5656 r8993  
    2020   !                                              !!* Namelist namagrif: AGRIF parameters 
    2121   LOGICAL , PUBLIC ::   ln_spc_dyn    = .FALSE.   !: 
    22    INTEGER , PUBLIC ::   nn_cln_update = 3         !: update frequency  
    2322   INTEGER , PUBLIC, PARAMETER ::   nn_sponge_len = 2  !: Sponge width (in number of parent grid points) 
    2423   REAL(wp), PUBLIC ::   rn_sponge_tra = 2800.     !: sponge coeff. for tracers 
    2524   REAL(wp), PUBLIC ::   rn_sponge_dyn = 2800.     !: sponge coeff. for dynamics 
    2625   LOGICAL , PUBLIC ::   ln_chk_bathy  = .FALSE.   !: check of parent bathymetry  
     26   LOGICAL , PUBLIC ::   lk_agrif_clp  = .FALSE.   !: Flag to retrieve clamped open boundaries 
    2727 
    2828   !                                              !!! OLD namelist names 
    29    INTEGER , PUBLIC ::   nbcline = 0               !: update counter 
    30    INTEGER , PUBLIC ::   nbclineupdate             !: update frequency  
    3129   REAL(wp), PUBLIC ::   visc_tra                  !: sponge coeff. for tracers 
    3230   REAL(wp), PUBLIC ::   visc_dyn                  !: sponge coeff. for dynamics 
     
    3533   LOGICAL , PUBLIC :: spongedoneU = .FALSE.       !: dynamics sponge layer indicator 
    3634   LOGICAL , PUBLIC :: lk_agrif_fstep = .TRUE.     !: if true: first step 
    37    LOGICAL , PUBLIC :: lk_agrif_doupd = .TRUE.     !: if true: send update from current grid 
    3835   LOGICAL , PUBLIC :: lk_agrif_debug = .FALSE.    !: if true: print debugging info 
    3936 
     
    6562   INTEGER :: e3t_id, e1u_id, e2v_id, sshn_id 
    6663   INTEGER :: scales_t_id 
    67 # if defined key_zdftke 
     64# if defined key_zdftke || defined key_zdfgls 
    6865   INTEGER :: avt_id, avm_id, en_id 
    6966# endif   
  • branches/2017/dev_METO_MERCATOR_2017/NEMOGCM/NEMO/NST_SRC/agrif_opa_interp.F90

    r7646 r8993  
    2424   USE agrif_oce 
    2525   USE phycst 
     26   USE dynspg_ts, ONLY: un_adv, vn_adv 
    2627   ! 
    2728   USE in_out_manager 
     
    3839   PUBLIC   interpunb, interpvnb, interpub2b, interpvb2b 
    3940   PUBLIC   interpe3t, interpumsk, interpvmsk 
    40 # if defined key_zdftke 
    41    PUBLIC   Agrif_tke, interpavm 
     41# if defined key_zdftke || defined key_zdfgls 
     42   PUBLIC   Agrif_avm, interpavm 
    4243# endif 
    4344 
     
    116117         ENDIF 
    117118         ! 
    118          DO jk=1,jpkm1                 ! Smooth 
    119             DO jj=j1,j2 
    120                ua(2,jj,jk) = 0.25_wp*(ua(1,jj,jk)+2._wp*ua(2,jj,jk)+ua(3,jj,jk)) 
    121                ua(2,jj,jk) = ua(2,jj,jk) * umask(2,jj,jk) 
    122             END DO 
    123          END DO 
     119         IF (.NOT.lk_agrif_clp) THEN 
     120            DO jk=1,jpkm1              ! Smooth 
     121               DO jj=j1,j2 
     122                  ua(2,jj,jk) = 0.25_wp*(ua(1,jj,jk)+2._wp*ua(2,jj,jk)+ua(3,jj,jk)) 
     123                  ua(2,jj,jk) = ua(2,jj,jk) * umask(2,jj,jk) 
     124               END DO 
     125            END DO 
     126         END IF 
    124127         ! 
    125128         zub(2,:) = 0._wp              ! Correct transport 
     
    185188         ENDIF 
    186189 
    187          DO jk = 1, jpkm1              ! Smooth 
    188             DO jj = j1, j2 
    189                ua(nlci-2,jj,jk) = 0.25_wp * umask(nlci-2,jj,jk)      & 
    190                   &             * ( ua(nlci-3,jj,jk) + 2._wp*ua(nlci-2,jj,jk) + ua(nlci-1,jj,jk) ) 
    191             END DO 
    192          END DO 
     190         IF (.NOT.lk_agrif_clp) THEN 
     191            DO jk = 1, jpkm1           ! Smooth 
     192               DO jj = j1, j2 
     193                  ua(nlci-2,jj,jk) = 0.25_wp * umask(nlci-2,jj,jk)      & 
     194                     &             * ( ua(nlci-3,jj,jk) + 2._wp*ua(nlci-2,jj,jk) + ua(nlci-1,jj,jk) ) 
     195               END DO 
     196            END DO 
     197         ENDIF 
    193198 
    194199         zub(nlci-2,:) = 0._wp        ! Correct transport 
     
    254259         ENDIF 
    255260         ! 
    256          DO jk = 1, jpkm1              ! Smooth 
    257             DO ji = i1, i2 
    258                va(ji,2,jk) = 0.25_wp * vmask(ji,2,jk)    & 
    259                   &        * ( va(ji,1,jk) + 2._wp*va(ji,2,jk) + va(ji,3,jk) ) 
    260             END DO 
    261          END DO 
     261         IF (.NOT.lk_agrif_clp) THEN 
     262            DO jk = 1, jpkm1              ! Smooth 
     263               DO ji = i1, i2 
     264                  va(ji,2,jk) = 0.25_wp * vmask(ji,2,jk)    & 
     265                     &        * ( va(ji,1,jk) + 2._wp*va(ji,2,jk) + va(ji,3,jk) ) 
     266               END DO 
     267            END DO 
     268         ENDIF 
    262269         ! 
    263270         zvb(:,2) = 0._wp              ! Correct transport 
     
    323330         ENDIF 
    324331         ! 
    325          DO jk = 1, jpkm1              ! Smooth 
    326             DO ji = i1, i2 
    327                va(ji,nlcj-2,jk) = 0.25_wp * vmask(ji,nlcj-2,jk)   & 
    328                   &             * ( va(ji,nlcj-3,jk) + 2._wp * va(ji,nlcj-2,jk) + va(ji,nlcj-1,jk) ) 
    329             END DO 
    330          END DO 
     332         IF (.NOT.lk_agrif_clp) THEN 
     333            DO jk = 1, jpkm1           ! Smooth 
     334               DO ji = i1, i2 
     335                  va(ji,nlcj-2,jk) = 0.25_wp * vmask(ji,nlcj-2,jk)   & 
     336                     &             * ( va(ji,nlcj-3,jk) + 2._wp * va(ji,nlcj-2,jk) + va(ji,nlcj-1,jk) ) 
     337               END DO 
     338            END DO 
     339         ENDIF 
    331340         ! 
    332341         zvb(:,nlcj-2) = 0._wp         ! Correct transport 
     
    449458      INTEGER :: ji, jj 
    450459      LOGICAL :: ll_int_cons 
    451       REAL(wp) :: zrhot, zt 
    452460      !!----------------------------------------------------------------------   
    453461      ! 
     
    456464      ll_int_cons = ln_bt_fw ! Assume conservative temporal integration in the forward case only 
    457465      ! 
    458       zrhot = Agrif_rhot() 
    459       ! 
    460       ! "Central" time index for interpolation: 
    461       IF( ln_bt_fw ) THEN 
    462          zt = REAL( Agrif_NbStepint()+0.5_wp, wp ) / zrhot 
    463       ELSE 
    464          zt = REAL( Agrif_NbStepint()       , wp ) / zrhot 
    465       ENDIF 
    466       ! 
    467       ! Linear interpolation of sea level 
    468       Agrif_SpecialValue    = 0._wp 
    469       Agrif_UseSpecialValue = .TRUE. 
    470       CALL Agrif_Bc_variable( sshn_id, calledweight=zt, procname=interpsshn ) 
    471       Agrif_UseSpecialValue = .FALSE. 
     466      ! Enforce volume conservation if no time refinement:   
     467      IF ( Agrif_rhot()==1 ) ll_int_cons=.TRUE.   
    472468      ! 
    473469      ! Interpolate barotropic fluxes 
    474       Agrif_SpecialValue=0. 
     470      Agrif_SpecialValue=0._wp 
    475471      Agrif_UseSpecialValue = ln_spc_dyn 
    476472      ! 
     
    491487         ubdy_n(:) = 0._wp   ;   vbdy_n(:) = 0._wp  
    492488         ubdy_s(:) = 0._wp   ;   vbdy_s(:) = 0._wp 
    493          CALL Agrif_Bc_variable( unb_id, calledweight=zt, procname=interpunb ) 
    494          CALL Agrif_Bc_variable( vnb_id, calledweight=zt, procname=interpvnb ) 
     489         CALL Agrif_Bc_variable( unb_id, procname=interpunb ) 
     490         CALL Agrif_Bc_variable( vnb_id, procname=interpvnb ) 
    495491      ENDIF 
    496492      Agrif_UseSpecialValue = .FALSE. 
     
    501497   SUBROUTINE Agrif_ssh( kt ) 
    502498      !!---------------------------------------------------------------------- 
    503       !!                  ***  ROUTINE Agrif_DYN  *** 
     499      !!                  ***  ROUTINE Agrif_ssh  *** 
    504500      !!----------------------------------------------------------------------   
    505501      INTEGER, INTENT(in) ::   kt 
    506502      !! 
     503      INTEGER :: ji, jj 
    507504      !!----------------------------------------------------------------------   
    508505      ! 
    509506      IF( Agrif_Root() )   RETURN 
     507      !       
     508      ! Linear interpolation in time of sea level 
     509      ! 
     510      Agrif_SpecialValue    = 0._wp 
     511      Agrif_UseSpecialValue = .TRUE. 
     512      CALL Agrif_Bc_variable(sshn_id, procname=interpsshn ) 
     513      Agrif_UseSpecialValue = .FALSE. 
    510514      ! 
    511515      IF((nbondi == -1).OR.(nbondi == 2)) THEN 
    512          ssha(2,:)=ssha(3,:) 
    513          sshn(2,:)=sshn(3,:) 
     516         DO jj=1,jpj 
     517            ssha(2,jj) = hbdy_w(jj) 
     518         END DO 
    514519      ENDIF 
    515520      ! 
    516521      IF((nbondi == 1).OR.(nbondi == 2)) THEN 
    517          ssha(nlci-1,:)=ssha(nlci-2,:) 
    518          sshn(nlci-1,:)=sshn(nlci-2,:) 
     522         DO jj=1,jpj 
     523            ssha(nlci-1,jj) = hbdy_e(jj) 
     524         END DO 
    519525      ENDIF 
    520526      ! 
    521527      IF((nbondj == -1).OR.(nbondj == 2)) THEN 
    522          ssha(:,2)=ssha(:,3) 
    523          sshn(:,2)=sshn(:,3) 
     528         DO ji=1,jpi 
     529            ssha(ji,2) = hbdy_s(ji) 
     530         END DO 
    524531      ENDIF 
    525532      ! 
    526533      IF((nbondj == 1).OR.(nbondj == 2)) THEN 
    527          ssha(:,nlcj-1)=ssha(:,nlcj-2) 
    528          sshn(:,nlcj-1)=sshn(:,nlcj-2) 
     534         DO ji=1,jpi 
     535            ssha(ji,nlcj-1) = hbdy_n(ji) 
     536         END DO 
    529537      ENDIF 
    530538      ! 
     
    541549      !!----------------------------------------------------------------------   
    542550      ! 
     551      ! 
     552      IF( Agrif_Root() )   RETURN 
     553      ! 
    543554      IF((nbondi == -1).OR.(nbondi == 2)) THEN 
    544555         DO jj = 1, jpj 
     
    567578   END SUBROUTINE Agrif_ssh_ts 
    568579 
    569 # if defined key_zdftke 
    570  
    571    SUBROUTINE Agrif_tke 
    572       !!---------------------------------------------------------------------- 
    573       !!                  ***  ROUTINE Agrif_tke  *** 
     580# if defined key_zdftke || defined key_zdfgls 
     581 
     582   SUBROUTINE Agrif_avm 
     583      !!---------------------------------------------------------------------- 
     584      !!                  ***  ROUTINE Agrif_avm  *** 
    574585      !!----------------------------------------------------------------------   
    575586      REAL(wp) ::   zalpha 
    576587      !!----------------------------------------------------------------------   
    577588      ! 
    578       zalpha = REAL( Agrif_NbStepint() + Agrif_IRhot() - 1, wp ) / REAL( Agrif_IRhot(), wp ) 
    579       IF( zalpha > 1. )   zalpha = 1. 
     589      IF( Agrif_Root() )   RETURN 
     590      ! 
     591!      zalpha = REAL( Agrif_NbStepint() + Agrif_IRhot() - 1, wp ) / REAL( Agrif_IRhot(), wp ) 
     592!      IF( zalpha > 1. )   zalpha = 1. 
     593      zalpha = 1._wp ! JC: proper time interpolation impossible   
     594                     ! => use last available value from parent  
    580595      ! 
    581596      Agrif_SpecialValue    = 0.e0 
     
    586601      Agrif_UseSpecialValue = .FALSE. 
    587602      ! 
    588    END SUBROUTINE Agrif_tke 
     603   END SUBROUTINE Agrif_avm 
    589604    
    590605# endif 
     
    609624         ptab(i1:i2,j1:j2,k1:k2,n1:n2) = tsn(i1:i2,j1:j2,k1:k2,n1:n2) 
    610625      ELSE 
     626         IF (lk_agrif_clp) THEN 
     627            DO jn = 1, jpts 
     628               DO jk = 1, jpkm1 
     629                  DO ji = i1,i2 
     630                     DO jj = j1,j2 
     631                        tsa(ji,jj,jk,jn) = ptab(ji,jj,jk,jn) 
     632                     END DO 
     633                  END DO 
     634               END DO 
     635            END DO            
     636            return 
     637         ENDIF 
    611638         ! 
    612639         western_side  = (nb == 1).AND.(ndir == 1) 
     
    781808      ! 
    782809      IF( before ) THEN  
    783          DO jk = k1, jpk 
     810         DO jk = 1, jpkm1 
    784811            ptab(i1:i2,j1:j2,jk) = e2u(i1:i2,j1:j2) * e3u_n(i1:i2,j1:j2,jk) * un(i1:i2,j1:j2,jk) 
    785812         END DO 
     
    788815         DO jk = 1, jpkm1 
    789816            DO jj=j1,j2 
    790                ua(i1:i2,jj,jk) = ptab(i1:i2,jj,jk) / ( zrhoy * e2u(i1:i2,jj) * e3u_n(i1:i2,jj,jk) ) 
     817               ua(i1:i2,jj,jk) = ptab(i1:i2,jj,jk) / ( zrhoy * e2u(i1:i2,jj) * e3u_a(i1:i2,jj,jk) ) 
    791818            END DO 
    792819         END DO 
     
    808835      !!---------------------------------------------------------------------- 
    809836      !       
    810       IF( before ) THEN       !interpv entre 1 et k2 et interpv2d en jpkp1 
    811          DO jk = k1, jpk 
     837      IF( before ) THEN   
     838         DO jk = 1, jpkm1 
    812839            ptab(i1:i2,j1:j2,jk) = e1v(i1:i2,j1:j2) * e3v_n(i1:i2,j1:j2,jk) * vn(i1:i2,j1:j2,jk) 
    813840         END DO 
     
    815842         zrhox= Agrif_Rhox() 
    816843         DO jk = 1, jpkm1 
    817             va(i1:i2,j1:j2,jk) = ptab(i1:i2,j1:j2,jk) / ( zrhox * e1v(i1:i2,j1:j2) * e3v_n(i1:i2,j1:j2,jk) ) 
     844            va(i1:i2,j1:j2,jk) = ptab(i1:i2,j1:j2,jk) / ( zrhox * e1v(i1:i2,j1:j2) * e3v_a(i1:i2,j1:j2,jk) ) 
    818845         END DO 
    819846      ENDIF 
     
    9781005      !!----------------------------------------------------------------------   
    9791006      IF( before ) THEN 
    980          ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * ub2_b(i1:i2,j1:j2) 
     1007         IF ( ln_bt_fw ) THEN 
     1008            ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * ub2_b(i1:i2,j1:j2) 
     1009         ELSE 
     1010            ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * un_adv(i1:i2,j1:j2) 
     1011         ENDIF 
    9811012      ELSE 
    9821013         western_side  = (nb == 1).AND.(ndir == 1) 
     
    10161047      ! 
    10171048      IF( before ) THEN 
    1018          ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vb2_b(i1:i2,j1:j2) 
     1049         IF ( ln_bt_fw ) THEN 
     1050            ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vb2_b(i1:i2,j1:j2) 
     1051         ELSE 
     1052            ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vn_adv(i1:i2,j1:j2) 
     1053         ENDIF 
    10191054      ELSE       
    10201055         western_side  = (nb == 1).AND.(ndir == 1) 
     
    11751210   END SUBROUTINE interpvmsk 
    11761211 
    1177 # if defined key_zdftke 
     1212# if defined key_zdftke || defined key_zdfgls 
    11781213 
    11791214   SUBROUTINE interpavm( ptab, i1, i2, j1, j2, k1, k2, before ) 
     
    11891224         ptab (i1:i2,j1:j2,k1:k2) = avm_k(i1:i2,j1:j2,k1:k2) 
    11901225      ELSE 
    1191          avm_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) 
     1226         avm  (i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) 
    11921227      ENDIF 
    11931228      ! 
    11941229   END SUBROUTINE interpavm 
    11951230 
    1196 # endif /* key_zdftke */ 
     1231# endif /* key_zdftke || key_zdfgls */ 
    11971232 
    11981233#else 
  • branches/2017/dev_METO_MERCATOR_2017/NEMOGCM/NEMO/NST_SRC/agrif_opa_sponge.F90

    r7646 r8993  
    3838      ! 
    3939#if defined SPONGE 
    40       timecoeff = REAL(Agrif_NbStepint(),wp)/Agrif_rhot() 
     40!!      timecoeff = REAL(Agrif_NbStepint(),wp)/Agrif_rhot() 
     41!! Assume persistence: 
     42      timecoeff = REAL(Agrif_rhot()-1,wp)/REAL(Agrif_rhot()) 
    4143 
    4244      CALL Agrif_Sponge 
     
    6163 
    6264#if defined SPONGE 
    63       timecoeff = REAL(Agrif_NbStepint(),wp)/Agrif_rhot() 
     65!!      timecoeff = REAL(Agrif_NbStepint(),wp)/Agrif_rhot() 
     66!! Assume persistence: 
     67      timecoeff = REAL(Agrif_rhot()-1,wp)/REAL(Agrif_rhot()) 
    6468 
    6569      Agrif_SpecialValue=0. 
     
    207211      ! 
    208212      IF( before ) THEN 
    209          tabres(i1:i2,j1:j2,k1:k2,n1:n2) = tsn(i1:i2,j1:j2,k1:k2,n1:n2) 
     213         tabres(i1:i2,j1:j2,k1:k2,n1:n2) = tsb(i1:i2,j1:j2,k1:k2,n1:n2) 
    210214      ELSE    
    211215         ! 
     
    276280      ! 
    277281      IF( before ) THEN 
    278          tabres = un(i1:i2,j1:j2,:) 
     282         tabres = ub(i1:i2,j1:j2,:) 
    279283      ELSE 
    280284         ubdiff(i1:i2,j1:j2,:) = (ub(i1:i2,j1:j2,:) - tabres(:,:,:))*umask(i1:i2,j1:j2,:) 
     
    373377 
    374378      IF( before ) THEN  
    375          tabres = vn(i1:i2,j1:j2,:) 
     379         tabres = vb(i1:i2,j1:j2,:) 
    376380      ELSE 
    377381         ! 
  • branches/2017/dev_METO_MERCATOR_2017/NEMOGCM/NEMO/NST_SRC/agrif_opa_update.F90

    r7646 r8993  
    11#define TWO_WAY        /* TWO WAY NESTING */ 
    2 #undef DECAL_FEEDBACK /* SEPARATION of INTERFACES*/ 
     2#undef DECAL_FEEDBACK  /* SEPARATION of INTERFACES*/ 
     3#undef VOL_REFLUX      /* VOLUME REFLUXING*/ 
    34  
    45MODULE agrif_opa_update 
     
    1213   USE wrk_nemo   
    1314   USE zdf_oce        ! vertical physics: ocean variables  
     15   USE domvvl         ! Need interpolation routines  
    1416 
    1517   IMPLICIT NONE 
    1618   PRIVATE 
    1719 
    18    PUBLIC Agrif_Update_Tra, Agrif_Update_Dyn,Update_Scales 
     20   PUBLIC Agrif_Update_Tra, Agrif_Update_Dyn, Update_Scales, Agrif_Update_vvl, Agrif_Update_ssh 
     21 
    1922# if defined key_zdftke 
    2023   PUBLIC Agrif_Update_Tke 
     
    2730CONTAINS 
    2831 
    29    RECURSIVE SUBROUTINE Agrif_Update_Tra( ) 
     32   SUBROUTINE Agrif_Update_Tra( ) 
    3033      !!--------------------------------------------- 
    3134      !!   *** ROUTINE Agrif_Update_Tra *** 
     
    3538      ! 
    3639#if defined TWO_WAY   
    37       IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update tracers  from grid Number',Agrif_Fixed(), 'nbcline', nbcline 
     40      IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update tracers  from grid Number',Agrif_Fixed() 
    3841 
    3942      Agrif_UseSpecialValueInUpdate = .TRUE. 
    4043      Agrif_SpecialValueFineGrid = 0. 
    4144      !  
    42       IF (MOD(nbcline,nbclineupdate) == 0) THEN 
    4345# if ! defined DECAL_FEEDBACK 
    44          CALL Agrif_Update_Variable(tsn_id, procname=updateTS) 
     46      CALL Agrif_Update_Variable(tsn_id, procname=updateTS) 
     47! near boundary update: 
     48!      CALL Agrif_Update_Variable(tsn_id,locupdate=(/0,2/), procname=updateTS) 
    4549# else 
    46          CALL Agrif_Update_Variable(tsn_id, locupdate=(/1,0/),procname=updateTS) 
     50      CALL Agrif_Update_Variable(tsn_id, locupdate=(/1,0/),procname=updateTS) 
     51! near boundary update: 
     52!      CALL Agrif_Update_Variable(tsn_id,locupdate=(/1,2/), procname=updateTS) 
    4753# endif 
    48       ELSE 
    49 # if ! defined DECAL_FEEDBACK 
    50          CALL Agrif_Update_Variable(tsn_id,locupdate=(/0,2/), procname=updateTS) 
    51 # else 
    52          CALL Agrif_Update_Variable(tsn_id,locupdate=(/1,2/), procname=updateTS) 
    53 # endif 
    54       ENDIF 
    5554      ! 
    5655      Agrif_UseSpecialValueInUpdate = .FALSE. 
    5756      ! 
    58       IF ( lk_agrif_doupd ) THEN ! Initialisation: do recursive update: 
    59          CALL Agrif_ChildGrid_To_ParentGrid() 
    60          CALL Agrif_Update_Tra() 
    61          CALL Agrif_ParentGrid_To_ChildGrid() 
    62       ENDIF 
    63       ! 
    6457#endif 
    6558      ! 
    6659   END SUBROUTINE Agrif_Update_Tra 
    6760 
    68  
    69    RECURSIVE SUBROUTINE Agrif_Update_Dyn( ) 
     61   SUBROUTINE Agrif_Update_Dyn( ) 
    7062      !!--------------------------------------------- 
    7163      !!   *** ROUTINE Agrif_Update_Dyn *** 
     
    7567      ! 
    7668#if defined TWO_WAY 
    77       IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update momentum from grid Number',Agrif_Fixed(), 'nbcline', nbcline 
     69      IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update momentum from grid Number',Agrif_Fixed() 
    7870 
    7971      Agrif_UseSpecialValueInUpdate = .FALSE. 
    8072      Agrif_SpecialValueFineGrid = 0. 
    8173      !      
    82       IF (mod(nbcline,nbclineupdate) == 0) THEN 
    8374# if ! defined DECAL_FEEDBACK 
    84          CALL Agrif_Update_Variable(un_update_id,procname = updateU) 
    85          CALL Agrif_Update_Variable(vn_update_id,procname = updateV) 
     75      CALL Agrif_Update_Variable(un_update_id,procname = updateU) 
     76      CALL Agrif_Update_Variable(vn_update_id,procname = updateV) 
     77! near boundary update: 
     78!      CALL Agrif_Update_Variable(un_update_id,locupdate=(/0,1/),procname = updateU) 
     79!      CALL Agrif_Update_Variable(vn_update_id,locupdate=(/0,1/),procname = updateV) 
    8680# else 
    87          CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateU) 
    88          CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updateV) 
     81      CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateU) 
     82      CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updateV) 
     83! near boundary update: 
     84!      CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,1/),locupdate2=(/1,1/),procname = updateU) 
     85!      CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/1,1/),locupdate2=(/0,1/),procname = updateV) 
    8986# endif 
    90       ELSE 
    91 # if ! defined DECAL_FEEDBACK 
    92          CALL Agrif_Update_Variable(un_update_id,locupdate=(/0,1/),procname = updateU) 
    93          CALL Agrif_Update_Variable(vn_update_id,locupdate=(/0,1/),procname = updateV)          
    94 # else 
    95          CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,1/),locupdate2=(/1,1/),procname = updateU) 
    96          CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/1,1/),locupdate2=(/0,1/),procname = updateV) 
    97 # endif 
    98       ENDIF 
    9987 
    10088# if ! defined DECAL_FEEDBACK 
     
    10593      CALL Agrif_Update_Variable(e2v_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updateV2d)   
    10694# endif 
    107  
     95      ! 
     96# if ! defined DECAL_FEEDBACK 
     97      ! Account for updated thicknesses at boundary edges 
     98      IF (.NOT.ln_linssh) THEN 
     99! For the time being calls below do not ensure reproducible results  
     100!         CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,0/),locupdate2=(/0,0/),procname = correct_u_bdy) 
     101!         CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/0,0/),locupdate2=(/0,0/),procname = correct_v_bdy) 
     102      ENDIF 
     103# endif 
     104      !  
    108105      IF ( ln_dynspg_ts.AND.ln_bt_fw ) THEN 
    109106         ! Update time integrated transports 
    110          IF (mod(nbcline,nbclineupdate) == 0) THEN 
    111107#  if ! defined DECAL_FEEDBACK 
    112             CALL Agrif_Update_Variable(ub2b_update_id,procname = updateub2b) 
    113             CALL Agrif_Update_Variable(vb2b_update_id,procname = updatevb2b) 
     108         CALL Agrif_Update_Variable(ub2b_update_id,procname = updateub2b) 
     109         CALL Agrif_Update_Variable(vb2b_update_id,procname = updatevb2b) 
    114110#  else 
    115             CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateub2b) 
    116             CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updatevb2b) 
     111         CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateub2b) 
     112         CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updatevb2b) 
    117113#  endif 
    118          ELSE 
    119 #  if ! defined DECAL_FEEDBACK 
    120             CALL Agrif_Update_Variable(ub2b_update_id,locupdate=(/0,1/),procname = updateub2b) 
    121             CALL Agrif_Update_Variable(vb2b_update_id,locupdate=(/0,1/),procname = updatevb2b) 
    122 #  else 
    123             CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/0,1/),locupdate2=(/1,1/),procname = updateub2b) 
    124             CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1,1/),locupdate2=(/0,1/),procname = updatevb2b) 
    125 #  endif 
    126          ENDIF 
    127114      END IF 
    128       ! 
    129       nbcline = nbcline + 1 
     115#endif 
     116      ! 
     117   END SUBROUTINE Agrif_Update_Dyn 
     118 
     119   SUBROUTINE Agrif_Update_ssh( ) 
     120      !!--------------------------------------------- 
     121      !!   *** ROUTINE Agrif_Update_ssh *** 
     122      !!--------------------------------------------- 
     123      !  
     124      IF (Agrif_Root()) RETURN 
     125      ! 
     126#if defined TWO_WAY 
    130127      ! 
    131128      Agrif_UseSpecialValueInUpdate = .TRUE. 
     
    136133      CALL Agrif_Update_Variable(sshn_id,locupdate=(/1,0/),procname = updateSSH) 
    137134# endif 
     135      ! 
    138136      Agrif_UseSpecialValueInUpdate = .FALSE. 
    139       !  
     137      ! 
     138#  if defined DECAL_FEEDBACK && defined VOL_REFLUX 
     139      IF ( ln_dynspg_ts.AND.ln_bt_fw ) THEN 
     140         ! Refluxing on ssh: 
     141         CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/0, 0/),locupdate2=(/1, 1/),procname = reflux_sshu) 
     142         CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1, 1/),locupdate2=(/0, 0/),procname = reflux_sshv) 
     143      END IF 
     144#  endif 
     145      ! 
    140146#endif 
    141147      ! 
    142       ! Do recursive update: 
    143       IF ( lk_agrif_doupd ) THEN ! Initialisation: do recursive update: 
    144          CALL Agrif_ChildGrid_To_ParentGrid() 
    145          CALL Agrif_Update_Dyn() 
    146          CALL Agrif_ParentGrid_To_ChildGrid() 
    147       ENDIF 
    148       ! 
    149    END SUBROUTINE Agrif_Update_Dyn 
     148   END SUBROUTINE Agrif_Update_ssh 
    150149 
    151150# if defined key_zdftke 
     
    156155      !!--------------------------------------------- 
    157156      !! 
    158       INTEGER, INTENT(in) :: kt 
     157      INTEGER, INTENT(in) :: kt  
     158      !  
     159      IF (Agrif_Root()) RETURN 
    159160      !        
    160161      IF( ( Agrif_NbStepint() .NE. 0 ) .AND. (kt /= 0) ) RETURN 
     
    176177# endif /* key_zdftke */ 
    177178 
     179   SUBROUTINE Agrif_Update_vvl( ) 
     180      !!--------------------------------------------- 
     181      !!   *** ROUTINE Agrif_Update_vvl *** 
     182      !!--------------------------------------------- 
     183      ! 
     184      IF (Agrif_Root()) RETURN 
     185      ! 
     186#if defined TWO_WAY   
     187      ! 
     188      IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update e3 from grid Number',Agrif_Fixed(), 'Step', Agrif_Nb_Step() 
     189      ! 
     190      Agrif_UseSpecialValueInUpdate = .TRUE. 
     191      Agrif_SpecialValueFineGrid = 0. 
     192      !  
     193      ! No interface separation here, update vertical grid at T points  
     194      ! everywhere over the overlapping regions (one account for refluxing in that case): 
     195      CALL Agrif_Update_Variable(e3t_id, procname=updatee3t)  
     196      ! 
     197      Agrif_UseSpecialValueInUpdate = .FALSE. 
     198      ! 
     199      CALL Agrif_ChildGrid_To_ParentGrid() 
     200      CALL dom_vvl_update_UVF 
     201      CALL Agrif_ParentGrid_To_ChildGrid() 
     202      ! 
     203#endif 
     204      ! 
     205   END SUBROUTINE Agrif_Update_vvl 
     206 
     207   SUBROUTINE dom_vvl_update_UVF 
     208      !!--------------------------------------------- 
     209      !!       *** ROUTINE dom_vvl_update_UVF *** 
     210      !!--------------------------------------------- 
     211      !! 
     212      INTEGER :: jk 
     213      REAL(wp):: zcoef 
     214      !!--------------------------------------------- 
     215 
     216      IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Finalize e3 on grid Number', & 
     217                  & Agrif_Fixed(), 'Step', Agrif_Nb_Step() 
     218 
     219      ! Save "old" scale factor (prior update) for subsequent asselin correction 
     220      ! of prognostic variables 
     221      ! ----------------------- 
     222      ! 
     223      e3u_a(:,:,:) = e3u_n(:,:,:) 
     224      e3v_a(:,:,:) = e3v_n(:,:,:) 
     225!      ua(:,:,:) = e3u_b(:,:,:) 
     226!      va(:,:,:) = e3v_b(:,:,:) 
     227      hu_a(:,:) = hu_n(:,:) 
     228      hv_a(:,:) = hv_n(:,:) 
     229 
     230      ! 1) NOW fields 
     231      !-------------- 
     232       
     233         ! Vertical scale factor interpolations 
     234         ! ------------------------------------ 
     235      CALL dom_vvl_interpol( e3t_n(:,:,:), e3u_n(:,:,:) ,  'U' ) 
     236      CALL dom_vvl_interpol( e3t_n(:,:,:), e3v_n(:,:,:) ,  'V' ) 
     237      CALL dom_vvl_interpol( e3u_n(:,:,:), e3f_n(:,:,:) ,  'F' ) 
     238 
     239      CALL dom_vvl_interpol( e3u_n(:,:,:), e3uw_n(:,:,:), 'UW' ) 
     240      CALL dom_vvl_interpol( e3v_n(:,:,:), e3vw_n(:,:,:), 'VW' ) 
     241 
     242         ! Update total depths: 
     243         ! -------------------- 
     244      hu_n(:,:) = 0._wp                        ! Ocean depth at U-points 
     245      hv_n(:,:) = 0._wp                        ! Ocean depth at V-points 
     246      DO jk = 1, jpkm1 
     247         hu_n(:,:) = hu_n(:,:) + e3u_n(:,:,jk) * umask(:,:,jk) 
     248         hv_n(:,:) = hv_n(:,:) + e3v_n(:,:,jk) * vmask(:,:,jk) 
     249      END DO 
     250      !                                        ! Inverse of the local depth 
     251      r1_hu_n(:,:) = ssumask(:,:) / ( hu_n(:,:) + 1._wp - ssumask(:,:) ) 
     252      r1_hv_n(:,:) = ssvmask(:,:) / ( hv_n(:,:) + 1._wp - ssvmask(:,:) ) 
     253 
     254 
     255      ! 2) BEFORE fields: 
     256      !------------------ 
     257!      IF (     (.NOT.(lk_agrif_fstep.AND.(neuler==0)).AND.(ln_dynspg_exp)) & 
     258!         & .OR.(.NOT.(lk_agrif_fstep.AND.(neuler==0)).AND.(ln_dynspg_ts    & 
     259!         & .AND.(.NOT.ln_bt_fw)))) THEN 
     260      IF (.NOT.(lk_agrif_fstep.AND.(neuler==0) )) THEN 
     261         ! 
     262         ! Vertical scale factor interpolations 
     263         ! ------------------------------------ 
     264         CALL dom_vvl_interpol( e3t_b(:,:,:), e3u_b(:,:,:),  'U'  ) 
     265         CALL dom_vvl_interpol( e3t_b(:,:,:), e3v_b(:,:,:),  'V'  ) 
     266 
     267         CALL dom_vvl_interpol( e3u_b(:,:,:), e3uw_b(:,:,:), 'UW' ) 
     268         CALL dom_vvl_interpol( e3v_b(:,:,:), e3vw_b(:,:,:), 'VW' ) 
     269 
     270         ! Update total depths: 
     271         ! -------------------- 
     272         hu_b(:,:) = 0._wp                     ! Ocean depth at U-points 
     273         hv_b(:,:) = 0._wp                     ! Ocean depth at V-points 
     274         DO jk = 1, jpkm1 
     275            hu_b(:,:) = hu_b(:,:) + e3u_b(:,:,jk) * umask(:,:,jk) 
     276            hv_b(:,:) = hv_b(:,:) + e3v_b(:,:,jk) * vmask(:,:,jk) 
     277         END DO 
     278         !                                     ! Inverse of the local depth 
     279         r1_hu_b(:,:) = ssumask(:,:) / ( hu_b(:,:) + 1._wp - ssumask(:,:) ) 
     280         r1_hv_b(:,:) = ssvmask(:,:) / ( hv_b(:,:) + 1._wp - ssvmask(:,:) ) 
     281      ENDIF 
     282      ! 
     283   END SUBROUTINE dom_vvl_update_UVF 
     284 
    178285   SUBROUTINE updateTS( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before ) 
    179286      !!--------------------------------------------- 
     
    185292      !! 
    186293      INTEGER :: ji,jj,jk,jn 
     294      REAL(wp) :: ztb, ztnu, ztno 
    187295      !!--------------------------------------------- 
    188296      ! 
     
    192300               DO jj=j1,j2 
    193301                  DO ji=i1,i2 
    194                      tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn) 
     302!> jc tmp 
     303                     tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn)  * e3t_n(ji,jj,jk) / e3t_0(ji,jj,jk) 
     304!                     tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn)  * e3t_n(ji,jj,jk) 
     305!< jc tmp 
    195306                  END DO 
    196307               END DO 
     
    198309         END DO 
    199310      ELSE 
     311!> jc tmp 
     312         DO jn = n1,n2 
     313            tabres(i1:i2,j1:j2,k1:k2,jn) =  tabres(i1:i2,j1:j2,k1:k2,jn) * e3t_0(i1:i2,j1:j2,k1:k2) & 
     314                                         & * tmask(i1:i2,j1:j2,k1:k2) 
     315         ENDDO 
     316!< jc tmp 
    200317         IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN 
    201318            ! Add asselin part 
     
    205322                     DO ji=i1,i2 
    206323                        IF( tabres(ji,jj,jk,jn) .NE. 0. ) THEN 
    207                            tsb(ji,jj,jk,jn) = tsb(ji,jj,jk,jn) &  
    208                                  & + atfp * ( tabres(ji,jj,jk,jn) & 
    209                                  &             - tsn(ji,jj,jk,jn) ) * tmask(ji,jj,jk) 
     324                           ztb  = tsb(ji,jj,jk,jn) * e3t_b(ji,jj,jk) ! fse3t_b prior update should be used 
     325                           ztnu = tabres(ji,jj,jk,jn) 
     326                           ztno = tsn(ji,jj,jk,jn) * e3t_a(ji,jj,jk) 
     327                           tsb(ji,jj,jk,jn) = ( ztb + atfp * ( ztnu - ztno) )  &  
     328                                     &        * tmask(ji,jj,jk) / e3t_b(ji,jj,jk) 
    210329                        ENDIF 
    211330                     ENDDO 
     
    219338                  DO ji=i1,i2 
    220339                     IF( tabres(ji,jj,jk,jn) .NE. 0. ) THEN  
    221                         tsn(ji,jj,jk,jn) = tabres(ji,jj,jk,jn) * tmask(ji,jj,jk) 
     340                        tsn(ji,jj,jk,jn) = tabres(ji,jj,jk,jn) / e3t_n(ji,jj,jk) 
    222341                     END IF 
    223342                  END DO 
     
    225344            END DO 
    226345         END DO 
     346         ! 
     347         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 
     348            tsb(i1:i2,j1:j2,k1:k2,n1:n2)  = tsn(i1:i2,j1:j2,k1:k2,n1:n2) 
     349         ENDIF 
     350         ! 
    227351      ENDIF 
    228352      !  
     
    238362      LOGICAL                               , INTENT(in   ) :: before 
    239363      ! 
    240       INTEGER  ::   ji, jj, jk 
    241       REAL(wp) ::   zrhoy 
     364      INTEGER  :: ji, jj, jk 
     365      REAL(wp) :: zrhoy, zub, zunu, zuno 
    242366      !!--------------------------------------------- 
    243367      !  
     
    251375            DO jj=j1,j2 
    252376               DO ji=i1,i2 
    253                   tabres(ji,jj,jk) = tabres(ji,jj,jk) * r1_e2u(ji,jj) / e3u_n(ji,jj,jk) 
     377                  tabres(ji,jj,jk) = tabres(ji,jj,jk) * r1_e2u(ji,jj)  
    254378                  ! 
    255379                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 
    256                      ub(ji,jj,jk) = ub(ji,jj,jk) &  
    257                            & + atfp * ( tabres(ji,jj,jk) - un(ji,jj,jk) ) * umask(ji,jj,jk) 
     380                     zub  = ub(ji,jj,jk) * e3u_b(ji,jj,jk)  ! fse3t_b prior update should be used 
     381                     zuno = un(ji,jj,jk) * e3u_a(ji,jj,jk) 
     382                     zunu = tabres(ji,jj,jk) 
     383                     ub(ji,jj,jk) = ( zub + atfp * ( zunu - zuno) ) &       
     384                                    & * umask(ji,jj,jk) / e3u_b(ji,jj,jk) 
    258385                  ENDIF 
    259386                  ! 
    260                   un(ji,jj,jk) = tabres(ji,jj,jk) * umask(ji,jj,jk) 
    261                END DO 
    262             END DO 
    263          END DO 
     387                  un(ji,jj,jk) = tabres(ji,jj,jk) * umask(ji,jj,jk) / e3u_n(ji,jj,jk) 
     388               END DO 
     389            END DO 
     390         END DO 
     391         ! 
     392         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 
     393            ub(i1:i2,j1:j2,k1:k2)  = un(i1:i2,j1:j2,k1:k2) 
     394         ENDIF 
     395         ! 
    264396      ENDIF 
    265397      !  
    266398   END SUBROUTINE updateu 
    267399 
    268  
    269    SUBROUTINE updatev( tabres, i1, i2, j1, j2, k1, k2, before ) 
     400   SUBROUTINE correct_u_bdy( tabres, i1, i2, j1, j2, k1, k2, before, nb, ndir ) 
     401      !!--------------------------------------------- 
     402      !!           *** ROUTINE correct_u_bdy *** 
     403      !!--------------------------------------------- 
     404      INTEGER                               , INTENT(in   ) :: i1, i2, j1, j2, k1, k2 
     405      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres 
     406      LOGICAL                               , INTENT(in   ) :: before 
     407      INTEGER                               , INTENT(in)    :: nb, ndir 
     408      !! 
     409      LOGICAL :: western_side, eastern_side  
     410      ! 
     411      INTEGER  :: jj, jk 
     412      REAL(wp) :: zcor 
     413      !!--------------------------------------------- 
     414      !  
     415      IF( .NOT.before ) THEN 
     416         ! 
     417         western_side  = (nb == 1).AND.(ndir == 1) 
     418         eastern_side  = (nb == 1).AND.(ndir == 2) 
     419         ! 
     420         IF (western_side) THEN 
     421            DO jj=j1,j2 
     422               zcor = un_b(i1-1,jj) * hu_a(i1-1,jj) * r1_hu_n(i1-1,jj) - un_b(i1-1,jj) 
     423               un_b(i1-1,jj) = un_b(i1-1,jj) + zcor 
     424               DO jk=1,jpkm1 
     425                  un(i1-1,jj,jk) = un(i1-1,jj,jk) + zcor * umask(i1-1,jj,jk) 
     426               END DO  
     427            END DO 
     428         ENDIF 
     429         ! 
     430         IF (eastern_side) THEN 
     431            DO jj=j1,j2 
     432               zcor = un_b(i2+1,jj) * hu_a(i2+1,jj) * r1_hu_n(i2+1,jj) - un_b(i2+1,jj) 
     433               un_b(i2+1,jj) = un_b(i2+1,jj) + zcor 
     434               DO jk=1,jpkm1 
     435                  un(i2+1,jj,jk) = un(i2+1,jj,jk) + zcor * umask(i2+1,jj,jk) 
     436               END DO  
     437            END DO 
     438         ENDIF 
     439         ! 
     440      ENDIF 
     441      !  
     442   END SUBROUTINE correct_u_bdy 
     443 
     444 
     445   SUBROUTINE updatev( tabres, i1, i2, j1, j2, k1, k2, before) 
    270446      !!--------------------------------------------- 
    271447      !!           *** ROUTINE updatev *** 
    272448      !!--------------------------------------------- 
    273       INTEGER :: i1,i2,j1,j2,k1,k2 
    274       INTEGER :: ji,jj,jk 
    275       REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: tabres 
    276       LOGICAL :: before 
    277       !! 
    278       REAL(wp) :: zrhox 
     449      INTEGER                               , INTENT(in   ) :: i1, i2, j1, j2, k1, k2 
     450      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres 
     451      LOGICAL                               , INTENT(in   ) :: before 
     452      ! 
     453      INTEGER  :: ji, jj, jk 
     454      REAL(wp) :: zrhox, zvb, zvnu, zvno 
    279455      !!---------------------------------------------       
    280456      ! 
     
    292468            DO jj=j1,j2 
    293469               DO ji=i1,i2 
    294                   tabres(ji,jj,jk) = tabres(ji,jj,jk) * r1_e1v(ji,jj) / e3v_n(ji,jj,jk) 
     470                  tabres(ji,jj,jk) = tabres(ji,jj,jk) * r1_e1v(ji,jj) 
    295471                  ! 
    296472                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 
    297                      vb(ji,jj,jk) = vb(ji,jj,jk) &  
    298                            & + atfp * ( tabres(ji,jj,jk) - vn(ji,jj,jk) ) * vmask(ji,jj,jk) 
     473                     zvb  = vb(ji,jj,jk) * e3v_b(ji,jj,jk) ! fse3t_b prior update should be used 
     474                     zvno = vn(ji,jj,jk) * e3v_a(ji,jj,jk) 
     475                     zvnu = tabres(ji,jj,jk) 
     476                     vb(ji,jj,jk) = ( zvb + atfp * ( zvnu - zvno) ) &       
     477                                    & * vmask(ji,jj,jk) / e3v_b(ji,jj,jk) 
    299478                  ENDIF 
    300479                  ! 
    301                   vn(ji,jj,jk) = tabres(ji,jj,jk) * vmask(ji,jj,jk) 
    302                END DO 
    303             END DO 
    304          END DO 
     480                  vn(ji,jj,jk) = tabres(ji,jj,jk) * vmask(ji,jj,jk) / e3v_n(ji,jj,jk) 
     481               END DO 
     482            END DO 
     483         END DO 
     484         ! 
     485         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 
     486            vb(i1:i2,j1:j2,k1:k2)  = vn(i1:i2,j1:j2,k1:k2) 
     487         ENDIF 
     488         ! 
    305489      ENDIF 
    306490      !  
    307491   END SUBROUTINE updatev 
     492 
     493   SUBROUTINE correct_v_bdy( tabres, i1, i2, j1, j2, k1, k2, before, nb, ndir ) 
     494      !!--------------------------------------------- 
     495      !!           *** ROUTINE correct_u_bdy *** 
     496      !!--------------------------------------------- 
     497      INTEGER                               , INTENT(in   ) :: i1, i2, j1, j2, k1, k2 
     498      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres 
     499      LOGICAL                               , INTENT(in   ) :: before 
     500      INTEGER                               , INTENT(in)    :: nb, ndir 
     501      !! 
     502      LOGICAL :: southern_side, northern_side  
     503      ! 
     504      INTEGER  :: ji, jk 
     505      REAL(wp) :: zcor 
     506      !!--------------------------------------------- 
     507      !  
     508      IF( .NOT.before ) THEN 
     509         ! 
     510         southern_side = (nb == 2).AND.(ndir == 1) 
     511         northern_side = (nb == 2).AND.(ndir == 2) 
     512         ! 
     513         IF (southern_side) THEN 
     514            DO ji=i1,i2 
     515               zcor = vn_b(ji,j1-1) * hv_a(ji,j1-1) * r1_hv_n(ji,j1-1) - vn_b(ji,j1-1) 
     516               vn_b(ji,j1-1) = vn_b(ji,j1-1) + zcor 
     517               DO jk=1,jpkm1 
     518                  vn(ji,j1-1,jk) = vn(ji,j1-1,jk) + zcor * vmask(ji,j1-1,jk) 
     519               END DO  
     520            END DO 
     521         ENDIF 
     522         ! 
     523         IF (northern_side) THEN 
     524            DO ji=i1,i2 
     525               zcor = vn_b(ji,j2+1) * hv_a(ji,j2+1) * r1_hv_n(ji,j2+1) - vn_b(ji,j2+1) 
     526               vn_b(ji,j2+1) = vn_b(ji,j2+1) + zcor 
     527               DO jk=1,jpkm1 
     528                  vn(ji,j2+1,jk) = vn(ji,j2+1,jk) + zcor * vmask(ji,j2+1,jk) 
     529               END DO  
     530            END DO 
     531         ENDIF 
     532         ! 
     533      ENDIF 
     534      !  
     535   END SUBROUTINE correct_v_bdy 
    308536 
    309537 
     
    316544      LOGICAL, INTENT(in) :: before 
    317545      !!  
    318       INTEGER :: ji, jj, jk 
     546      INTEGER  :: ji, jj, jk 
    319547      REAL(wp) :: zrhoy 
    320548      REAL(wp) :: zcorr 
     
    331559         DO jj=j1,j2 
    332560            DO ji=i1,i2 
    333                tabres(ji,jj) =  tabres(ji,jj) * r1_hu_n(ji,jj) * r1_e2u(ji,jj)   
     561               tabres(ji,jj) =  tabres(ji,jj) * r1_e2u(ji,jj)   
    334562               !     
    335563               ! Update "now" 3d velocities: 
     
    338566                  spgu(ji,jj) = spgu(ji,jj) + e3u_n(ji,jj,jk) * un(ji,jj,jk) 
    339567               END DO 
    340                spgu(ji,jj) = spgu(ji,jj) * r1_hu_n(ji,jj) 
    341568               ! 
    342                zcorr = tabres(ji,jj) - spgu(ji,jj) 
     569               zcorr = (tabres(ji,jj) - spgu(ji,jj)) * r1_hu_n(ji,jj) 
    343570               DO jk=1,jpkm1               
    344571                  un(ji,jj,jk) = un(ji,jj,jk) + zcorr * umask(ji,jj,jk)            
     
    348575               IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN 
    349576                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 
    350                      zcorr = tabres(ji,jj) - un_b(ji,jj) 
     577                     zcorr = (tabres(ji,jj) - un_b(ji,jj) * hu_a(ji,jj)) * r1_hu_b(ji,jj) 
    351578                     ub_b(ji,jj) = ub_b(ji,jj) + atfp * zcorr * umask(ji,jj,1) 
    352579                  END IF 
    353                ENDIF              
    354                un_b(ji,jj) = tabres(ji,jj) * umask(ji,jj,1) 
     580               ENDIF     
     581               un_b(ji,jj) = tabres(ji,jj) * r1_hu_n(ji,jj) * umask(ji,jj,1) 
    355582               !        
    356583               ! Correct "before" velocities to hold correct bt component: 
     
    359586                  spgu(ji,jj) = spgu(ji,jj) + e3u_b(ji,jj,jk) * ub(ji,jj,jk) 
    360587               END DO 
    361                spgu(ji,jj) = spgu(ji,jj) * r1_hu_b(ji,jj) 
    362588               ! 
    363                zcorr = ub_b(ji,jj) - spgu(ji,jj) 
     589               zcorr = ub_b(ji,jj) - spgu(ji,jj) * r1_hu_b(ji,jj) 
    364590               DO jk=1,jpkm1               
    365591                  ub(ji,jj,jk) = ub(ji,jj,jk) + zcorr * umask(ji,jj,jk)            
     
    368594            END DO 
    369595         END DO 
     596         ! 
     597         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 
     598            ub_b(i1:i2,j1:j2)  = un_b(i1:i2,j1:j2) 
     599         ENDIF 
    370600      ENDIF 
    371601      ! 
     
    381611      LOGICAL, INTENT(in) :: before 
    382612      !!  
    383       INTEGER :: ji, jj, jk 
     613      INTEGER  :: ji, jj, jk 
    384614      REAL(wp) :: zrhox 
    385615      REAL(wp) :: zcorr 
     
    396626         DO jj=j1,j2 
    397627            DO ji=i1,i2 
    398                tabres(ji,jj) =  tabres(ji,jj) * r1_hv_n(ji,jj) * r1_e1v(ji,jj)   
     628               tabres(ji,jj) =  tabres(ji,jj) * r1_e1v(ji,jj)   
    399629               !     
    400630               ! Update "now" 3d velocities: 
     
    403633                  spgv(ji,jj) = spgv(ji,jj) + e3v_n(ji,jj,jk) * vn(ji,jj,jk) 
    404634               END DO 
    405                spgv(ji,jj) = spgv(ji,jj) * r1_hv_n(ji,jj) 
    406635               ! 
    407                zcorr = tabres(ji,jj) - spgv(ji,jj) 
     636               zcorr = (tabres(ji,jj) - spgv(ji,jj)) * r1_hv_n(ji,jj) 
    408637               DO jk=1,jpkm1               
    409638                  vn(ji,jj,jk) = vn(ji,jj,jk) + zcorr * vmask(ji,jj,jk)            
     
    413642               IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN 
    414643                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 
    415                      zcorr = tabres(ji,jj) - vn_b(ji,jj) 
     644                     zcorr = (tabres(ji,jj) - vn_b(ji,jj) * hv_a(ji,jj)) * r1_hv_b(ji,jj) 
    416645                     vb_b(ji,jj) = vb_b(ji,jj) + atfp * zcorr * vmask(ji,jj,1) 
    417646                  END IF 
    418647               ENDIF               
    419                vn_b(ji,jj) = tabres(ji,jj) * vmask(ji,jj,1) 
     648               vn_b(ji,jj) = tabres(ji,jj) * r1_hv_n(ji,jj) * vmask(ji,jj,1) 
    420649               !        
    421650               ! Correct "before" velocities to hold correct bt component: 
     
    424653                  spgv(ji,jj) = spgv(ji,jj) + e3v_b(ji,jj,jk) * vb(ji,jj,jk) 
    425654               END DO 
    426                spgv(ji,jj) = spgv(ji,jj) * r1_hv_b(ji,jj) 
    427655               ! 
    428                zcorr = vb_b(ji,jj) - spgv(ji,jj) 
     656               zcorr = vb_b(ji,jj) - spgv(ji,jj) * r1_hv_b(ji,jj) 
    429657               DO jk=1,jpkm1               
    430658                  vb(ji,jj,jk) = vb(ji,jj,jk) + zcorr * vmask(ji,jj,jk)            
     
    433661            END DO 
    434662         END DO 
     663         ! 
     664         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 
     665            vb_b(i1:i2,j1:j2)  = vn_b(i1:i2,j1:j2) 
     666         ENDIF 
     667         ! 
    435668      ENDIF 
    436669      !  
     
    456689         END DO 
    457690      ELSE 
    458          IF( .NOT.ln_dynspg_ts .OR. ( ln_dynspg_ts .AND. .NOT.ln_bt_fw ) ) THEN 
    459             IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN 
    460                DO jj=j1,j2 
    461                   DO ji=i1,i2 
    462                      sshb(ji,jj) =   sshb(ji,jj) & 
    463                            & + atfp * ( tabres(ji,jj) - sshn(ji,jj) ) * tmask(ji,jj,1) 
    464                   END DO 
    465                END DO 
    466             ENDIF 
     691         IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN 
     692            DO jj=j1,j2 
     693               DO ji=i1,i2 
     694                  sshb(ji,jj) =   sshb(ji,jj) & 
     695                        & + atfp * ( tabres(ji,jj) - sshn(ji,jj) ) * tmask(ji,jj,1) 
     696               END DO 
     697            END DO 
    467698         ENDIF 
    468699         ! 
     
    472703            END DO 
    473704         END DO 
     705         ! 
     706         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 
     707            sshb(i1:i2,j1:j2)  = sshn(i1:i2,j1:j2) 
     708         ENDIF 
     709         ! 
     710 
    474711      ENDIF 
    475712      ! 
     
    486723      !! 
    487724      INTEGER :: ji, jj 
    488       REAL(wp) :: zrhoy 
     725      REAL(wp) :: zrhoy, za1, zcor 
    489726      !!--------------------------------------------- 
    490727      ! 
     
    498735         tabres = zrhoy * tabres 
    499736      ELSE 
     737         !  
     738         tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e2u(i1:i2,j1:j2) 
     739         ! 
     740         za1 = 1._wp / REAL(Agrif_rhot(), wp) 
    500741         DO jj=j1,j2 
    501742            DO ji=i1,i2 
    502                ub2_b(ji,jj) = tabres(ji,jj) / e2u(ji,jj) 
     743               zcor=tabres(ji,jj) - ub2_b(ji,jj) 
     744               ! Update time integrated fluxes also in case of multiply nested grids: 
     745               ub2_i_b(ji,jj) = ub2_i_b(ji,jj) + za1 * zcor  
     746               ! Update corrective fluxes: 
     747               un_bf(ji,jj)  = un_bf(ji,jj) + zcor 
     748               ! Update half step back fluxes: 
     749               ub2_b(ji,jj) = tabres(ji,jj) 
    503750            END DO 
    504751         END DO 
     
    507754   END SUBROUTINE updateub2b 
    508755 
     756   SUBROUTINE reflux_sshu( tabres, i1, i2, j1, j2, before, nb, ndir ) 
     757      !!--------------------------------------------- 
     758      !!          *** ROUTINE reflux_sshu *** 
     759      !!--------------------------------------------- 
     760      INTEGER, INTENT(in) :: i1, i2, j1, j2 
     761      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres 
     762      LOGICAL, INTENT(in) :: before 
     763      INTEGER, INTENT(in) :: nb, ndir 
     764      !! 
     765      LOGICAL :: western_side, eastern_side  
     766      INTEGER :: ji, jj 
     767      REAL(wp) :: zrhoy, za1, zcor 
     768      !!--------------------------------------------- 
     769      ! 
     770      IF (before) THEN 
     771         zrhoy = Agrif_Rhoy() 
     772         DO jj=j1,j2 
     773            DO ji=i1,i2 
     774               tabres(ji,jj) = ub2_i_b(ji,jj) * e2u(ji,jj) 
     775            END DO 
     776         END DO 
     777         tabres = zrhoy * tabres 
     778      ELSE 
     779         !  
     780         tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e2u(i1:i2,j1:j2) 
     781         ! 
     782         western_side  = (nb == 1).AND.(ndir == 1) 
     783         eastern_side  = (nb == 1).AND.(ndir == 2) 
     784         ! 
     785         IF (western_side) THEN 
     786            DO jj=j1,j2 
     787               zcor = rdt * r1_e1e2t(i1  ,jj) * e2u(i1,jj) * (ub2_b(i1,jj)-tabres(i1,jj))  
     788               sshn(i1  ,jj) = sshn(i1  ,jj) + zcor 
     789               IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(i1  ,jj) = sshb(i1  ,jj) + atfp * zcor 
     790            END DO 
     791         ENDIF 
     792         IF (eastern_side) THEN 
     793            DO jj=j1,j2 
     794               zcor = - rdt * r1_e1e2t(i2+1,jj) * e2u(i2,jj) * (ub2_b(i2,jj)-tabres(i2,jj)) 
     795               sshn(i2+1,jj) = sshn(i2+1,jj) + zcor 
     796               IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(i2+1,jj) = sshb(i2+1,jj) + atfp * zcor 
     797            END DO 
     798         ENDIF 
     799         ! 
     800      ENDIF 
     801      ! 
     802   END SUBROUTINE reflux_sshu 
    509803 
    510804   SUBROUTINE updatevb2b( tabres, i1, i2, j1, j2, before ) 
     
    517811      !! 
    518812      INTEGER :: ji, jj 
    519       REAL(wp) :: zrhox 
     813      REAL(wp) :: zrhox, za1, zcor 
    520814      !!--------------------------------------------- 
    521815      ! 
     
    529823         tabres = zrhox * tabres 
    530824      ELSE 
     825         !  
     826         tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e1v(i1:i2,j1:j2) 
     827         ! 
     828         za1 = 1._wp / REAL(Agrif_rhot(), wp) 
    531829         DO jj=j1,j2 
    532830            DO ji=i1,i2 
    533                vb2_b(ji,jj) = tabres(ji,jj) / e1v(ji,jj) 
     831               zcor=tabres(ji,jj) - vb2_b(ji,jj) 
     832               ! Update time integrated fluxes also in case of multiply nested grids: 
     833               vb2_i_b(ji,jj) = vb2_i_b(ji,jj) + za1 * zcor  
     834               ! Update corrective fluxes: 
     835               vn_bf(ji,jj)  = vn_bf(ji,jj) + zcor 
     836               ! Update half step back fluxes: 
     837               vb2_b(ji,jj) = tabres(ji,jj) 
    534838            END DO 
    535839         END DO 
     
    538842   END SUBROUTINE updatevb2b 
    539843 
     844   SUBROUTINE reflux_sshv( tabres, i1, i2, j1, j2, before, nb, ndir ) 
     845      !!--------------------------------------------- 
     846      !!          *** ROUTINE reflux_sshv *** 
     847      !!--------------------------------------------- 
     848      INTEGER, INTENT(in) :: i1, i2, j1, j2 
     849      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres 
     850      LOGICAL, INTENT(in) :: before 
     851      INTEGER, INTENT(in) :: nb, ndir 
     852      !! 
     853      LOGICAL :: southern_side, northern_side  
     854      INTEGER :: ji, jj 
     855      REAL(wp) :: zrhox, za1, zcor 
     856      !!--------------------------------------------- 
     857      ! 
     858      IF (before) THEN 
     859         zrhox = Agrif_Rhox() 
     860         DO jj=j1,j2 
     861            DO ji=i1,i2 
     862               tabres(ji,jj) = vb2_i_b(ji,jj) * e1v(ji,jj)  
     863            END DO 
     864         END DO 
     865         tabres = zrhox * tabres 
     866      ELSE 
     867         !  
     868         tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e1v(i1:i2,j1:j2) 
     869         ! 
     870         southern_side = (nb == 2).AND.(ndir == 1) 
     871         northern_side = (nb == 2).AND.(ndir == 2) 
     872         ! 
     873         IF (southern_side) THEN 
     874            DO ji=i1,i2 
     875               zcor = rdt * r1_e1e2t(ji,j1  ) * e1v(ji,j1  ) * (vb2_b(ji,j1)-tabres(ji,j1)) 
     876               sshn(ji,j1  ) = sshn(ji,j1  ) + zcor 
     877               IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(ji,j1  ) = sshb(ji,j1) + atfp * zcor 
     878            END DO 
     879         ENDIF 
     880         IF (northern_side) THEN                
     881            DO ji=i1,i2 
     882               zcor = - rdt * r1_e1e2t(ji,j2+1) * e1v(ji,j2  ) * (vb2_b(ji,j2)-tabres(ji,j2)) 
     883               sshn(ji,j2+1) = sshn(ji,j2+1) + zcor 
     884               IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(ji,j2+1) = sshb(ji,j2+1) + atfp * zcor 
     885            END DO 
     886         ENDIF 
     887         !  
     888      ENDIF 
     889      ! 
     890   END SUBROUTINE reflux_sshv 
    540891 
    541892   SUBROUTINE update_scales( tabres, i1, i2, j1, j2, k1, k2, n1,n2, before ) 
     
    644995# endif /* key_zdftke */  
    645996 
     997   SUBROUTINE updatee3t(ptab_dum, i1, i2, j1, j2, k1, k2, before ) 
     998      !!--------------------------------------------- 
     999      !!           *** ROUTINE updatee3t *** 
     1000      !!--------------------------------------------- 
     1001      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: ptab_dum 
     1002      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2 
     1003      LOGICAL, INTENT(in) :: before 
     1004      ! 
     1005      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptab 
     1006      INTEGER :: ji,jj,jk 
     1007      REAL(wp) :: zcoef 
     1008      !!--------------------------------------------- 
     1009      ! 
     1010      IF (.NOT.before) THEN 
     1011         ! 
     1012         ALLOCATE(ptab(i1:i2,j1:j2,1:jpk))  
     1013         ! 
     1014         ! Update e3t from ssh (z* case only) 
     1015         DO jk = 1, jpkm1 
     1016            DO jj=j1,j2 
     1017               DO ji=i1,i2 
     1018                  ptab(ji,jj,jk) = e3t_0(ji,jj,jk) * (1._wp + sshn(ji,jj) & 
     1019                                     & *ssmask(ji,jj)/(ht_0(ji,jj)-1._wp + ssmask(ji,jj))) 
     1020               END DO 
     1021            END DO 
     1022         END DO 
     1023         ! 
     1024         ! 1) Updates at BEFORE time step: 
     1025         ! ------------------------------- 
     1026         ! 
     1027         ! Save "old" scale factor (prior update) for subsequent asselin correction 
     1028         ! of prognostic variables 
     1029         e3t_a(i1:i2,j1:j2,1:jpkm1) = e3t_n(i1:i2,j1:j2,1:jpkm1) 
     1030 
     1031         ! One should also save e3t_b, but lacking of workspace... 
     1032!         hdivn(i1:i2,j1:j2,1:jpkm1)   = e3t_b(i1:i2,j1:j2,1:jpkm1) 
     1033 
     1034         IF (.NOT.(lk_agrif_fstep.AND.(neuler==0) )) THEN 
     1035            DO jk = 1, jpkm1 
     1036               DO jj=j1,j2 
     1037                  DO ji=i1,i2 
     1038                     e3t_b(ji,jj,jk) =  e3t_b(ji,jj,jk) & 
     1039                           & + atfp * ( ptab(ji,jj,jk) - e3t_n(ji,jj,jk) ) 
     1040                  END DO 
     1041               END DO 
     1042            END DO 
     1043            ! 
     1044            e3w_b  (i1:i2,j1:j2,1) = e3w_0(i1:i2,j1:j2,1) + e3t_b(i1:i2,j1:j2,1) - e3t_0(i1:i2,j1:j2,1) 
     1045            gdepw_b(i1:i2,j1:j2,1) = 0.0_wp 
     1046            gdept_b(i1:i2,j1:j2,1) = 0.5_wp * e3w_b(i1:i2,j1:j2,1) 
     1047            ! 
     1048            DO jk = 2, jpk 
     1049               DO jj = j1,j2 
     1050                  DO ji = i1,i2             
     1051                     zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk)) 
     1052                     e3w_b(ji,jj,jk)  = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) *        &  
     1053                     &                                        ( e3t_b(ji,jj,jk-1) - e3t_0(ji,jj,jk-1) )  & 
     1054                     &                                  +            0.5_wp * tmask(ji,jj,jk)   *        & 
     1055                     &                                        ( e3t_b(ji,jj,jk  ) - e3t_0(ji,jj,jk  ) ) 
     1056                     gdepw_b(ji,jj,jk) = gdepw_b(ji,jj,jk-1) + e3t_b(ji,jj,jk-1) 
     1057                     gdept_b(ji,jj,jk) =      zcoef  * ( gdepw_b(ji,jj,jk  ) + 0.5 * e3w_b(ji,jj,jk))  & 
     1058                         &               + (1-zcoef) * ( gdept_b(ji,jj,jk-1) +       e3w_b(ji,jj,jk))  
     1059                  END DO 
     1060               END DO 
     1061            END DO 
     1062            ! 
     1063         ENDIF         
     1064         ! 
     1065         ! 2) Updates at NOW time step: 
     1066         ! ---------------------------- 
     1067         ! 
     1068         ! Update vertical scale factor at T-points: 
     1069         e3t_n(i1:i2,j1:j2,1:jpkm1) = ptab(i1:i2,j1:j2,1:jpkm1) 
     1070         ! 
     1071         ! Update total depth: 
     1072         ht_n(i1:i2,j1:j2) = 0._wp 
     1073         DO jk = 1, jpkm1 
     1074            ht_n(i1:i2,j1:j2) = ht_n(i1:i2,j1:j2) + e3t_n(i1:i2,j1:j2,jk) * tmask(i1:i2,j1:j2,jk) 
     1075         END DO 
     1076         ! 
     1077         ! Update vertical scale factor at W-points and depths: 
     1078         e3w_n (i1:i2,j1:j2,1) = e3w_0(i1:i2,j1:j2,1) + e3t_n(i1:i2,j1:j2,1) - e3t_0(i1:i2,j1:j2,1) 
     1079         gdept_n(i1:i2,j1:j2,1) = 0.5_wp * e3w_n(i1:i2,j1:j2,1) 
     1080         gdepw_n(i1:i2,j1:j2,1) = 0.0_wp 
     1081         gde3w_n(i1:i2,j1:j2,1) = gdept_n(i1:i2,j1:j2,1) - (ht_n(i1:i2,j1:j2)-ht_0(i1:i2,j1:j2)) ! Last term in the rhs is ssh 
     1082         ! 
     1083         DO jk = 2, jpk 
     1084            DO jj = j1,j2 
     1085               DO ji = i1,i2             
     1086               zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk)) 
     1087               e3w_n(ji,jj,jk)  = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) * ( e3t_n(ji,jj,jk-1) - e3t_0(ji,jj,jk-1) )   & 
     1088               &                                  +            0.5_wp * tmask(ji,jj,jk)   * ( e3t_n(ji,jj,jk  ) - e3t_0(ji,jj,jk  ) ) 
     1089               gdepw_n(ji,jj,jk) = gdepw_n(ji,jj,jk-1) + e3t_n(ji,jj,jk-1) 
     1090               gdept_n(ji,jj,jk) =      zcoef  * ( gdepw_n(ji,jj,jk  ) + 0.5 * e3w_n(ji,jj,jk))  & 
     1091                   &               + (1-zcoef) * ( gdept_n(ji,jj,jk-1) +       e3w_n(ji,jj,jk))  
     1092               gde3w_n(ji,jj,jk) = gdept_n(ji,jj,jk) - (ht_n(ji,jj)-ht_0(ji,jj)) ! Last term in the rhs is ssh 
     1093               END DO 
     1094            END DO 
     1095         END DO 
     1096         ! 
     1097         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 
     1098            e3t_b (i1:i2,j1:j2,1:jpk)  = e3t_n (i1:i2,j1:j2,1:jpk) 
     1099            e3w_b (i1:i2,j1:j2,1:jpk)  = e3w_n (i1:i2,j1:j2,1:jpk) 
     1100            gdepw_b(i1:i2,j1:j2,1:jpk) = gdepw_n(i1:i2,j1:j2,1:jpk) 
     1101            gdept_b(i1:i2,j1:j2,1:jpk) = gdept_n(i1:i2,j1:j2,1:jpk) 
     1102         ENDIF 
     1103         ! 
     1104         DEALLOCATE(ptab) 
     1105      ENDIF 
     1106      ! 
     1107   END SUBROUTINE updatee3t 
     1108 
    6461109#else 
    6471110CONTAINS 
  • branches/2017/dev_METO_MERCATOR_2017/NEMOGCM/NEMO/NST_SRC/agrif_top_sponge.F90

    r6140 r8993  
    4646      ! 
    4747#if defined SPONGE_TOP 
    48       timecoeff = REAL( Agrif_NbStepint(), wp ) / Agrif_rhot() 
     48!!      timecoeff = REAL( Agrif_NbStepint(), wp ) / Agrif_rhot() 
     49!! Assume persistence  
     50      timecoeff = REAL(Agrif_rhot()-1,wp)/REAL(Agrif_rhot()) 
    4951      CALL Agrif_sponge 
    5052      Agrif_SpecialValue    = 0._wp 
     
    7375      ! 
    7476      IF( before ) THEN 
    75          tabres(i1:i2,j1:j2,k1:k2,n1:n2) = trn(i1:i2,j1:j2,k1:k2,n1:n2) 
     77         tabres(i1:i2,j1:j2,k1:k2,n1:n2) = trb(i1:i2,j1:j2,k1:k2,n1:n2) 
    7678      ELSE       
    7779!!gm line below use of :,:  versus i1:i2,j1:j2  ....   strange, not wrong.    ===>> to be corrected 
  • branches/2017/dev_METO_MERCATOR_2017/NEMOGCM/NEMO/NST_SRC/agrif_top_update.F90

    r6140 r8993  
    2424   PUBLIC Agrif_Update_Trc 
    2525 
    26    INTEGER, PUBLIC ::   nbcline_trc = 0   !: ??? 
    27  
    2826   !!---------------------------------------------------------------------- 
    2927   !! NEMO/NST 3.7 , NEMO Consortium (2015) 
     
    3331CONTAINS 
    3432 
    35    SUBROUTINE Agrif_Update_Trc( kt ) 
     33   SUBROUTINE Agrif_Update_Trc( ) 
    3634      !!---------------------------------------------------------------------- 
    3735      !!                   *** ROUTINE Agrif_Update_Trc *** 
    3836      !!---------------------------------------------------------------------- 
    39       INTEGER, INTENT(in) ::   kt 
    40       !!---------------------------------------------------------------------- 
    4137      !  
    42       IF((Agrif_NbStepint() .NE. (Agrif_irhot()-1)).AND.(kt /= 0)) RETURN 
     38      IF (Agrif_Root()) RETURN  
     39      ! 
    4340#if defined TWO_WAY    
    4441      Agrif_UseSpecialValueInUpdate = .TRUE. 
    4542      Agrif_SpecialValueFineGrid    = 0._wp 
    4643      !  
    47       IF( MOD(nbcline_trc,nbclineupdate) == 0 ) THEN 
    4844# if ! defined DECAL_FEEDBACK 
    49          CALL Agrif_Update_Variable(trn_id, procname=updateTRC ) 
     45      CALL Agrif_Update_Variable(trn_id, procname=updateTRC ) 
     46!      CALL Agrif_Update_Variable( trn_id, locupdate=(/0,2/), procname=updateTRC ) 
    5047# else 
    51          CALL Agrif_Update_Variable(trn_id, locupdate=(/1,0/),procname=updateTRC ) 
     48      CALL Agrif_Update_Variable(trn_id, locupdate=(/1,0/),procname=updateTRC ) 
     49!      CALL Agrif_Update_Variable( trn_id, locupdate=(/1,2/), procname=updateTRC ) 
    5250# endif 
    53       ELSE 
    54 # if ! defined DECAL_FEEDBACK 
    55          CALL Agrif_Update_Variable( trn_id, locupdate=(/0,2/), procname=updateTRC ) 
    56 # else 
    57          CALL Agrif_Update_Variable( trn_id, locupdate=(/1,2/), procname=updateTRC ) 
    58 # endif 
    59       ENDIF 
    6051      ! 
    6152      Agrif_UseSpecialValueInUpdate = .FALSE. 
    62       nbcline_trc = nbcline_trc + 1 
     53      ! 
    6354#endif 
    6455      ! 
     
    6657 
    6758 
    68    SUBROUTINE updateTRC( ptab, i1, i2, j1, j2, k1, k2, n1, n2, before ) 
     59   SUBROUTINE updateTRC( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before ) 
    6960      !!---------------------------------------------------------------------- 
    70       !!                      *** ROUTINE updateT *** 
     61      !!                      *** ROUTINE updateTRC *** 
    7162      !!---------------------------------------------------------------------- 
    7263      INTEGER                                    , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2, n1, n2 
    73       REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) ::   ptab 
     64      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) ::   tabres 
    7465      LOGICAL                                    , INTENT(in   ) ::   before 
    7566      !! 
    76       INTEGER ::   ji, jj, jk, jn 
     67      INTEGER :: ji,jj,jk,jn 
     68      REAL(wp) :: ztb, ztnu, ztno 
    7769      !!---------------------------------------------------------------------- 
    7870      ! 
    79       IF( before ) THEN 
    80          ptab(i1:i2,j1:j2,k1:k2,n1:n2) = trn(i1:i2,j1:j2,k1:k2,n1:n2) 
    81       ELSE 
    82          IF( .NOT.(lk_agrif_fstep.AND.(neuler==0)) ) THEN 
    83             ! Add asselin part 
    84             DO jn = n1,n2 
    85                DO jk = k1, k2 
    86                   DO jj = j1, j2 
    87                      DO ji = i1, i2 
    88                         IF( ptab(ji,jj,jk,jn) /= 0._wp ) THEN 
    89                            trb(ji,jj,jk,jn) = trb(ji,jj,jk,jn)             &  
    90                               &             + atfp * ( ptab(ji,jj,jk,jn)   & 
    91                                  &                    - trn(ji,jj,jk,jn) ) * tmask(ji,jj,jk) 
    92                         ENDIF 
    93                      END DO 
     71      ! 
     72      IF (before) THEN 
     73         DO jn = n1,n2 
     74            DO jk=k1,k2 
     75               DO jj=j1,j2 
     76                  DO ji=i1,i2 
     77!> jc tmp 
     78                     tabres(ji,jj,jk,jn) = trn(ji,jj,jk,jn)  * e3t_n(ji,jj,jk) / e3t_0(ji,jj,jk) 
     79!                     tabres(ji,jj,jk,jn) = trn(ji,jj,jk,jn)  * e3t_n(ji,jj,jk) 
     80!< jc tmp 
    9481                  END DO 
    9582               END DO 
    9683            END DO 
     84         END DO 
     85      ELSE 
     86!> jc tmp 
     87         DO jn = n1,n2 
     88            tabres(i1:i2,j1:j2,k1:k2,jn) =  tabres(i1:i2,j1:j2,k1:k2,jn) * e3t_0(i1:i2,j1:j2,k1:k2) & 
     89                                         & * tmask(i1:i2,j1:j2,k1:k2) 
     90         ENDDO 
     91!< jc tmp 
     92         IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN 
     93            ! Add asselin part 
     94            DO jn = n1,n2 
     95               DO jk=k1,k2 
     96                  DO jj=j1,j2 
     97                     DO ji=i1,i2 
     98                        IF( tabres(ji,jj,jk,jn) .NE. 0. ) THEN 
     99                           ztb  = trb(ji,jj,jk,jn) * e3t_b(ji,jj,jk) ! fse3t_b prior update should be used 
     100                           ztnu = tabres(ji,jj,jk,jn) 
     101                           ztno = trn(ji,jj,jk,jn) * e3t_a(ji,jj,jk) 
     102                           trb(ji,jj,jk,jn) = ( ztb + atfp * ( ztnu - ztno) )  &  
     103                                     &        * tmask(ji,jj,jk) / e3t_b(ji,jj,jk) 
     104                        ENDIF 
     105                     ENDDO 
     106                  ENDDO 
     107               ENDDO 
     108            ENDDO 
    97109         ENDIF 
    98          DO jn = n1, n2 
    99             DO jk = k1, k2 
    100                DO jj = j1, j2 
    101                   DO ji = i1, i2 
    102                      IF( ptab(ji,jj,jk,jn) /= 0._wp ) THEN  
    103                         trn(ji,jj,jk,jn) = ptab(ji,jj,jk,jn) * tmask(ji,jj,jk) 
     110         DO jn = n1,n2 
     111            DO jk=k1,k2 
     112               DO jj=j1,j2 
     113                  DO ji=i1,i2 
     114                     IF( tabres(ji,jj,jk,jn) .NE. 0. ) THEN  
     115                        trn(ji,jj,jk,jn) = tabres(ji,jj,jk,jn) / e3t_n(ji,jj,jk) 
    104116                     END IF 
    105117                  END DO 
     
    107119            END DO 
    108120         END DO 
     121         ! 
     122         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 
     123            trb(i1:i2,j1:j2,k1:k2,n1:n2)  = trn(i1:i2,j1:j2,k1:k2,n1:n2) 
     124         ENDIF 
     125         ! 
    109126      ENDIF 
    110127      !  
  • branches/2017/dev_METO_MERCATOR_2017/NEMOGCM/NEMO/NST_SRC/agrif_user.F90

    r7761 r8993  
     1#undef UPD_HIGH   /* MIX HIGH UPDATE */ 
    12#if defined key_agrif 
    23!!---------------------------------------------------------------------- 
     
    8889# endif 
    8990   ! 
     91   IF ( Agrif_Level().EQ.Agrif_MaxLevel() ) CALL agrif_Update_ini() 
     92 
     93   Agrif_UseSpecialValueInUpdate = .FALSE.      
     94 
    9095END SUBROUTINE Agrif_initvalues 
    9196 
     
    144149   CALL Agrif_Set_bc(e2v_id,(/0,0/)) 
    145150 
    146    ! 5. Update type 
     151   ! 4. Update type 
    147152   !---------------  
     153# if defined UPD_HIGH 
     154   CALL Agrif_Set_Updatetype(e1u_id,update1 = Agrif_Update_Average, update2=Agrif_Update_Full_Weighting) 
     155   CALL Agrif_Set_Updatetype(e2v_id,update1 = Agrif_Update_Full_Weighting, update2=Agrif_Update_Average) 
     156#else 
    148157   CALL Agrif_Set_Updatetype(e1u_id,update1 = Agrif_Update_Copy, update2=Agrif_Update_Average) 
    149158   CALL Agrif_Set_Updatetype(e2v_id,update1 = Agrif_Update_Average, update2=Agrif_Update_Copy) 
    150  
    151 ! High order updates 
    152 !   CALL Agrif_Set_Updatetype(e1u_id,update1 = Agrif_Update_Average,            update2=Agrif_Update_Full_Weighting) 
    153 !   CALL Agrif_Set_Updatetype(e2v_id,update1 = Agrif_Update_Full_Weighting,     update2=Agrif_Update_Average) 
    154     ! 
     159#endif 
     160 
    155161END SUBROUTINE agrif_declare_var_dom 
    156162 
     
    175181   ! 
    176182   LOGICAL :: check_namelist 
    177    CHARACTER(len=15) :: cl_check1, cl_check2, cl_check3 
     183   CHARACTER(len=15) :: cl_check1, cl_check2, cl_check3, cl_check4  
    178184   !!---------------------------------------------------------------------- 
    179185 
     
    205211   Agrif_UseSpecialValue = .TRUE. 
    206212   CALL Agrif_Bc_variable(sshn_id,calledweight=1., procname=interpsshn ) 
     213   hbdy_w(:) = 0.e0 ; hbdy_e(:) = 0.e0 ; hbdy_n(:) = 0.e0 ; hbdy_s(:) = 0.e0 
     214   ssha(:,:) = 0.e0 
    207215 
    208216   IF ( ln_dynspg_ts ) THEN 
     
    212220      CALL Agrif_Bc_variable(ub2b_interp_id,calledweight=1.,procname=interpub2b) 
    213221      CALL Agrif_Bc_variable(vb2b_interp_id,calledweight=1.,procname=interpvb2b) 
    214       ubdy_w(:) = 0.e0 ; vbdy_w(:) = 0.e0 ; hbdy_w(:) =0.e0 
    215       ubdy_e(:) = 0.e0 ; vbdy_e(:) = 0.e0 ; hbdy_e(:) =0.e0  
    216       ubdy_n(:) = 0.e0 ; vbdy_n(:) = 0.e0 ; hbdy_n(:) =0.e0  
    217       ubdy_s(:) = 0.e0 ; vbdy_s(:) = 0.e0 ; hbdy_s(:) =0.e0 
     222      ubdy_w(:) = 0.e0 ; vbdy_w(:) = 0.e0 
     223      ubdy_e(:) = 0.e0 ; vbdy_e(:) = 0.e0 
     224      ubdy_n(:) = 0.e0 ; vbdy_n(:) = 0.e0 
     225      ubdy_s(:) = 0.e0 ; vbdy_s(:) = 0.e0 
    218226   ENDIF 
    219227 
     
    234242         WRITE(cl_check2,*)  NINT(rdt) 
    235243         WRITE(cl_check3,*)  NINT(Agrif_Parent(rdt)/Agrif_Rhot()) 
    236          CALL ctl_stop( 'incompatible time step between ocean grids',   & 
     244         CALL ctl_stop( 'Incompatible time step between ocean grids',   & 
    237245               &               'parent grid value : '//cl_check1    ,   &  
    238246               &               'child  grid value : '//cl_check2    ,   &  
     
    245253         WRITE(cl_check1,*)  (Agrif_Parent(nit000)-1)*Agrif_IRhot() + 1 
    246254         WRITE(cl_check2,*)   Agrif_Parent(nitend)   *Agrif_IRhot() 
    247          CALL ctl_warn( 'incompatible run length between grids'               ,   & 
    248                &              ' nit000 on fine grid will be change to : '//cl_check1,   & 
    249                &              ' nitend on fine grid will be change to : '//cl_check2    ) 
     255         CALL ctl_warn( 'Incompatible run length between grids'                      ,   & 
     256               &               'nit000 on fine grid will be changed to : '//cl_check1,   & 
     257               &               'nitend on fine grid will be changed to : '//cl_check2    ) 
    250258         nit000 = (Agrif_Parent(nit000)-1)*Agrif_IRhot() + 1 
    251259         nitend =  Agrif_Parent(nitend)   *Agrif_IRhot() 
    252260      ENDIF 
    253261 
    254       ! Check coordinates 
    255      !SF  IF( ln_zps ) THEN 
    256      !SF     ! check parameters for partial steps  
    257      !SF     IF( Agrif_Parent(e3zps_min) .NE. e3zps_min ) THEN 
    258      !SF        WRITE(*,*) 'incompatible e3zps_min between grids' 
    259      !SF        WRITE(*,*) 'parent grid :',Agrif_Parent(e3zps_min) 
    260      !SF        WRITE(*,*) 'child grid  :',e3zps_min 
    261      !SF        WRITE(*,*) 'those values should be identical' 
    262      !SF        STOP 
    263      !SF     ENDIF 
    264      !SF     IF( Agrif_Parent(e3zps_rat) /= e3zps_rat ) THEN 
    265      !SF        WRITE(*,*) 'incompatible e3zps_rat between grids' 
    266      !SF        WRITE(*,*) 'parent grid :',Agrif_Parent(e3zps_rat) 
    267      !SF        WRITE(*,*) 'child grid  :',e3zps_rat 
    268      !SF        WRITE(*,*) 'those values should be identical'                   
    269      !SF        STOP 
    270      !SF     ENDIF 
    271      !SF  ENDIF 
    272  
    273262      ! Check free surface scheme 
    274263      IF ( ( Agrif_Parent(ln_dynspg_ts ).AND.ln_dynspg_exp ).OR.& 
    275264         & ( Agrif_Parent(ln_dynspg_exp).AND.ln_dynspg_ts ) ) THEN 
    276          WRITE(*,*) 'incompatible free surface scheme between grids' 
    277          WRITE(*,*) 'parent grid ln_dynspg_ts  :', Agrif_Parent(ln_dynspg_ts ) 
    278          WRITE(*,*) 'parent grid ln_dynspg_exp :', Agrif_Parent(ln_dynspg_exp) 
    279          WRITE(*,*) 'child grid  ln_dynspg_ts  :', ln_dynspg_ts 
    280          WRITE(*,*) 'child grid  ln_dynspg_exp :', ln_dynspg_exp 
    281          WRITE(*,*) 'those logicals should be identical'                   
     265         WRITE(cl_check1,*)  Agrif_Parent( ln_dynspg_ts ) 
     266         WRITE(cl_check2,*)  ln_dynspg_ts 
     267         WRITE(cl_check3,*)  Agrif_Parent( ln_dynspg_exp ) 
     268         WRITE(cl_check4,*)  ln_dynspg_exp 
     269         CALL ctl_stop( 'Incompatible free surface scheme between grids' ,  & 
     270               &               'parent grid ln_dynspg_ts  :'//cl_check1  ,  &  
     271               &               'child  grid ln_dynspg_ts  :'//cl_check2  ,  & 
     272               &               'parent grid ln_dynspg_exp :'//cl_check3  ,  & 
     273               &               'child  grid ln_dynspg_exp :'//cl_check4  ,  & 
     274               &               'those logicals should be identical' )                  
     275         STOP 
     276      ENDIF 
     277 
     278      ! Check if identical linear free surface option 
     279      IF ( ( Agrif_Parent(ln_linssh ).AND.(.NOT.ln_linssh )).OR.& 
     280         & ( (.NOT.Agrif_Parent(ln_linssh)).AND.ln_linssh ) ) THEN 
     281         WRITE(cl_check1,*)  Agrif_Parent(ln_linssh ) 
     282         WRITE(cl_check2,*)  ln_linssh 
     283         CALL ctl_stop( 'Incompatible linearized fs option between grids',  & 
     284               &               'parent grid ln_linssh  :'//cl_check1     ,  & 
     285               &               'child  grid ln_linssh  :'//cl_check2     ,  & 
     286               &               'those logicals should be identical' )                   
    282287         STOP 
    283288      ENDIF 
     
    306311   ENDIF 
    307312   !  
    308    ! Do update at initialisation because not done before writing restarts 
    309    ! This would indeed change boundary conditions values at initial time 
    310    ! hence produce restartability issues. 
    311    ! Note that update below is recursive (with lk_agrif_doupd=T): 
    312    !  
    313 ! JC: I am not sure if Agrif_MaxLevel() is the "relative" 
    314 !     or the absolute maximum nesting level...TBC                         
    315    IF ( Agrif_Level().EQ.Agrif_MaxLevel() ) THEN  
    316       ! NB: Do tracers first, dynamics after because nbcline incremented in dynamics 
    317       CALL Agrif_Update_tra() 
    318       CALL Agrif_Update_dyn() 
    319    ENDIF 
    320    ! 
     313END SUBROUTINE Agrif_InitValues_cont 
     314 
     315RECURSIVE SUBROUTINE Agrif_Update_ini( ) 
     316   !!---------------------------------------------------------------------- 
     317   !!                 *** ROUTINE agrif_Update_ini *** 
     318   !! 
     319   !! ** Purpose :: Recursive update done at initialization 
     320   !!---------------------------------------------------------------------- 
     321   USE dom_oce 
     322   USE agrif_opa_update 
     323#if defined key_top 
     324   USE agrif_top_update 
     325#endif 
     326   ! 
     327   IMPLICIT NONE 
     328   !!---------------------------------------------------------------------- 
     329   ! 
     330   IF (Agrif_Root()) RETURN 
     331   ! 
     332                       CALL Agrif_Update_ssh() 
     333   IF (.NOT.ln_linssh) CALL Agrif_Update_vvl() 
     334                       CALL Agrif_Update_tra() 
     335#if defined key_top 
     336                       CALL Agrif_Update_Trc() 
     337#endif 
     338                       CALL Agrif_Update_dyn() 
    321339# if defined key_zdftke 
    322    CALL Agrif_Update_tke(0) 
    323 # endif 
    324    ! 
    325    Agrif_UseSpecialValueInUpdate = .FALSE. 
    326    nbcline = 0 
    327    lk_agrif_doupd = .FALSE. 
    328    ! 
    329 END SUBROUTINE Agrif_InitValues_cont 
    330  
     340! JC remove update because this precludes from perfect restartability 
     341!!                       CALL Agrif_Update_tke(0) 
     342# endif 
     343 
     344   CALL Agrif_ChildGrid_To_ParentGrid() 
     345   CALL Agrif_Update_ini() 
     346   CALL Agrif_ParentGrid_To_ChildGrid() 
     347 
     348END SUBROUTINE agrif_update_ini 
    331349 
    332350SUBROUTINE agrif_declare_var 
     
    371389   CALL agrif_declare_variable((/2,2/),(/3,3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),sshn_id) 
    372390 
    373 # if defined key_zdftke 
    374    CALL agrif_declare_variable((/2,2,0/),(/3,3,0/),(/'x','y','N'/),(/1,1,1/),(/jpi,jpj,jpk/), en_id) 
    375    CALL agrif_declare_variable((/2,2,0/),(/3,3,0/),(/'x','y','N'/),(/1,1,1/),(/jpi,jpj,jpk/),avt_id) 
    376    CALL agrif_declare_variable((/2,2,0/),(/3,3,0/),(/'x','y','N'/),(/1,1,1/),(/jpi,jpj,jpk/),avm_id) 
     391# if defined key_zdftke || defined key_zdfgls 
     392   CALL agrif_declare_variable((/2,2,0/),(/3,3,0/),(/'x','y','N'/),(/1,1,1/),(/nlci,nlcj,jpk/), en_id) 
     393   CALL agrif_declare_variable((/2,2,0/),(/3,3,0/),(/'x','y','N'/),(/1,1,1/),(/nlci,nlcj,jpk/),avt_id) 
     394   CALL agrif_declare_variable((/2,2,0/),(/3,3,0/),(/'x','y','N'/),(/1,1,1/),(/nlci,nlcj,jpk/),avm_id) 
    377395# endif 
    378396 
     
    400418   CALL Agrif_Set_bcinterp(vmsk_id,interp=AGRIF_constant) 
    401419 
    402 # if defined key_zdftke 
     420# if defined key_zdftke || defined key_zdfgls 
    403421   CALL Agrif_Set_bcinterp(avm_id ,interp=AGRIF_linear) 
    404422# endif 
     
    411429   CALL Agrif_Set_bc(vn_interp_id,(/0,1/)) 
    412430 
    413 !   CALL Agrif_Set_bc(tsn_sponge_id,(/-3*Agrif_irhox(),0/)) 
    414 !   CALL Agrif_Set_bc(un_sponge_id,(/-2*Agrif_irhox()-1,0/)) 
    415 !   CALL Agrif_Set_bc(vn_sponge_id,(/-2*Agrif_irhox()-1,0/)) 
    416431   CALL Agrif_Set_bc(tsn_sponge_id,(/-nn_sponge_len*Agrif_irhox()-1,0/)) 
    417432   CALL Agrif_Set_bc(un_sponge_id ,(/-nn_sponge_len*Agrif_irhox()-1,0/)) 
     
    428443   CALL Agrif_Set_bc(vmsk_id,(/0,0/)) 
    429444 
     445# if defined key_zdftke || defined key_zdfgls 
     446   CALL Agrif_Set_bc(avm_id ,(/0,1/)) 
     447# endif 
     448 
     449   ! 4. Update type 
     450   !---------------  
     451   CALL Agrif_Set_Updatetype(scales_t_id, update = AGRIF_Update_Average) 
     452 
     453# if defined UPD_HIGH 
     454   CALL Agrif_Set_Updatetype(tsn_id, update = Agrif_Update_Full_Weighting) 
     455   CALL Agrif_Set_Updatetype(un_update_id,update1 = Agrif_Update_Average, update2 = Agrif_Update_Full_Weighting) 
     456   CALL Agrif_Set_Updatetype(vn_update_id,update1 = Agrif_Update_Full_Weighting, update2 = Agrif_Update_Average) 
     457 
     458   CALL Agrif_Set_Updatetype(ub2b_update_id,update1 = Agrif_Update_Average, update2 = Agrif_Update_Full_Weighting) 
     459   CALL Agrif_Set_Updatetype(vb2b_update_id,update1 = Agrif_Update_Full_Weighting, update2 = Agrif_Update_Average) 
     460   CALL Agrif_Set_Updatetype(sshn_id, update = Agrif_Update_Full_Weighting) 
     461   CALL Agrif_Set_Updatetype(e3t_id, update = Agrif_Update_Full_Weighting) 
     462 
    430463# if defined key_zdftke 
    431    CALL Agrif_Set_bc(avm_id ,(/0,1/)) 
    432 # endif 
    433  
    434    ! 5. Update type 
    435    !---------------  
     464   CALL Agrif_Set_Updatetype( en_id, update = AGRIF_Update_Full_Weighting) 
     465   CALL Agrif_Set_Updatetype(avt_id, update = AGRIF_Update_Full_Weighting) 
     466   CALL Agrif_Set_Updatetype(avm_id, update = AGRIF_Update_Full_Weighting) 
     467# endif 
     468 
     469#else 
    436470   CALL Agrif_Set_Updatetype(tsn_id, update = AGRIF_Update_Average) 
    437  
    438    CALL Agrif_Set_Updatetype(scales_t_id, update = AGRIF_Update_Average) 
    439  
    440471   CALL Agrif_Set_Updatetype(un_update_id,update1 = Agrif_Update_Copy, update2 = Agrif_Update_Average) 
    441472   CALL Agrif_Set_Updatetype(vn_update_id,update1 = Agrif_Update_Average, update2 = Agrif_Update_Copy) 
    442473 
    443    CALL Agrif_Set_Updatetype(sshn_id, update = AGRIF_Update_Average) 
    444  
    445474   CALL Agrif_Set_Updatetype(ub2b_update_id,update1 = Agrif_Update_Copy, update2 = Agrif_Update_Average) 
    446475   CALL Agrif_Set_Updatetype(vb2b_update_id,update1 = Agrif_Update_Average, update2 = Agrif_Update_Copy) 
     476   CALL Agrif_Set_Updatetype(sshn_id, update = AGRIF_Update_Average) 
     477   CALL Agrif_Set_Updatetype(e3t_id, update = AGRIF_Update_Average) 
    447478 
    448479# if defined key_zdftke 
     
    452483# endif 
    453484 
    454 ! High order updates 
    455 !   CALL Agrif_Set_Updatetype(tsn_id, update = Agrif_Update_Full_Weighting) 
    456 !   CALL Agrif_Set_Updatetype(un_update_id,update1 = Agrif_Update_Average, update2 = Agrif_Update_Full_Weighting) 
    457 !   CALL Agrif_Set_Updatetype(vn_update_id,update1 = Agrif_Update_Full_Weighting, update2 = Agrif_Update_Average) 
    458 ! 
    459 !   CALL Agrif_Set_Updatetype(ub2b_update_id,update1 = Agrif_Update_Average, update2 = Agrif_Update_Full_Weighting) 
    460 !   CALL Agrif_Set_Updatetype(vb2b_update_id,update1 = Agrif_Update_Full_Weighting, update2 = Agrif_Update_Average) 
    461 !   CALL Agrif_Set_Updatetype(sshn_id, update = Agrif_Update_Full_Weighting) 
    462  
     485#endif 
    463486   ! 
    464487END SUBROUTINE agrif_declare_var 
     
    594617      CALL ctl_stop('rhot * nn_fsbc(parent) /= N * nn_fsbc(child), therefore nn_fsbc(child) should be set to 1 or nn_fsbc(parent)') 
    595618   ENDIF 
    596  
    597    ! stop if update frequency is different from nn_fsbc 
    598    IF( nbclineupdate > nn_fsbc )  CALL ctl_stop('With ice model on child grid, nn_cln_update should be set to 1 or nn_fsbc') 
    599  
    600  
    601619   ! First Interpolations (using "after" ice subtime step => lim_nbstep=1) 
    602620   !---------------------------------------------------------------------- 
     
    733751      ENDIF 
    734752 
    735       ! Check coordinates 
    736       IF( ln_zps ) THEN 
    737          ! check parameters for partial steps  
    738          IF( Agrif_Parent(e3zps_min) .NE. e3zps_min ) THEN 
    739             WRITE(*,*) 'incompatible e3zps_min between grids' 
    740             WRITE(*,*) 'parent grid :',Agrif_Parent(e3zps_min) 
    741             WRITE(*,*) 'child grid  :',e3zps_min 
    742             WRITE(*,*) 'those values should be identical' 
    743             STOP 
    744          ENDIF 
    745          IF( Agrif_Parent(e3zps_rat) .NE. e3zps_rat ) THEN 
    746             WRITE(*,*) 'incompatible e3zps_rat between grids' 
    747             WRITE(*,*) 'parent grid :',Agrif_Parent(e3zps_rat) 
    748             WRITE(*,*) 'child grid  :',e3zps_rat 
    749             WRITE(*,*) 'those values should be identical'                   
    750             STOP 
    751          ENDIF 
    752753      ENDIF 
    753754      ! Check passive tracer cell 
     
    756757      ENDIF 
    757758   ENDIF 
    758  
    759    CALL Agrif_Update_trc(0) 
    760    ! 
    761    Agrif_UseSpecialValueInUpdate = .FALSE. 
    762    nbcline_trc = 0 
    763759   ! 
    764760END SUBROUTINE Agrif_InitValues_cont_top 
     
    792788   !----------------------------- 
    793789   CALL Agrif_Set_bc(trn_id,(/0,1/)) 
    794 !   CALL Agrif_Set_bc(trn_sponge_id,(/-3*Agrif_irhox(),0/)) 
    795790   CALL Agrif_Set_bc(trn_sponge_id,(/-nn_sponge_len*Agrif_irhox()-1,0/)) 
    796791 
    797    ! 5. Update type 
     792   ! 4. Update type 
    798793   !---------------  
     794# if defined UPD_HIGH 
     795   CALL Agrif_Set_Updatetype(trn_id, update = Agrif_Update_Full_Weighting) 
     796#else 
    799797   CALL Agrif_Set_Updatetype(trn_id, update = AGRIF_Update_Average) 
    800  
    801 !   Higher order update 
    802 !   CALL Agrif_Set_Updatetype(tsn_id, update = Agrif_Update_Full_Weighting) 
    803  
     798#endif 
    804799   ! 
    805800END SUBROUTINE agrif_declare_var_top 
     
    832827   INTEGER  ::   ios                 ! Local integer output status for namelist read 
    833828   INTEGER  ::   iminspon 
    834    NAMELIST/namagrif/ nn_cln_update, rn_sponge_tra, rn_sponge_dyn, ln_spc_dyn, ln_chk_bathy 
     829   NAMELIST/namagrif/ rn_sponge_tra, rn_sponge_dyn, ln_spc_dyn, ln_chk_bathy 
    835830   !!-------------------------------------------------------------------------------------- 
    836831   ! 
     
    849844      WRITE(numout,*) '~~~~~~~~~~~~~~~' 
    850845      WRITE(numout,*) '   Namelist namagrif : set AGRIF parameters' 
    851       WRITE(numout,*) '      baroclinic update frequency       nn_cln_update = ', nn_cln_update 
    852846      WRITE(numout,*) '      sponge coefficient for tracers    rn_sponge_tra = ', rn_sponge_tra, ' s' 
    853847      WRITE(numout,*) '      sponge coefficient for dynamics   rn_sponge_tra = ', rn_sponge_dyn, ' s' 
     
    858852   ! 
    859853   ! convert DOCTOR namelist name into OLD names 
    860    nbclineupdate = nn_cln_update 
    861854   visc_tra      = rn_sponge_tra 
    862855   visc_dyn      = rn_sponge_dyn 
     
    878871SUBROUTINE Agrif_InvLoc( indloc, nprocloc, i, indglob ) 
    879872   !!---------------------------------------------------------------------- 
    880    !!                     *** ROUTINE Agrif_detect *** 
     873   !!                     *** ROUTINE Agrif_InvLoc *** 
    881874   !!---------------------------------------------------------------------- 
    882875   USE dom_oce 
Note: See TracChangeset for help on using the changeset viewer.