- Timestamp:
- 2015-10-31T08:40:45+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2015/dev_r5836_NOC3_vvl_by_default/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90
r5643 r5845 55 55 PUBLIC dyn_nxt ! routine called by step.F90 56 56 57 !! * Substitutions58 # include "domzgr_substitute.h90"59 57 !!---------------------------------------------------------------------- 60 58 !! NEMO/OPA 3.3 , NEMO Consortium (2010) … … 146 144 ELSE ! applied on thickness weighted velocity 147 145 DO jk = 1, jpkm1 148 ua(:,:,jk) = ( ub(:,:,jk) * fse3u_b(:,:,jk) & 149 & + z2dt * ua(:,:,jk) * fse3u_n(:,:,jk) ) & 150 & / fse3u_a(:,:,jk) * umask(:,:,jk) 151 va(:,:,jk) = ( vb(:,:,jk) * fse3v_b(:,:,jk) & 152 & + z2dt * va(:,:,jk) * fse3v_n(:,:,jk) ) & 153 & / fse3v_a(:,:,jk) * vmask(:,:,jk) 146 ua(:,:,jk) = ( ub(:,:,jk) * e3u_b(:,:,jk) & 147 & + z2dt * ua(:,:,jk) * e3u_n(:,:,jk) ) / e3u_a(:,:,jk) * umask(:,:,jk) 148 va(:,:,jk) = ( vb(:,:,jk) * e3v_b(:,:,jk) & 149 & + z2dt * va(:,:,jk) * e3v_n(:,:,jk) ) / e3v_a(:,:,jk) * vmask(:,:,jk) 154 150 END DO 155 151 ENDIF … … 160 156 ! Ensure below that barotropic velocities match time splitting estimate 161 157 ! Compute actual transport and replace it with ts estimate at "after" time step 162 zue(:,:) = fse3u_a(:,:,1) * ua(:,:,1) * umask(:,:,1)163 zve(:,:) = fse3v_a(:,:,1) * va(:,:,1) * vmask(:,:,1)158 zue(:,:) = e3u_a(:,:,1) * ua(:,:,1) * umask(:,:,1) 159 zve(:,:) = e3v_a(:,:,1) * va(:,:,1) * vmask(:,:,1) 164 160 DO jk = 2, jpkm1 165 zue(:,:) = zue(:,:) + fse3u_a(:,:,jk) * ua(:,:,jk) * umask(:,:,jk)166 zve(:,:) = zve(:,:) + fse3v_a(:,:,jk) * va(:,:,jk) * vmask(:,:,jk)161 zue(:,:) = zue(:,:) + e3u_a(:,:,jk) * ua(:,:,jk) * umask(:,:,jk) 162 zve(:,:) = zve(:,:) + e3v_a(:,:,jk) * va(:,:,jk) * vmask(:,:,jk) 167 163 END DO 168 164 DO jk = 1, jpkm1 169 ua(:,:,jk) = ( ua(:,:,jk) - zue(:,:) * hur_a(:,:) + ua_b(:,:) ) * umask(:,:,jk)170 va(:,:,jk) = ( va(:,:,jk) - zve(:,:) * hvr_a(:,:) + va_b(:,:) ) * vmask(:,:,jk)165 ua(:,:,jk) = ( ua(:,:,jk) - zue(:,:) * r1_hu_a(:,:) + ua_b(:,:) ) * umask(:,:,jk) 166 va(:,:,jk) = ( va(:,:,jk) - zve(:,:) * r1_hv_a(:,:) + va_b(:,:) ) * vmask(:,:,jk) 171 167 END DO 172 168 … … 231 227 IF (lk_vvl) THEN 232 228 DO jk = 1, jpkm1 233 fse3t_b(:,:,jk) = fse3t_n(:,:,jk)234 fse3u_b(:,:,jk) = fse3u_n(:,:,jk)235 fse3v_b(:,:,jk) = fse3v_n(:,:,jk)229 e3t_b(:,:,jk) = e3t_n(:,:,jk) 230 e3u_b(:,:,jk) = e3u_n(:,:,jk) 231 e3v_b(:,:,jk) = e3v_n(:,:,jk) 236 232 ENDDO 237 233 ENDIF … … 261 257 IF (lk_dynspg_ts.AND.ln_bt_fw) THEN 262 258 ! No asselin filtering on thicknesses if forward time splitting 263 fse3t_b(:,:,:) = fse3t_n(:,:,:)259 e3t_b(:,:,:) = e3t_n(:,:,:) 264 260 ELSE 265 fse3t_b(:,:,:) = fse3t_n(:,:,:) + atfp * ( fse3t_b(:,:,:) - 2._wp * fse3t_n(:,:,:) + fse3t_a(:,:,:) )261 e3t_b(:,:,:) = e3t_n(:,:,:) + atfp * ( e3t_b(:,:,:) - 2._wp * e3t_n(:,:,:) + e3t_a(:,:,:) ) 266 262 ! Add volume filter correction: compatibility with tracer advection scheme 267 263 ! => time filter + conservation correction (only at the first level) 268 264 IF ( nn_isf == 0) THEN ! if no ice shelf melting 269 fse3t_b(:,:,1) = fse3t_b(:,:,1) - atfp * rdt * r1_rau0 * ( emp_b(:,:) - emp(:,:) &265 e3t_b(:,:,1) = e3t_b(:,:,1) - atfp * rdt * r1_rau0 * ( emp_b(:,:) - emp(:,:) & 270 266 & -rnf_b(:,:) + rnf(:,:) ) * tmask(:,:,1) 271 267 ELSE ! if ice shelf melting … … 273 269 DO ji = 1,jpi 274 270 jk = mikt(ji,jj) 275 fse3t_b(ji,jj,jk) = fse3t_b(ji,jj,jk) - atfp * rdt * r1_rau0 &276 & 277 & 278 & 271 e3t_b(ji,jj,jk) = e3t_b(ji,jj,jk) - atfp * rdt * r1_rau0 & 272 & * ( (emp_b(ji,jj) - emp(ji,jj) ) & 273 & - (rnf_b(ji,jj) - rnf(ji,jj) ) & 274 & + (fwfisf_b(ji,jj) - fwfisf(ji,jj)) ) * tmask(ji,jj,jk) 279 275 END DO 280 276 END DO … … 285 281 ! Before scale factor at (u/v)-points 286 282 ! ----------------------------------- 287 CALL dom_vvl_interpol( fse3t_b(:,:,:), fse3u_b(:,:,:), 'U' )288 CALL dom_vvl_interpol( fse3t_b(:,:,:), fse3v_b(:,:,:), 'V' )283 CALL dom_vvl_interpol( e3t_b(:,:,:), e3u_b(:,:,:), 'U' ) 284 CALL dom_vvl_interpol( e3t_b(:,:,:), e3v_b(:,:,:), 'V' ) 289 285 ! Leap-Frog - Asselin filter and swap: applied on velocity 290 286 ! ----------------------------------- … … 306 302 ! Temporary filtered scale factor at (u/v)-points (will become before scale factor) 307 303 !------------------------------------------------ 308 CALL dom_vvl_interpol( fse3t_b(:,:,:), ze3u_f, 'U' )309 CALL dom_vvl_interpol( fse3t_b(:,:,:), ze3v_f, 'V' )304 CALL dom_vvl_interpol( e3t_b(:,:,:), ze3u_f, 'U' ) 305 CALL dom_vvl_interpol( e3t_b(:,:,:), ze3v_f, 'V' ) 310 306 ! Leap-Frog - Asselin filter and swap: applied on thickness weighted velocity 311 307 ! ----------------------------------- =========================== … … 313 309 DO jj = 1, jpj 314 310 DO ji = 1, jpi 315 zue3a = ua(ji,jj,jk) * fse3u_a(ji,jj,jk)316 zve3a = va(ji,jj,jk) * fse3v_a(ji,jj,jk)317 zue3n = un(ji,jj,jk) * fse3u_n(ji,jj,jk)318 zve3n = vn(ji,jj,jk) * fse3v_n(ji,jj,jk)319 zue3b = ub(ji,jj,jk) * fse3u_b(ji,jj,jk)320 zve3b = vb(ji,jj,jk) * fse3v_b(ji,jj,jk)311 zue3a = ua(ji,jj,jk) * e3u_a(ji,jj,jk) 312 zve3a = va(ji,jj,jk) * e3v_a(ji,jj,jk) 313 zue3n = un(ji,jj,jk) * e3u_n(ji,jj,jk) 314 zve3n = vn(ji,jj,jk) * e3v_n(ji,jj,jk) 315 zue3b = ub(ji,jj,jk) * e3u_b(ji,jj,jk) 316 zve3b = vb(ji,jj,jk) * e3v_b(ji,jj,jk) 321 317 ! 322 318 zuf = ( zue3n + atfp * ( zue3b - 2._wp * zue3n + zue3a ) ) / ze3u_f(ji,jj,jk) … … 330 326 END DO 331 327 END DO 332 fse3u_b(:,:,1:jpkm1) = ze3u_f(:,:,1:jpkm1) ! e3u_b <-- filtered scale factor333 fse3v_b(:,:,1:jpkm1) = ze3v_f(:,:,1:jpkm1)328 e3u_b(:,:,1:jpkm1) = ze3u_f(:,:,1:jpkm1) ! e3u_b <-- filtered scale factor 329 e3v_b(:,:,1:jpkm1) = ze3v_f(:,:,1:jpkm1) 334 330 ENDIF 335 331 ! … … 339 335 ! Revert "before" velocities to time split estimate 340 336 ! Doing it here also means that asselin filter contribution is removed 341 zue(:,:) = fse3u_b(:,:,1) * ub(:,:,1) * umask(:,:,1)342 zve(:,:) = fse3v_b(:,:,1) * vb(:,:,1) * vmask(:,:,1)337 zue(:,:) = e3u_b(:,:,1) * ub(:,:,1) * umask(:,:,1) 338 zve(:,:) = e3v_b(:,:,1) * vb(:,:,1) * vmask(:,:,1) 343 339 DO jk = 2, jpkm1 344 zue(:,:) = zue(:,:) + fse3u_b(:,:,jk) * ub(:,:,jk) * umask(:,:,jk)345 zve(:,:) = zve(:,:) + fse3v_b(:,:,jk) * vb(:,:,jk) * vmask(:,:,jk)340 zue(:,:) = zue(:,:) + e3u_b(:,:,jk) * ub(:,:,jk) * umask(:,:,jk) 341 zve(:,:) = zve(:,:) + e3v_b(:,:,jk) * vb(:,:,jk) * vmask(:,:,jk) 346 342 END DO 347 343 DO jk = 1, jpkm1 348 ub(:,:,jk) = ub(:,:,jk) - (zue(:,:) * hur(:,:) - un_b(:,:)) * umask(:,:,jk)349 vb(:,:,jk) = vb(:,:,jk) - (zve(:,:) * hvr(:,:) - vn_b(:,:)) * vmask(:,:,jk)344 ub(:,:,jk) = ub(:,:,jk) - (zue(:,:) * r1_hu_n(:,:) - un_b(:,:)) * umask(:,:,jk) 345 vb(:,:,jk) = vb(:,:,jk) - (zve(:,:) * r1_hv_n(:,:) - vn_b(:,:)) * vmask(:,:,jk) 350 346 END DO 351 347 ENDIF … … 359 355 ! 360 356 IF (lk_vvl) THEN 361 hu_b(:,:) = 0.362 hv_b(:,:) = 0.363 DO jk = 1, jpkm1364 hu_b(:,:) = hu_b(:,:) + fse3u_b(:,:,jk) * umask(:,:,jk)365 hv_b(:,:) = hv_b(:,:) + fse3v_b(:,:,jk) * vmask(:,:,jk)357 hu_b(:,:) = e3u_b(:,:,1) * umask(:,:,1) 358 hv_b(:,:) = e3v_b(:,:,1) * vmask(:,:,1) 359 DO jk = 2, jpkm1 360 hu_b(:,:) = hu_b(:,:) + e3u_b(:,:,jk) * umask(:,:,jk) 361 hv_b(:,:) = hv_b(:,:) + e3v_b(:,:,jk) * vmask(:,:,jk) 366 362 END DO 367 hur_b(:,:) = umask_i(:,:) / ( hu_b(:,:) + 1._wp - umask_i(:,:) ) 368 hvr_b(:,:) = vmask_i(:,:) / ( hv_b(:,:) + 1._wp - vmask_i(:,:) )369 ENDIF370 !371 un_b(:,:) = 0._wp ; vn_b(:,:) = 0._wp372 u b_b(:,:) = 0._wp ; vb_b(:,:) = 0._wp373 !363 !!gm don't understand the use of umask_i .... 364 r1_hu_b(:,:) = umask_i(:,:) / ( hu_b(:,:) + 1._wp - umask_i(:,:) ) 365 r1_hv_b(:,:) = vmask_i(:,:) / ( hv_b(:,:) + 1._wp - vmask_i(:,:) ) 366 ENDIF 367 ! 368 un_b(:,:) = 0._wp ; vn_b(:,:) = 0._wp 369 ub_b(:,:) = 0._wp ; vb_b(:,:) = 0._wp 374 370 DO jk = 1, jpkm1 375 371 DO jj = 1, jpj 376 372 DO ji = 1, jpi 377 un_b(ji,jj) = un_b(ji,jj) + fse3u_a(ji,jj,jk) * un(ji,jj,jk) * umask(ji,jj,jk)378 vn_b(ji,jj) = vn_b(ji,jj) + fse3v_a(ji,jj,jk) * vn(ji,jj,jk) * vmask(ji,jj,jk)373 un_b(ji,jj) = un_b(ji,jj) + e3u_a(ji,jj,jk) * un(ji,jj,jk) * umask(ji,jj,jk) 374 vn_b(ji,jj) = vn_b(ji,jj) + e3v_a(ji,jj,jk) * vn(ji,jj,jk) * vmask(ji,jj,jk) 379 375 ! 380 ub_b(ji,jj) = ub_b(ji,jj) + fse3u_b(ji,jj,jk) * ub(ji,jj,jk) * umask(ji,jj,jk)381 vb_b(ji,jj) = vb_b(ji,jj) + fse3v_b(ji,jj,jk) * vb(ji,jj,jk) * vmask(ji,jj,jk)376 ub_b(ji,jj) = ub_b(ji,jj) + e3u_b(ji,jj,jk) * ub(ji,jj,jk) * umask(ji,jj,jk) 377 vb_b(ji,jj) = vb_b(ji,jj) + e3v_b(ji,jj,jk) * vb(ji,jj,jk) * vmask(ji,jj,jk) 382 378 END DO 383 379 END DO 384 380 END DO 385 ! 386 ! 387 un_b(:,:) = un_b(:,:) * hur_a(:,:) 388 vn_b(:,:) = vn_b(:,:) * hvr_a(:,:) 389 ub_b(:,:) = ub_b(:,:) * hur_b(:,:) 390 vb_b(:,:) = vb_b(:,:) * hvr_b(:,:) 391 ! 392 ! 393 381 un_b(:,:) = un_b(:,:) * r1_hu_a(:,:) 382 vn_b(:,:) = vn_b(:,:) * r1_hv_a(:,:) 383 ub_b(:,:) = ub_b(:,:) * r1_hu_b(:,:) 384 vb_b(:,:) = vb_b(:,:) * r1_hv_b(:,:) 385 ! 394 386 IF( l_trddyn ) THEN ! 3D output: asselin filter trends on momentum 395 387 zua(:,:,:) = ( ub(:,:,:) - zua(:,:,:) ) * z1_2dt
Note: See TracChangeset
for help on using the changeset viewer.
