- Timestamp:
- 2018-06-30T12:51:02+02:00 (6 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2018/dev_r9838_ENHANCE04_MLF/src/OCE/DYN/dynnxt.F90
r9598 r9863 64 64 CONTAINS 65 65 66 SUBROUTINE dyn_nxt 66 SUBROUTINE dyn_nxt( kt ) 67 67 !!---------------------------------------------------------------------- 68 68 !! *** ROUTINE dyn_nxt *** … … 92 92 !! un,vn now horizontal velocity of next time-step 93 93 !!---------------------------------------------------------------------- 94 INTEGER, INTENT( in ) :: kt 94 INTEGER, INTENT( in ) :: kt ! ocean time-step index 95 95 ! 96 96 INTEGER :: ji, jj, jk ! dummy loop indices 97 97 INTEGER :: ikt ! local integers 98 REAL(wp) :: zue3a, zue3n, zue3b, zuf, zcoef 99 REAL(wp) :: zve3a, zve3n, zve3b, zvf , z1_2dt! - -98 REAL(wp) :: zue3a, zue3n, zue3b, zuf, zcoef ! local scalars 99 REAL(wp) :: zve3a, zve3n, zve3b, zvf ! - - 100 100 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zue, zve 101 101 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ze3u_f, ze3v_f, zua, zva … … 132 132 ! so that asselin contribution is removed at the same time 133 133 DO jk = 1, jpkm1 134 un(:,:,jk) = ( un(:,:,jk) - un_adv(:,:)*r1_hu_n(:,:) + un_b(:,:) ) *umask(:,:,jk)135 vn(:,:,jk) = ( vn(:,:,jk) - vn_adv(:,:)*r1_hv_n(:,:) + vn_b(:,:) ) *vmask(:,:,jk)134 un(:,:,jk) = ( un(:,:,jk) - un_adv(:,:)*r1_hu_n(:,:) + un_b(:,:) ) * umask(:,:,jk) 135 vn(:,:,jk) = ( vn(:,:,jk) - vn_adv(:,:)*r1_hv_n(:,:) + vn_b(:,:) ) * vmask(:,:,jk) 136 136 END DO 137 137 ENDIF … … 152 152 !!$ Do we need a call to bdy_vol here?? 153 153 ! 154 IF( l_trddyn ) THEN ! prepare the atf trend computation + some diagnostics155 z1_2dt = 1._wp / (2. * rdt) ! Euler or leap-frog time step156 IF( neuler == 0 .AND. kt == nit000 ) z1_2dt = 1._wp / rdt157 !158 ! ! Kinetic energy and Conversion159 IF( ln_KE_trd ) CALL trd_dyn( ua, va, jpdyn_ken, kt )160 !161 IF( ln_dyn_trd ) THEN ! 3D output: total momentum trends162 zua(:,:,:) = ( ua(:,:,:) - ub(:,:,:) ) * z1_2dt163 zva(:,:,:) = ( va(:,:,:) - vb(:,:,:) ) * z1_2dt164 CALL iom_put( "utrd_tot", zua ) ! total momentum trends, except the asselin timefilter154 IF( l_trddyn ) THEN !* prepare the atf trend computation + some diagnostics 155 IF( ln_KE_trd ) CALL trd_dyn( ua, va, jpdyn_ken, kt ) ! Kinetic energy and Conversion 156 IF( ln_dyn_trd ) THEN ! 3D output: total momentum trends 157 IF( ln_dynadv_vec ) THEN 158 zua(:,:,:) = ( ua(:,:,:) - ub(:,:,:) ) * r1_2dt 159 zva(:,:,:) = ( va(:,:,:) - vb(:,:,:) ) * r1_2dt 160 ELSE 161 zua(:,:,:) = ( e3u_a(:,:,:)*ua(:,:,:) - e3u_b(:,:,:)*ub(:,:,:) ) / e3u_n(:,:,:) * r1_2dt 162 zva(:,:,:) = ( e3v_a(:,:,:)*va(:,:,:) - e3v_b(:,:,:)*vb(:,:,:) ) / e3v_n(:,:,:) * r1_2dt 163 ENDIF 164 CALL iom_put( "utrd_tot", zua ) ! total momentum trends, except the asselin filter 165 165 CALL iom_put( "vtrd_tot", zva ) 166 166 ENDIF 167 ! 168 zua(:,:,:) = un(:,:,:) ! save the now velocity before the asselin filter 169 zva(:,:,:) = vn(:,:,:) ! (caution: there will be a shift by 1 timestep in the 170 ! ! computation of the asselin filter trends) 167 zua(:,:,:) = un(:,:,:) ! save the now velocity before the asselin filter 168 zva(:,:,:) = vn(:,:,:) ! (caution: the Asselin filter trends computation will be shifted by 1 timestep) 171 169 ENDIF 172 170 173 171 ! Time filter and swap of dynamics arrays 174 ! --------------------------------------- ---175 IF( neuler == 0 .AND. kt == nit000 ) THEN !* Euler at first time-step: only swap172 ! --------------------------------------- 173 IF( l_1st_euler ) THEN !== Euler at 1st time-step ==! (swap only) 176 174 DO jk = 1, jpkm1 177 175 un(:,:,jk) = ua(:,:,jk) ! un <-- ua 178 176 vn(:,:,jk) = va(:,:,jk) 179 177 END DO 180 IF( .NOT.ln_linssh ) THEN ! e3._b <-- e3._n 181 !!gm BUG ???? I don't understand why it is not : e3._n <-- e3._a 178 IF( .NOT.ln_linssh ) THEN ! e3._n <-- e3._a 182 179 DO jk = 1, jpkm1 183 ! e3t_b(:,:,jk) = e3t_n(:,:,jk)184 ! e3u_b(:,:,jk) = e3u_n(:,:,jk)185 ! e3v_b(:,:,jk) = e3v_n(:,:,jk)186 !187 180 e3t_n(:,:,jk) = e3t_a(:,:,jk) 188 181 e3u_n(:,:,jk) = e3u_a(:,:,jk) 189 182 e3v_n(:,:,jk) = e3v_a(:,:,jk) 190 183 END DO 191 !!gm BUG end 192 ENDIF 193 ! 194 195 ELSE !* Leap-Frog : Asselin filter and swap 184 ENDIF 185 ! 186 ELSE !== Leap-Frog ==! (Asselin filter and swap) 187 ! 196 188 ! ! =============! 197 189 IF( ln_linssh ) THEN ! Fixed volume ! … … 213 205 ELSE ! Variable volume ! 214 206 ! ! ================! 215 ! Before scale factor at t-points 216 ! (used as a now filtered scale factor until the swap) 217 ! ---------------------------------------------------- 207 ! Before scale factor at t-points (used as a now filtered scale factor until the swap) 218 208 DO jk = 1, jpkm1 219 209 e3t_b(:,:,jk) = e3t_n(:,:,jk) + atfp * ( e3t_b(:,:,jk) - 2._wp * e3t_n(:,:,jk) + e3t_a(:,:,jk) ) 220 210 END DO 221 211 ! Add volume filter correction: compatibility with tracer advection scheme 222 ! => time filter + conservation correction (only at the first level)212 ! => time filter + conservation correction (only at the first level) 223 213 zcoef = atfp * rdt * r1_rau0 224 214 … … 232 222 IF( jk <= nk_rnf(ji,jj) ) THEN 233 223 e3t_b(ji,jj,jk) = e3t_b(ji,jj,jk) - zcoef * ( - rnf_b(ji,jj) + rnf(ji,jj) ) & 234 &* ( e3t_n(ji,jj,jk) / h_rnf(ji,jj) ) * tmask(ji,jj,jk)224 & * ( e3t_n(ji,jj,jk) / h_rnf(ji,jj) ) * tmask(ji,jj,jk) 235 225 ENDIF 236 END DO237 END DO238 END DO226 END DO 227 END DO 228 END DO 239 229 ELSE 240 230 e3t_b(:,:,1) = e3t_b(:,:,1) - zcoef * ( -rnf_b(:,:) + rnf(:,:))*tmask(:,:,1) 241 231 ENDIF 242 END 232 ENDIF 243 233 244 234 IF ( ln_isf ) THEN ! if ice shelf melting … … 253 243 END DO 254 244 END DO 255 END 245 ENDIF 256 246 ! 257 247 IF( ln_dynadv_vec ) THEN ! Asselin filter applied on velocity … … 322 312 ENDIF 323 313 ! 324 ENDIF ! neuler =/0314 ENDIF ! end Leap-Frog time stepping 325 315 ! 326 316 ! Set "now" and "before" barotropic velocities for next time step: 327 ! JC: Would be more clever to swap variables than to make a full vertical 328 ! integration 317 ! JC: Would be more clever to swap variables than to make a full vertical integration 329 318 ! 330 319 ! … … 360 349 ENDIF 361 350 IF( l_trddyn ) THEN ! 3D output: asselin filter trends on momentum 362 zua(:,:,:) = ( ub(:,:,:) - zua(:,:,:) ) * z1_2dt363 zva(:,:,:) = ( vb(:,:,:) - zva(:,:,:) ) * z1_2dt351 zua(:,:,:) = ( ub(:,:,:) - zua(:,:,:) ) * r1_2dt 352 zva(:,:,:) = ( vb(:,:,:) - zva(:,:,:) ) * r1_2dt 364 353 CALL trd_dyn( zua, zva, jpdyn_atf, kt ) 365 354 ENDIF
Note: See TracChangeset
for help on using the changeset viewer.