Changeset 13603
- Timestamp:
- 2020-10-14T18:04:16+02:00 (4 years ago)
- File:
-
- 1 moved
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/dev_r13327_KERNEL-06_2_techene_e3/src/SWE/stpMLF.F90
r13602 r13603 1 MODULE stp LF1 MODULE stpMLF 2 2 !!====================================================================== 3 !! *** MODULE st ep***3 !! *** MODULE stpMLF *** 4 4 !! Time-stepping : manager of the shallow water equation time stepping 5 !! Modified Leap Frog scheme 5 6 !!====================================================================== 6 7 !! History : NEMO ! 2020-03 (A. Nasser, G. Madec) Original code from 4.0.2 7 !!---------------------------------------------------------------------- 8 9 !!---------------------------------------------------------------------- 10 !! stp : Shallow Water time-stepping 11 !!---------------------------------------------------------------------- 12 USE step_oce ! time stepping definition modules 13 USE phycst ! physical constants 14 USE usrdef_nam 8 !! - ! 2020-10 (S. Techene, G. Madec) cleanning 9 !!---------------------------------------------------------------------- 10 11 !!---------------------------------------------------------------------- 12 !! stp_MLF : MLF Shallow Water Eq. time-stepping 13 !!---------------------------------------------------------------------- 14 USE stp_oce ! modules used in nemo_init and stp_MLF 15 15 ! 16 USE iom ! xIOs server 17 USE domqco 18 16 USE domqco ! quasi-eulerian coordinate 17 USE phycst ! physical constants 18 USE usrdef_nam ! user defined namelist parameters 19 !!st USE usrdef_sbc ! user defined surface boundary cond 20 19 21 IMPLICIT NONE 20 22 PRIVATE 21 23 22 PUBLIC stp_ LF ! called by nemogcm.F9024 PUBLIC stp_MLF ! called by nemogcm.F90 23 25 24 !!---------------------------------------------------------------------- 25 !! time level indices 26 !!---------------------------------------------------------------------- 27 INTEGER, PUBLIC :: Nbb, Nnn, Naa, Nrhs !! used by nemo_init 26 ! !** time level indices **! 27 INTEGER, PUBLIC :: Nbb, Nnn, Naa, Nrhs !: used by nemo_init 28 28 29 29 !! * Substitutions … … 37 37 CONTAINS 38 38 39 #if defined key_agrif 40 RECURSIVE SUBROUTINE stp_LF( ) 41 INTEGER :: kstp ! ocean time-step index 42 #else 43 SUBROUTINE stp_LF( kstp ) 44 INTEGER, INTENT(in) :: kstp ! ocean time-step index 45 #endif 39 SUBROUTINE stp_MLF( kstp ) 46 40 !!---------------------------------------------------------------------- 47 !! *** ROUTINE stp ***41 !! *** ROUTINE stp_MLF *** 48 42 !! 49 43 !! ** Purpose : - Time stepping of shallow water (SHW) (momentum and ssh eqs.) … … 56 50 !! -6- Outputs and diagnostics 57 51 !!---------------------------------------------------------------------- 58 INTEGER :: ji, jj, jk ! dummy loop indice 59 INTEGER :: indic ! error indicator if < 0 60 !!gm kcall can be removed, I guess 61 INTEGER :: kcall ! optional integer argument (dom_vvl_sf_nxt) 62 REAL(wp):: z1_2rho0 ! local scalars 63 64 REAL(wp) :: zue3a, zue3n, zue3b ! local scalars 65 REAL(wp) :: zve3a, zve3n, zve3b ! - - 66 REAL(wp) :: ze3t_tf, ze3u_tf, ze3v_tf, zua, zva 52 INTEGER, INTENT(in) :: kstp ! ocean time-step index 53 ! 54 INTEGER :: ji, jj, jk ! dummy loop indices 55 INTEGER :: indic ! error indicator if < 0 56 REAL(wp):: z1_2rho0 ! local scalars 57 REAL(wp):: zrhs_u, zue3a, zue3n, zue3b, zua ! local scalars 58 REAL(wp):: zrhs_v, zve3a, zve3n, zve3b, zva ! - - 67 59 !! --------------------------------------------------------------------- 68 #if defined key_agrif 69 kstp = nit000 + Agrif_Nb_Step() 70 Kbb_a = Nbb; Kmm_a = Nnn; Krhs_a = Nrhs ! agrif_oce module copies of time level indices 71 IF( lk_agrif_debug ) THEN 72 IF( Agrif_Root() .and. lwp) WRITE(*,*) '---' 73 IF(lwp) WRITE(*,*) 'Grid Number', Agrif_Fixed(),' time step ', kstp, 'int tstep', Agrif_NbStepint() 74 ENDIF 75 IF( kstp == nit000 + 1 ) lk_agrif_fstep = .FALSE. 76 # if defined key_iomput 77 IF( Agrif_Nbstepint() == 0 ) CALL iom_swap( cxios_context ) 78 # endif 79 #endif 80 ! 81 IF( ln_timing ) CALL timing_start('stp_LF') 60 ! 61 IF( ln_timing ) CALL timing_start('stp_MLF') 82 62 ! 83 63 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> … … 85 65 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 86 66 ! 87 IF( l_1st_euler ) THEN 88 rDt = 67 IF( l_1st_euler ) THEN ! start or restart with Euler 1st time-step 68 rDt = rn_Dt 89 69 r1_Dt = 1._wp / rDt 90 70 ENDIF 91 71 92 IF 72 IF( kstp == nit000 ) ww(:,:,:) = 0._wp ! initialize vertical velocity one for all to zero 93 73 94 74 ! … … 97 77 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 98 78 indic = 0 ! reset to no error condition 99 79 100 80 IF( kstp == nit000 ) THEN ! initialize IOM context (must be done after nemo_init for AGRIF+XIOS+OASIS) 101 81 CALL iom_init( cxios_context, ld_closedef=.FALSE. ) ! for model grid (including passible AGRIF zoom) 102 IF( lk_diamlr ) CALL dia_mlr_iom_init ! with additional setup for multiple-linear-regression analysis103 82 CALL iom_init_closedef 104 IF( ln_crs ) CALL iom_init( TRIM(cxios_context)//"_crs" ) ! for coarse grid105 83 ENDIF 106 84 IF( kstp /= nit000 ) CALL day( kstp ) ! Calendar (day was already called at nit000 in day_init) 107 85 CALL iom_setkt( kstp - nit000 + 1, cxios_context ) ! tell IOM we are at time step kstp 108 IF( ln_crs ) CALL iom_setkt( kstp - nit000 + 1, TRIM(cxios_context)//"_crs" ) ! tell IOM we are at time step kstp 109 110 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 111 ! Update external forcing (tides, open boundaries, ice shelf interaction and surface boundary condition (including sea-ice) 112 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 113 IF( ln_tide ) CALL tide_update( kstp ) ! update tide potential 114 IF( ln_apr_dyn ) CALL sbc_apr ( kstp ) ! atmospheric pressure (NB: call before bdy_dta which needs ssh_ib) 115 IF( ln_bdy ) CALL bdy_dta ( kstp, Nnn ) ! update dynamic & tracer data at open boundaries 116 CALL sbc ( kstp, Nbb, Nnn ) ! Sea Boundary Condition 117 118 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 119 ! Ocean physics update 120 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 121 ! LATERAL PHYSICS 122 ! ! eddy diffusivity coeff. 86 87 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 88 ! Update external forcing (SWE: surface boundary condition only) 89 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 90 91 CALL sbc ( kstp, Nbb, Nnn ) ! Sea Boundary Condition 92 93 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 94 ! Ocean physics update (SWE: eddy viscosity only) 95 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 96 123 97 IF( l_ldfdyn_time ) CALL ldf_dyn( kstp, Nbb ) ! eddy viscosity coeff. 124 98 125 99 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 126 ! Ocean dynamics : hdiv, ssh, e3, u, v, w 127 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 128 129 CALL ssh_nxt ( kstp, Nbb, Nnn, ssh, Naa ) ! after ssh (includes call to div_hor) 100 ! RHS of horizontal velocity 101 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 102 130 103 uu(:,:,:,Nrhs) = 0._wp ! set dynamics trends to zero 131 104 vv(:,:,:,Nrhs) = 0._wp 132 105 133 IF( ln_bdy ) CALL bdy_dyn3d_dmp ( kstp, Nbb, uu, vv, Nrhs ) ! bdy damping trends 134 135 #if defined key_agrif 136 IF(.NOT. Agrif_Root()) & 137 & CALL Agrif_Sponge_dyn ! momentum sponge 138 #endif 139 CALL dyn_adv( kstp, Nbb, Nnn , uu, vv, Nrhs ) ! advection (VF or FF) ==> RHS 140 141 CALL dyn_vor( kstp, Nnn , uu, vv, Nrhs ) ! vorticity ==> RHS 142 143 CALL dyn_ldf( kstp, Nbb, Nnn , uu, vv, Nrhs ) ! lateral mixing 144 145 IF( .NOT.ln_linssh ) CALL dom_qco_r3c ( ssh(:,:,Naa), r3t(:,:,Naa), r3u(:,:,Naa), r3v(:,:,Naa), r3f(:,:) ) ! "after" ssh./h._0 ratio explicit 146 !IF( .NOT.ln_linssh ) CALL dom_vvl_sf_nxt_st( kstp, Nbb, Nnn, Naa ) ! after vertical scale factors 147 !!an - calcul du gradient de pression horizontal (explicit) 106 CALL dyn_adv( kstp, Nbb, Nnn, uu, vv, Nrhs ) ! advection (VF or FF) ==> RHS 107 CALL dyn_vor( kstp, Nnn, uu, vv, Nrhs ) ! vorticity ==> RHS 108 CALL dyn_ldf( kstp, Nbb, Nnn, uu, vv, Nrhs ) ! lateral mixing ==> RHS 109 110 z1_2rho0 = 0.5_wp * r1_rho0 111 ! 148 112 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 149 uu(ji,jj,jk,Nrhs) = uu(ji,jj,jk,Nrhs) - grav * ( ssh(ji+1,jj,Nnn) - ssh(ji,jj,Nnn) ) * r1_e1u(ji,jj) 150 vv(ji,jj,jk,Nrhs) = vv(ji,jj,jk,Nrhs) - grav * ( ssh(ji,jj+1,Nnn) - ssh(ji,jj,Nnn) ) * r1_e2v(ji,jj) 113 ! ! horizontal pressure gradient 114 zrhs_u = - grav * ( ssh(ji+1,jj,Nnn) - ssh(ji,jj,Nnn) ) * r1_e1u(ji,jj) 115 zrhs_v = - grav * ( ssh(ji,jj+1,Nnn) - ssh(ji,jj,Nnn) ) * r1_e2v(ji,jj) 116 ! ! wind stress and layer friction 117 zrhs_u = zrhs_u + z1_2rho0 * ( utau_b(ji,jj) + utau(ji,jj) ) / e3u(ji,jj,jk,Nnn) & 118 & - rn_rfr * uu(ji,jj,jk,Nbb) 119 zrhs_v = zrhs_v + z1_2rho0 * ( vtau_b(ji,jj) + vtau(ji,jj) ) / e3v(ji,jj,jk,Nnn) & 120 & - rn_rfr * vv(ji,jj,jk,Nbb) 121 ! ! ==> RHS 122 uu(ji,jj,jk,Nrhs) = uu(ji,jj,jk,Nrhs) + zrhs_u 123 vv(ji,jj,jk,Nrhs) = vv(ji,jj,jk,Nrhs) + zrhs_v 151 124 END_3D 152 ! 153 ! add wind stress forcing and layer linear friction to the RHS 154 z1_2rho0 = 0.5_wp * r1_rho0 155 DO_3D( 0, 0, 0, 0,1,jpkm1) 156 uu(ji,jj,jk,Nrhs) = uu(ji,jj,jk,Nrhs) + z1_2rho0 * ( utau_b(ji,jj) + utau(ji,jj) ) / e3u(ji,jj,jk,Nnn) & 157 & - rn_rfr * uu(ji,jj,jk,Nbb) 158 vv(ji,jj,jk,Nrhs) = vv(ji,jj,jk,Nrhs) + z1_2rho0 * ( vtau_b(ji,jj) + vtau(ji,jj) ) / e3v(ji,jj,jk,Nnn) & 159 & - rn_rfr * vv(ji,jj,jk,Nbb) 160 END_3D 161 !!an 162 163 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 164 ! Leap-Frog time splitting + Robert-Asselin time filter on u,v,e3 165 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 166 !!st ssh_atf : add ssh filtering up there 167 IF ( .NOT.( l_1st_euler ) ) THEN ! Only do time filtering for leapfrog timesteps 168 ! ! filtering "now" field 125 126 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 127 ! Time stepping of ssh Eq. (and update r3_Naa) 128 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 129 ! ! Leap Frog time stepping ==> ssh_Naa and r3_Naa 130 CALL ssh_nxt( kstp, Nbb, Nnn, ssh , Naa ) ! after ssh 131 ! ! after ssh/h_0 ratio explicit 132 CALL dom_qco_r3c( ssh(:,:,Naa), r3t(:,:,Naa), r3u(:,:,Naa), r3v(:,:,Naa), r3f(:,:) ) 133 ! ! Asselin filter ==> ssh_Nnn filtered 134 IF ( .NOT.( l_1st_euler ) ) THEN ! Time filtering of now ssh 169 135 ssh(:,:,Nnn) = ssh(:,:,Nnn) + rn_atfp * ( ssh(:,:,Nbb) - 2._wp * ssh(:,:,Nnn) + ssh(:,:,Naa) ) 170 136 ENDIF 171 !!st ssh_atf 172 173 !! what about IF( .NOT.ln_linssh ) ? 174 !!an futur module dyn_nxt (a la place de dyn_atf) 175 137 138 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 139 ! Time stepping of dynamics (u,v) 140 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 141 176 142 IF( ln_dynadv_vec ) THEN ! vector invariant form : applied on velocity 177 IF( l_1st_euler ) THEN 178 DO_3D( 0, 0, 0, 0,1,jpkm1)143 IF( l_1st_euler ) THEN ! Euler time stepping (no Asselin filter) 144 DO_3D( 0,0, 0,0, 1,jpkm1) 179 145 uu(ji,jj,jk,Naa) = uu(ji,jj,jk,Nbb) + rDt * uu(ji,jj,jk,Nrhs) * umask(ji,jj,jk) 180 146 vv(ji,jj,jk,Naa) = vv(ji,jj,jk,Nbb) + rDt * vv(ji,jj,jk,Nrhs) * vmask(ji,jj,jk) 181 147 END_3D 182 ELSE 183 DO_3D( 1, 1, 1, 1,1,jpkm1)148 ELSE ! Leap Frog time stepping + Asselin filter 149 DO_3D( 0,0, 0,0, 1,jpkm1) 184 150 zua = uu(ji,jj,jk,Nbb) + rDt * uu(ji,jj,jk,Nrhs) * umask(ji,jj,jk) 185 151 zva = vv(ji,jj,jk,Nbb) + rDt * vv(ji,jj,jk,Nrhs) * vmask(ji,jj,jk) … … 191 157 vv(ji,jj,jk,Naa) = zva 192 158 END_3D 193 CALL dom_qco_r3c( ssh(:,:,Nnn), r3t(:,:,Nnn), r3u(:,:,Nnn), r3v(:,:,Nnn) ) ! "now" ssh/h_0 ratio from filtrered ssh 194 #if ! defined key_qco 195 DO jk = 1, jpkm1 196 e3t(:,:,jk,Nnn) = e3t_0(:,:,jk) * ( 1._wp + r3t(:,:,Nnn) ) 197 e3u(:,:,jk,Nnn) = e3u_0(:,:,jk) * ( 1._wp + r3u(:,:,Nnn) ) 198 e3v(:,:,jk,Nnn) = e3v_0(:,:,jk) * ( 1._wp + r3v(:,:,Nnn) ) 199 END DO 200 #endif 159 ! ! Update r3_Nnn 160 CALL dom_qco_r3c( ssh(:,:,Nnn), r3t(:,:,Nnn), r3u(:,:,Nnn), r3v(:,:,Nnn) ) ! now ssh/h_0 ratio from filtrered ssh 201 161 ENDIF 202 162 ! 203 163 ELSE ! flux form : applied on thickness weighted velocity 204 IF( l_1st_euler ) THEN 205 DO_3D( 0, 0, 0, 0,1,jpkm1)164 IF( l_1st_euler ) THEN ! Euler time stepping (no Asselin filter) 165 DO_3D( 0,0, 0,0, 1,jpkm1) 206 166 zue3b = e3u(ji,jj,jk,Nbb) * uu(ji,jj,jk,Nbb) 207 167 zve3b = e3v(ji,jj,jk,Nbb) * vv(ji,jj,jk,Nbb) 208 ! ! LFtime stepping209 zue3a = zue3b + rDt * e3u(ji,jj,jk,N rhs) * uu(ji,jj,jk,Nrhs) * umask(ji,jj,jk)210 zve3a = zve3b + rDt * e3v(ji,jj,jk,N rhs) * vv(ji,jj,jk,Nrhs) * vmask(ji,jj,jk)168 ! ! Euler time stepping 169 zue3a = zue3b + rDt * e3u(ji,jj,jk,Nnn) * uu(ji,jj,jk,Nrhs) * umask(ji,jj,jk) 170 zve3a = zve3b + rDt * e3v(ji,jj,jk,Nnn) * vv(ji,jj,jk,Nrhs) * vmask(ji,jj,jk) 211 171 ! 212 172 uu(ji,jj,jk,Naa) = zue3a / e3u(ji,jj,jk,Naa) … … 214 174 END_3D 215 175 ELSE ! Leap Frog time stepping + Asselin filter 216 CALL dom_qco_r3c( ssh(:,:,Nnn), r3t_f(:,:), r3u_f(:,:), r3v_f(:,:) ) ! "now"ssh/h_0 ratio from filtrered ssh217 DO_3D( 1, 1, 1, 1,1,jpkm1)218 zue3n = e3u(ji,jj,jk,Nnn) * uu(ji,jj,jk,Nnn)219 zve3n = e3v(ji,jj,jk,Nnn) * vv(ji,jj,jk,Nnn)220 zue3b = e3u(ji,jj,jk,Nbb) * uu(ji,jj,jk,Nbb)221 zve3b = e3v(ji,jj,jk,Nbb) * vv(ji,jj,jk,Nbb)176 CALL dom_qco_r3c( ssh(:,:,Nnn), r3t_f(:,:), r3u_f(:,:), r3v_f(:,:) ) ! now ssh/h_0 ratio from filtrered ssh 177 DO_3D( 0,0, 0,0, 1,jpkm1) 178 zue3n = ( 1._wp + r3u(ji,jj,Nnn) ) * uu(ji,jj,jk,Nnn) 179 zve3n = ( 1._wp + r3v(ji,jj,Nnn) ) * vv(ji,jj,jk,Nnn) 180 zue3b = ( 1._wp + r3u(ji,jj,Nbb) ) * uu(ji,jj,jk,Nbb) 181 zve3b = ( 1._wp + r3v(ji,jj,Nbb) ) * vv(ji,jj,jk,Nbb) 222 182 ! ! LF time stepping 223 zue3a = zue3b + rDt * e3u(ji,jj,jk,Nrhs) * uu(ji,jj,jk,Nrhs) * umask(ji,jj,jk) 224 zve3a = zve3b + rDt * e3v(ji,jj,jk,Nrhs) * vv(ji,jj,jk,Nrhs) * vmask(ji,jj,jk) 225 ! ! Asselin time filter on e3u/v/t 226 ze3u_tf = e3u(ji,jj,jk,Nnn) + rn_atfp * ( e3u(ji,jj,jk,Nbb) - 2._wp * e3u(ji,jj,jk,Nnn) + e3u(ji,jj,jk,Naa) ) 227 ze3v_tf = e3v(ji,jj,jk,Nnn) + rn_atfp * ( e3v(ji,jj,jk,Nbb) - 2._wp * e3v(ji,jj,jk,Nnn) + e3v(ji,jj,jk,Naa) ) 228 ze3t_tf = e3t(ji,jj,jk,Nnn) + rn_atfp * ( e3t(ji,jj,jk,Nbb) - 2._wp * e3t(ji,jj,jk,Nnn) + e3t(ji,jj,jk,Naa) ) 183 zue3a = zue3b + rDt * ( 1._wp + r3u(ji,jj,Nnn) ) * uu(ji,jj,jk,Nrhs) * umask(ji,jj,jk) 184 zve3a = zve3b + rDt * ( 1._wp + r3v(ji,jj,Nnn) ) * vv(ji,jj,jk,Nrhs) * vmask(ji,jj,jk) 229 185 ! ! Asselin time filter on u,v (Nnn) 230 uu(ji,jj,jk,Nnn) = ( zue3n + rn_atfp * ( zue3b - 2._wp * zue3n + zue3a ) ) / ( e3u_0(ji,jj,jk) * ( 1._wp + r3u_f(ji,jj) ))231 vv(ji,jj,jk,Nnn) = ( zve3n + rn_atfp * ( zve3b - 2._wp * zve3n + zve3a ) ) / ( e3v_0(ji,jj,jk) * ( 1._wp + r3v_f(ji,jj) ))186 uu(ji,jj,jk,Nnn) = ( zue3n + rn_atfp * ( zue3b - 2._wp * zue3n + zue3a ) ) / ( 1._wp + r3u_f(ji,jj) ) 187 vv(ji,jj,jk,Nnn) = ( zve3n + rn_atfp * ( zve3b - 2._wp * zve3n + zve3a ) ) / ( 1._wp + r3v_f(ji,jj) ) 232 188 ! 233 uu(ji,jj,jk,Naa) = zue3a / e3u(ji,jj,jk,Naa)234 vv(ji,jj,jk,Naa) = zve3a / e3v(ji,jj,jk,Naa)189 uu(ji,jj,jk,Naa) = zue3a / ( 1._wp + r3u(ji,jj,Naa) ) 190 vv(ji,jj,jk,Naa) = zve3a / ( 1._wp + r3v(ji,jj,Naa) ) 235 191 END_3D 192 ! ! Update r3_Nnn with time filtered values 236 193 r3t(:,:,Nnn) = r3t_f(:,:) 237 194 r3u(:,:,Nnn) = r3u_f(:,:) 238 195 r3v(:,:,Nnn) = r3v_f(:,:) 239 #if ! defined key_qco240 DO jk = 1, jpkm1241 e3t(:,:,jk,Nnn) = e3t_0(:,:,jk) * ( 1._wp + r3t(:,:,Nnn) )242 e3u(:,:,jk,Nnn) = e3u_0(:,:,jk) * ( 1._wp + r3u(:,:,Nnn) )243 e3v(:,:,jk,Nnn) = e3v_0(:,:,jk) * ( 1._wp + r3v(:,:,Nnn) )244 END DO245 #endif246 196 ENDIF 247 197 ENDIF 248 198 249 250 CALL lbc_lnk_multi( 'stp', uu(:,:,:,Nnn), 'U', -1., vv(:,:,:,Nnn), 'V', -1., & !* local domain boundaries 251 & uu(:,:,:,Naa), 'U', -1., vv(:,:,:,Naa), 'V', -1. ) 252 253 !!an 199 CALL lbc_lnk_multi( 'stp_MLF', uu(:,:,:,Nnn), 'U', -1., vv(:,:,:,Nnn), 'V', -1., & !* local domain boundaries 200 & uu(:,:,:,Naa), 'U', -1., vv(:,:,:,Naa), 'V', -1. ) 254 201 255 202 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 256 203 ! Set boundary conditions, time filter and swap time levels 257 204 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 258 259 !!an TO BE ADDED : dyn_nxt260 !! CALL dyn_atf ( kstp, Nbb, Nnn, Naa, uu, vv, e3t, e3u, e3v ) ! time filtering of "now" velocities and scale factors261 !!an TO BE ADDED : a simplifier262 ! CALL ssh_atf ( kstp, Nbb, Nnn, Naa, ssh ) ! time filtering of "now" sea surface height263 !!st copied above264 !! IF ( .NOT.( l_1st_euler ) ) THEN ! Only do time filtering for leapfrog timesteps265 !! ! ! filtering "now" field266 !! ssh(:,:,Nnn) = ssh(:,:,Nnn) + rn_atfp * ( ssh(:,:,Nbb) - 2 * ssh(:,:,Nnn) + ssh(:,:,Naa) )267 !! ENDIF268 !!st269 !!an270 271 205 272 206 ! Swap time levels … … 275 209 Nnn = Naa 276 210 Naa = Nrhs 277 ! 278 ! CALL dom_vvl_sf_update_st( kstp, Nbb, Nnn, Naa ) ! recompute vertical scale factors 211 279 212 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 280 213 ! diagnostics and outputs 281 214 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 282 IF( ln_floats ) CALL flo_stp ( kstp, Nbb, Nnn ) ! drifting Floats 283 IF( ln_diacfl ) CALL dia_cfl ( kstp, Nnn ) ! Courant number diagnostics 284 285 CALL dia_wri ( kstp, Nnn ) ! ocean model: outputs 286 287 ! 288 IF( lrst_oce ) CALL rst_write ( kstp, Nbb, Nnn ) ! write output ocean restart file 289 290 291 #if defined key_agrif 292 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 293 ! AGRIF 294 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 295 Kbb_a = Nbb; Kmm_a = Nnn; Krhs_a = Nrhs ! agrif_oce module copies of time level indices 296 CALL Agrif_Integrate_ChildGrids( stp_LF ) ! allows to finish all the Child Grids before updating 297 298 IF( Agrif_NbStepint() == 0 ) THEN 299 CALL Agrif_update_all( ) ! Update all components 300 ENDIF 301 #endif 302 IF( ln_diaobs ) CALL dia_obs ( kstp, Nnn ) ! obs-minus-model (assimilation) diagnostics (call after dynamics update) 215 216 IF( ln_diacfl ) CALL dia_cfl ( kstp, Nnn ) ! Courant number diagnostics 217 CALL dia_wri ( kstp, Nnn ) ! ocean model: outputs 218 ! 219 IF( lrst_oce ) CALL rst_write( kstp, Nbb, Nnn ) ! write output ocean restart file 303 220 304 221 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 305 222 ! Control 306 223 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 307 CALL stp_ctl ( kstp, Nbb, Nnn, indic ) 308 309 224 CALL stp_ctl ( kstp, Nnn ) 225 310 226 IF( kstp == nit000 ) THEN ! 1st time step only 311 227 CALL iom_close( numror ) ! close input ocean restart file … … 316 232 ! 317 233 #if defined key_iomput 318 IF( kstp == nitend .OR. indic < 0 ) THEN 319 CALL iom_context_finalize( cxios_context ) ! needed for XIOS+AGRIF 320 IF(lrxios) CALL iom_context_finalize( crxios_context ) 234 IF( kstp == nitend .OR. indic < 0 ) THEN 235 !!st : cxios_context needed ? because opened earlier ??? 236 CALL iom_context_finalize( cxios_context ) ! needed for XIOS+AGRIF 237 !!st : crxios_context not needed associated with coarsening ! 238 IF(lrxios) CALL iom_context_finalize( crxios_context ) 321 239 ENDIF 322 240 #endif 323 241 ! 324 242 IF( l_1st_euler ) THEN ! recover Leap-frog timestep 325 rDt = 2._wp * rn_Dt243 rDt = 2._wp * rn_Dt 326 244 r1_Dt = 1._wp / rDt 327 245 l_1st_euler = .FALSE. 328 246 ENDIF 329 247 ! 330 IF( ln_timing ) CALL timing_stop('stp ')331 ! 332 END SUBROUTINE stp_ LF333 ! 248 IF( ln_timing ) CALL timing_stop('stp_MLF') 249 ! 250 END SUBROUTINE stp_MLF 251 334 252 !!====================================================================== 335 END MODULE stp LF253 END MODULE stpMLF
Note: See TracChangeset
for help on using the changeset viewer.
