- Timestamp:
- 2019-11-22T15:29:17+01:00 (4 years ago)
- Location:
- NEMO/branches/2019/dev_r11943_MERGE_2019/src
- Files:
-
- 2 edited
-
. (modified) (1 prop)
-
NST/agrif_oce_update.F90 (modified) (54 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r11943_MERGE_2019/src
- Property svn:mergeinfo deleted
-
NEMO/branches/2019/dev_r11943_MERGE_2019/src/NST/agrif_oce_update.F90
r10068 r11949 230 230 ! ----------------------- 231 231 ! 232 e3u _a(:,:,:) = e3u_n(:,:,:)233 e3v _a(:,:,:) = e3v_n(:,:,:)234 ! u a(:,:,:) = e3u_b(:,:,:)235 ! v a(:,:,:) = e3v_b(:,:,:)236 hu _a(:,:) = hu_n(:,:)237 hv _a(:,:) = hv_n(:,:)232 e3u(:,:,:,Krhs_a) = e3u(:,:,:,Kmm_a) 233 e3v(:,:,:,Krhs_a) = e3v(:,:,:,Kmm_a) 234 ! uu(:,:,:,Krhs_a) = e3u(:,:,:,Kbb_a) 235 ! vv(:,:,:,Krhs_a) = e3v(:,:,:,Kbb_a) 236 hu(:,:,Krhs_a) = hu(:,:,Kmm_a) 237 hv(:,:,Krhs_a) = hv(:,:,Kmm_a) 238 238 239 239 ! 1) NOW fields … … 242 242 ! Vertical scale factor interpolations 243 243 ! ------------------------------------ 244 CALL dom_vvl_interpol( e3t _n(:,:,:), e3u_n(:,:,:) , 'U' )245 CALL dom_vvl_interpol( e3t _n(:,:,:), e3v_n(:,:,:) , 'V' )246 CALL dom_vvl_interpol( e3u _n(:,:,:), e3f_n(:,:,:) , 'F' )247 248 CALL dom_vvl_interpol( e3u _n(:,:,:), e3uw_n(:,:,:), 'UW' )249 CALL dom_vvl_interpol( e3v _n(:,:,:), e3vw_n(:,:,:), 'VW' )244 CALL dom_vvl_interpol( e3t(:,:,:,Kmm_a), e3u(:,:,:,Kmm_a) , 'U' ) 245 CALL dom_vvl_interpol( e3t(:,:,:,Kmm_a), e3v(:,:,:,Kmm_a) , 'V' ) 246 CALL dom_vvl_interpol( e3u(:,:,:,Kmm_a), e3f(:,:,:) , 'F' ) 247 248 CALL dom_vvl_interpol( e3u(:,:,:,Kmm_a), e3uw(:,:,:,Kmm_a), 'UW' ) 249 CALL dom_vvl_interpol( e3v(:,:,:,Kmm_a), e3vw(:,:,:,Kmm_a), 'VW' ) 250 250 251 251 ! Update total depths: 252 252 ! -------------------- 253 hu _n(:,:) = 0._wp ! Ocean depth at U-points254 hv _n(:,:) = 0._wp ! Ocean depth at V-points253 hu(:,:,Kmm_a) = 0._wp ! Ocean depth at U-points 254 hv(:,:,Kmm_a) = 0._wp ! Ocean depth at V-points 255 255 DO jk = 1, jpkm1 256 hu _n(:,:) = hu_n(:,:) + e3u_n(:,:,jk) * umask(:,:,jk)257 hv _n(:,:) = hv_n(:,:) + e3v_n(:,:,jk) * vmask(:,:,jk)256 hu(:,:,Kmm_a) = hu(:,:,Kmm_a) + e3u(:,:,jk,Kmm_a) * umask(:,:,jk) 257 hv(:,:,Kmm_a) = hv(:,:,Kmm_a) + e3v(:,:,jk,Kmm_a) * vmask(:,:,jk) 258 258 END DO 259 259 ! ! Inverse of the local depth 260 r1_hu _n(:,:) = ssumask(:,:) / ( hu_n(:,:) + 1._wp - ssumask(:,:) )261 r1_hv _n(:,:) = ssvmask(:,:) / ( hv_n(:,:) + 1._wp - ssvmask(:,:) )260 r1_hu(:,:,Kmm_a) = ssumask(:,:) / ( hu(:,:,Kmm_a) + 1._wp - ssumask(:,:) ) 261 r1_hv(:,:,Kmm_a) = ssvmask(:,:) / ( hv(:,:,Kmm_a) + 1._wp - ssvmask(:,:) ) 262 262 263 263 … … 268 268 ! Vertical scale factor interpolations 269 269 ! ------------------------------------ 270 CALL dom_vvl_interpol( e3t _b(:,:,:), e3u_b(:,:,:), 'U' )271 CALL dom_vvl_interpol( e3t _b(:,:,:), e3v_b(:,:,:), 'V' )272 273 CALL dom_vvl_interpol( e3u _b(:,:,:), e3uw_b(:,:,:), 'UW' )274 CALL dom_vvl_interpol( e3v _b(:,:,:), e3vw_b(:,:,:), 'VW' )270 CALL dom_vvl_interpol( e3t(:,:,:,Kbb_a), e3u(:,:,:,Kbb_a), 'U' ) 271 CALL dom_vvl_interpol( e3t(:,:,:,Kbb_a), e3v(:,:,:,Kbb_a), 'V' ) 272 273 CALL dom_vvl_interpol( e3u(:,:,:,Kbb_a), e3uw(:,:,:,Kbb_a), 'UW' ) 274 CALL dom_vvl_interpol( e3v(:,:,:,Kbb_a), e3vw(:,:,:,Kbb_a), 'VW' ) 275 275 276 276 ! Update total depths: 277 277 ! -------------------- 278 hu _b(:,:) = 0._wp ! Ocean depth at U-points279 hv _b(:,:) = 0._wp ! Ocean depth at V-points278 hu(:,:,Kbb_a) = 0._wp ! Ocean depth at U-points 279 hv(:,:,Kbb_a) = 0._wp ! Ocean depth at V-points 280 280 DO jk = 1, jpkm1 281 hu _b(:,:) = hu_b(:,:) + e3u_b(:,:,jk) * umask(:,:,jk)282 hv _b(:,:) = hv_b(:,:) + e3v_b(:,:,jk) * vmask(:,:,jk)281 hu(:,:,Kbb_a) = hu(:,:,Kbb_a) + e3u(:,:,jk,Kbb_a) * umask(:,:,jk) 282 hv(:,:,Kbb_a) = hv(:,:,Kbb_a) + e3v(:,:,jk,Kbb_a) * vmask(:,:,jk) 283 283 END DO 284 284 ! ! Inverse of the local depth 285 r1_hu _b(:,:) = ssumask(:,:) / ( hu_b(:,:) + 1._wp - ssumask(:,:) )286 r1_hv _b(:,:) = ssvmask(:,:) / ( hv_b(:,:) + 1._wp - ssvmask(:,:) )285 r1_hu(:,:,Kbb_a) = ssumask(:,:) / ( hu(:,:,Kbb_a) + 1._wp - ssumask(:,:) ) 286 r1_hv(:,:,Kbb_a) = ssvmask(:,:) / ( hv(:,:,Kbb_a) + 1._wp - ssvmask(:,:) ) 287 287 ENDIF 288 288 ! … … 315 315 DO jj=j1,j2 316 316 DO ji=i1,i2 317 tabres(ji,jj,jk,jn) = ( tsn(ji,jj,jk,jn) * e3t_n(ji,jj,jk) )&318 * tmask(ji,jj,jk) + (tmask(ji,jj,jk)-1)*999._wp317 tabres(ji,jj,jk,jn) = ( ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Kmm_a) ) & 318 & * tmask(ji,jj,jk) + (tmask(ji,jj,jk) - 1._wp)*999._wp 319 319 END DO 320 320 END DO … … 324 324 DO jj=j1,j2 325 325 DO ji=i1,i2 326 tabres(ji,jj,jk,n2) = tmask(ji,jj,jk) * e3t _n(ji,jj,jk) &327 + (tmask(ji,jj,jk)-1)*999._wp326 tabres(ji,jj,jk,n2) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm_a) & 327 & + (tmask(ji,jj,jk) - 1._wp)*999._wp 328 328 END DO 329 329 END DO 330 330 END DO 331 331 ELSE 332 tabres_child(:,:,:,:) = 0. 332 tabres_child(:,:,:,:) = 0._wp 333 333 AGRIF_SpecialValue = 0._wp 334 334 DO jj=j1,j2 … … 345 345 IF (tmask(ji,jj,jk) < -900) EXIT ! TODO: Will not work with ISF 346 346 N_out = N_out + 1 347 h_out(N_out) = e3t _n(ji,jj,jk)347 h_out(N_out) = e3t(ji,jj,jk,Kmm_a) 348 348 ENDDO 349 349 IF (N_in > 0) THEN !Remove this? … … 356 356 ENDIF 357 357 DO jn=n1,n2-1 358 CALL reconstructandremap(tabin(1:N_in,jn),h_in(1:N_in),tabres_child(ji,jj,1:N_out,jn),h_out(1:N_out),N_in,N_out) 358 CALL reconstructandremap( tabin(1:N_in,jn), h_in(1:N_in), tabres_child(ji,jj,1:N_out,jn), & 359 & h_out(1:N_out) , N_in , N_out ) 359 360 ENDDO 360 361 ENDIF … … 369 370 DO ji=i1,i2 370 371 IF( tabres_child(ji,jj,jk,jn) .NE. 0. ) THEN 371 tsb(ji,jj,jk,jn) = tsb(ji,jj,jk,jn) & 372 & + atfp * ( tabres_child(ji,jj,jk,jn) & 373 & - tsn(ji,jj,jk,jn) ) * tmask(ji,jj,jk) 372 ts(ji,jj,jk,jn,Kbb_a) = ts(ji,jj,jk,jn,Kbb_a) & 373 & + atfp * ( tabres_child(ji,jj,jk,jn) & 374 & - ts(ji,jj,jk,jn,Kmm_a) & 375 & ) * tmask(ji,jj,jk) 374 376 ENDIF 375 377 ENDDO … … 383 385 DO ji=i1,i2 384 386 IF( tabres_child(ji,jj,jk,jn) .NE. 0. ) THEN 385 ts n(ji,jj,jk,jn) = tabres_child(ji,jj,jk,jn) * tmask(ji,jj,jk)387 ts(ji,jj,jk,jn,Kmm_a) = tabres_child(ji,jj,jk,jn) * tmask(ji,jj,jk) 386 388 END IF 387 389 END DO … … 413 415 DO ji=i1,i2 414 416 !> jc tmp 415 tabres(ji,jj,jk,jn) = ts n(ji,jj,jk,jn) * e3t_n(ji,jj,jk) / e3t_0(ji,jj,jk)416 ! tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn) * e3t_n(ji,jj,jk)417 tabres(ji,jj,jk,jn) = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Kmm_a) / e3t_0(ji,jj,jk) 418 ! tabres(ji,jj,jk,jn) = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Kmm_a) 417 419 !< jc tmp 418 420 END DO … … 434 436 DO ji = i1, i2 435 437 IF( tabres(ji,jj,jk,jn) /= 0._wp ) THEN 436 ztb = ts b(ji,jj,jk,jn) * e3t_b(ji,jj,jk) ! fse3t_b prior update should be used438 ztb = ts(ji,jj,jk,jn,Kbb_a) * e3t(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used 437 439 ztnu = tabres(ji,jj,jk,jn) 438 ztno = ts n(ji,jj,jk,jn) * e3t_a(ji,jj,jk)439 ts b(ji,jj,jk,jn) = ( ztb + atfp * ( ztnu - ztno) ) &440 & * tmask(ji,jj,jk) / e3t_b(ji,jj,jk)440 ztno = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Krhs_a) 441 ts(ji,jj,jk,jn,Kbb_a) = ( ztb + atfp * ( ztnu - ztno) ) & 442 & * tmask(ji,jj,jk) / e3t(ji,jj,jk,Kbb_a) 441 443 ENDIF 442 444 END DO … … 450 452 DO ji=i1,i2 451 453 IF( tabres(ji,jj,jk,jn) /= 0._wp ) THEN 452 ts n(ji,jj,jk,jn) = tabres(ji,jj,jk,jn) / e3t_n(ji,jj,jk)454 ts(ji,jj,jk,jn,Kmm_a) = tabres(ji,jj,jk,jn) / e3t(ji,jj,jk,Kmm_a) 453 455 END IF 454 456 END DO … … 458 460 ! 459 461 IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 460 ts b(i1:i2,j1:j2,k1:k2,1:jpts) = tsn(i1:i2,j1:j2,k1:k2,1:jpts)462 ts(i1:i2,j1:j2,k1:k2,1:jpts,Kbb_a) = ts(i1:i2,j1:j2,k1:k2,1:jpts,Kmm_a) 461 463 ENDIF 462 464 ! … … 495 497 DO jj=j1,j2 496 498 DO ji=i1,i2 497 tabres(ji,jj,jk,1) = zrhoy * e2u(ji,jj) * e3u _n(ji,jj,jk) * umask(ji,jj,jk) * un(ji,jj,jk) &498 + (umask(ji,jj,jk)-1)*999._wp499 tabres(ji,jj,jk,2) = zrhoy * umask(ji,jj,jk) * e2u(ji,jj) * e3u _n(ji,jj,jk) &500 + (umask(ji,jj,jk)-1)*999._wp499 tabres(ji,jj,jk,1) = zrhoy * e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) * umask(ji,jj,jk) * uu(ji,jj,jk,Kmm_a) & 500 & + (umask(ji,jj,jk)-1._wp)*999._wp 501 tabres(ji,jj,jk,2) = zrhoy * umask(ji,jj,jk) * e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) & 502 & + (umask(ji,jj,jk)-1._wp)*999._wp 501 503 END DO 502 504 END DO … … 513 515 IF( tabres(ji,jj,jk,2) < -900) EXIT 514 516 N_in = N_in + 1 515 tabin(jk) = tabres(ji,jj,jk,1)/tabres(ji,jj,jk,2)517 tabin(jk) = tabres(ji,jj,jk,1)/tabres(ji,jj,jk,2) 516 518 h_in(N_in) = tabres(ji,jj,jk,2)/e2u(ji,jj) 517 519 ENDDO … … 520 522 IF (umask(ji,jj,jk) == 0) EXIT 521 523 N_out = N_out + 1 522 h_out(N_out) = e3u _n(ji,jj,jk)524 h_out(N_out) = e3u(ji,jj,jk,Kmm_a) 523 525 ENDDO 524 526 IF (N_in * N_out > 0) THEN 525 h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in))527 h_diff = SUM( h_out(1:N_out) ) - SUM( h_in(1:N_in) ) 526 528 IF (h_diff < -1.e-4) THEN 527 529 !Even if bathy at T points match it's possible for the U points to be deeper in the child grid. … … 540 542 ENDDO 541 543 ENDIF 542 CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),tabres_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out) 543 tabres_child(ji,jj,N_out) = tabres_child(ji,jj,N_out) + excess/(e2u(ji,jj)*h_out(N_out)) 544 CALL reconstructandremap( tabin(1:N_in) , h_in(1:N_in), tabres_child(ji,jj,1:N_out), & 545 & h_out(1:N_out), N_in , N_out ) 546 tabres_child(ji,jj,N_out) = tabres_child(ji,jj,N_out) + excess/( e2u(ji,jj)*h_out(N_out) ) 544 547 ENDIF 545 548 ENDDO … … 550 553 DO ji=i1,i2 551 554 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 552 u b(ji,jj,jk) = ub(ji,jj,jk)&553 & + atfp * ( tabres_child(ji,jj,jk) - un(ji,jj,jk) ) * umask(ji,jj,jk)555 uu(ji,jj,jk,Kbb_a) = uu(ji,jj,jk,Kbb_a) & 556 & + atfp * ( tabres_child(ji,jj,jk) - uu(ji,jj,jk,Kmm_a) ) * umask(ji,jj,jk) 554 557 ENDIF 555 558 ! 556 u n(ji,jj,jk) = tabres_child(ji,jj,jk) * umask(ji,jj,jk)559 uu(ji,jj,jk,Kmm_a) = tabres_child(ji,jj,jk) * umask(ji,jj,jk) 557 560 END DO 558 561 END DO … … 579 582 zrhoy = Agrif_Rhoy() 580 583 DO jk = k1, k2 581 tabres(i1:i2,j1:j2,jk,1) = zrhoy * e2u(i1:i2,j1:j2) * e3u_n(i1:i2,j1:j2,jk) * un(i1:i2,j1:j2,jk) 584 tabres(i1:i2,j1:j2,jk,1) = zrhoy * e2u(i1:i2,j1:j2) * e3u(i1:i2,j1:j2,jk,Kmm_a) & 585 & * uu(i1:i2,j1:j2,jk,Kmm_a) 582 586 END DO 583 587 ELSE … … 588 592 ! 589 593 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 590 zub = u b(ji,jj,jk) * e3u_b(ji,jj,jk) ! fse3t_b prior update should be used591 zuno = u n(ji,jj,jk) * e3u_a(ji,jj,jk)594 zub = uu(ji,jj,jk,Kbb_a) * e3u(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used 595 zuno = uu(ji,jj,jk,Kmm_a) * e3u(ji,jj,jk,Krhs_a) 592 596 zunu = tabres(ji,jj,jk,1) 593 u b(ji,jj,jk) = ( zub + atfp * ( zunu - zuno) )&594 & * umask(ji,jj,jk) / e3u_b(ji,jj,jk)597 uu(ji,jj,jk,Kbb_a) = ( zub + atfp * ( zunu - zuno ) ) & 598 & * umask(ji,jj,jk) / e3u(ji,jj,jk,Kbb_a) 595 599 ENDIF 596 600 ! 597 u n(ji,jj,jk) = tabres(ji,jj,jk,1) * umask(ji,jj,jk) / e3u_n(ji,jj,jk)601 uu(ji,jj,jk,Kmm_a) = tabres(ji,jj,jk,1) * umask(ji,jj,jk) / e3u(ji,jj,jk,Kmm_a) 598 602 END DO 599 603 END DO … … 601 605 ! 602 606 IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 603 u b(i1:i2,j1:j2,k1:k2) = un(i1:i2,j1:j2,k1:k2)607 uu(i1:i2,j1:j2,k1:k2,Kbb_a) = uu(i1:i2,j1:j2,k1:k2,Kmm_a) 604 608 ENDIF 605 609 ! … … 632 636 IF (western_side) THEN 633 637 DO jj=j1,j2 634 zcor = u n_b(i1-1,jj) * hu_a(i1-1,jj) * r1_hu_n(i1-1,jj) - un_b(i1-1,jj)635 u n_b(i1-1,jj) = un_b(i1-1,jj) + zcor638 zcor = uu_b(i1-1,jj,Kmm_a) * hu(i1-1,jj,Krhs_a) * r1_hu(i1-1,jj,Kmm_a) - uu_b(i1-1,jj,Kmm_a) 639 uu_b(i1-1,jj,Kmm_a) = uu_b(i1-1,jj,Kmm_a) + zcor 636 640 DO jk=1,jpkm1 637 u n(i1-1,jj,jk) = un(i1-1,jj,jk) + zcor * umask(i1-1,jj,jk)641 uu(i1-1,jj,jk,Kmm_a) = uu(i1-1,jj,jk,Kmm_a) + zcor * umask(i1-1,jj,jk) 638 642 END DO 639 643 END DO … … 642 646 IF (eastern_side) THEN 643 647 DO jj=j1,j2 644 zcor = u n_b(i2+1,jj) * hu_a(i2+1,jj) * r1_hu_n(i2+1,jj) - un_b(i2+1,jj)645 u n_b(i2+1,jj) = un_b(i2+1,jj) + zcor648 zcor = uu_b(i2+1,jj,Kmm_a) * hu(i2+1,jj,Krhs_a) * r1_hu(i2+1,jj,Kmm_a) - uu_b(i2+1,jj,Kmm_a) 649 uu_b(i2+1,jj,Kmm_a) = uu_b(i2+1,jj,Kmm_a) + zcor 646 650 DO jk=1,jpkm1 647 u n(i2+1,jj,jk) = un(i2+1,jj,jk) + zcor * umask(i2+1,jj,jk)651 uu(i2+1,jj,jk,Kmm_a) = uu(i2+1,jj,jk,Kmm_a) + zcor * umask(i2+1,jj,jk) 648 652 END DO 649 653 END DO … … 682 686 DO jj=j1,j2 683 687 DO ji=i1,i2 684 tabres(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v_n(ji,jj,jk) * vmask(ji,jj,jk) * vn(ji,jj,jk) & 685 + (vmask(ji,jj,jk)-1)*999._wp 686 tabres(ji,jj,jk,2) = vmask(ji,jj,jk) * zrhox * e1v(ji,jj) * e3v_n(ji,jj,jk) & 688 tabres(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) & 689 & * vmask(ji,jj,jk) * vv(ji,jj,jk,Kmm_a) & 690 & + (vmask(ji,jj,jk)-1)*999._wp 691 tabres(ji,jj,jk,2) = vmask(ji,jj,jk) * zrhox * e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) & 687 692 + (vmask(ji,jj,jk)-1)*999._wp 688 693 END DO … … 705 710 IF (vmask(ji,jj,jk) == 0) EXIT 706 711 N_out = N_out + 1 707 h_out(N_out) = e3v _n(ji,jj,jk)712 h_out(N_out) = e3v(ji,jj,jk,Kmm_a) 708 713 ENDDO 709 714 IF (N_in * N_out > 0) THEN 710 h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in))715 h_diff = SUM( h_out(1:N_out) ) - SUM( h_in(1:N_in) ) 711 716 IF (h_diff < -1.e-4) then 712 717 !Even if bathy at T points match it's possible for the U points to be deeper in the child grid. … … 714 719 excess = 0._wp 715 720 DO jk=N_in,1,-1 716 thick = MIN( -1*h_diff, h_in(jk))721 thick = MIN( -1._wp*h_diff, h_in(jk) ) 717 722 excess = excess + tabin(jk)*thick*e2u(ji,jj) 718 tabin(jk) = tabin(jk)*(1. - thick/h_in(jk))723 tabin(jk) = tabin(jk)*(1._wp - thick/h_in(jk)) 719 724 h_diff = h_diff + thick 720 725 IF ( h_diff == 0) THEN … … 725 730 ENDDO 726 731 ENDIF 727 CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),tabres_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out) 732 CALL reconstructandremap( tabin(1:N_in) , h_in(1:N_in), tabres_child(ji,jj,1:N_out), & 733 & h_out(1:N_out), N_in , N_out ) 728 734 tabres_child(ji,jj,N_out) = tabres_child(ji,jj,N_out) + excess/(e1v(ji,jj)*h_out(N_out)) 729 735 ENDIF … … 736 742 ! 737 743 IF( .NOT.(lk_agrif_fstep.AND.(neuler==0)) ) THEN ! Add asselin part 738 vb(ji,jj,jk) = vb(ji,jj,jk) & 739 & + atfp * ( tabres_child(ji,jj,jk) - vn(ji,jj,jk) ) * vmask(ji,jj,jk) 744 vv(ji,jj,jk,Kbb_a) = vv(ji,jj,jk,Kbb_a) & 745 & + atfp * ( tabres_child(ji,jj,jk) - vv(ji,jj,jk,Kmm_a) ) & 746 & * vmask(ji,jj,jk) 740 747 ENDIF 741 748 ! 742 v n(ji,jj,jk) = tabres_child(ji,jj,jk) * vmask(ji,jj,jk)749 vv(ji,jj,jk,Kmm_a) = tabres_child(ji,jj,jk) * vmask(ji,jj,jk) 743 750 END DO 744 751 END DO … … 767 774 DO jj=j1,j2 768 775 DO ji=i1,i2 769 tabres(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v _n(ji,jj,jk) * vn(ji,jj,jk)776 tabres(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) * vv(ji,jj,jk,Kmm_a) 770 777 END DO 771 778 END DO … … 778 785 ! 779 786 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 780 zvb = v b(ji,jj,jk) * e3v_b(ji,jj,jk) ! fse3t_b prior update should be used781 zvno = v n(ji,jj,jk) * e3v_a(ji,jj,jk)787 zvb = vv(ji,jj,jk,Kbb_a) * e3v(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used 788 zvno = vv(ji,jj,jk,Kmm_a) * e3v(ji,jj,jk,Krhs_a) 782 789 zvnu = tabres(ji,jj,jk,1) 783 v b(ji,jj,jk) = ( zvb + atfp * ( zvnu - zvno) )&784 & * vmask(ji,jj,jk) / e3v_b(ji,jj,jk)790 vv(ji,jj,jk,Kbb_a) = ( zvb + atfp * ( zvnu - zvno) ) & 791 & * vmask(ji,jj,jk) / e3v(ji,jj,jk,Kbb_a) 785 792 ENDIF 786 793 ! 787 v n(ji,jj,jk) = tabres(ji,jj,jk,1) * vmask(ji,jj,jk) / e3v_n(ji,jj,jk)794 vv(ji,jj,jk,Kmm_a) = tabres(ji,jj,jk,1) * vmask(ji,jj,jk) / e3v(ji,jj,jk,Kmm_a) 788 795 END DO 789 796 END DO … … 791 798 ! 792 799 IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 793 v b(i1:i2,j1:j2,k1:k2) = vn(i1:i2,j1:j2,k1:k2)800 vv(i1:i2,j1:j2,k1:k2,Kbb_a) = vv(i1:i2,j1:j2,k1:k2,Kmm_a) 794 801 ENDIF 795 802 ! … … 822 829 IF (southern_side) THEN 823 830 DO ji=i1,i2 824 zcor = v n_b(ji,j1-1) * hv_a(ji,j1-1) * r1_hv_n(ji,j1-1) - vn_b(ji,j1-1)825 v n_b(ji,j1-1) = vn_b(ji,j1-1) + zcor831 zcor = vv_b(ji,j1-1,Kmm_a) * hv(ji,j1-1,Krhs_a) * r1_hv(ji,j1-1,Kmm_a) - vv_b(ji,j1-1,Kmm_a) 832 vv_b(ji,j1-1,Kmm_a) = vv_b(ji,j1-1,Kmm_a) + zcor 826 833 DO jk=1,jpkm1 827 v n(ji,j1-1,jk) = vn(ji,j1-1,jk) + zcor * vmask(ji,j1-1,jk)834 vv(ji,j1-1,jk,Kmm_a) = vv(ji,j1-1,jk,Kmm_a) + zcor * vmask(ji,j1-1,jk) 828 835 END DO 829 836 END DO … … 832 839 IF (northern_side) THEN 833 840 DO ji=i1,i2 834 zcor = v n_b(ji,j2+1) * hv_a(ji,j2+1) * r1_hv_n(ji,j2+1) - vn_b(ji,j2+1)835 v n_b(ji,j2+1) = vn_b(ji,j2+1) + zcor841 zcor = vv_b(ji,j2+1,Kmm_a) * hv(ji,j2+1,Krhs_a) * r1_hv(ji,j2+1,Kmm_a) - vv_b(ji,j2+1,Kmm_a) 842 vv_b(ji,j2+1,Kmm_a) = vv_b(ji,j2+1,Kmm_a) + zcor 836 843 DO jk=1,jpkm1 837 v n(ji,j2+1,jk) = vn(ji,j2+1,jk) + zcor * vmask(ji,j2+1,jk)844 vv(ji,j2+1,jk,Kmm_a) = vv(ji,j2+1,jk,Kmm_a) + zcor * vmask(ji,j2+1,jk) 838 845 END DO 839 846 END DO … … 862 869 DO jj=j1,j2 863 870 DO ji=i1,i2 864 tabres(ji,jj) = zrhoy * u n_b(ji,jj) * hu_n(ji,jj) * e2u(ji,jj)871 tabres(ji,jj) = zrhoy * uu_b(ji,jj,Kmm_a) * hu(ji,jj,Kmm_a) * e2u(ji,jj) 865 872 END DO 866 873 END DO … … 873 880 spgu(ji,jj) = 0._wp 874 881 DO jk=1,jpkm1 875 spgu(ji,jj) = spgu(ji,jj) + e3u _n(ji,jj,jk) * un(ji,jj,jk)882 spgu(ji,jj) = spgu(ji,jj) + e3u(ji,jj,jk,Kmm_a) * uu(ji,jj,jk,Kmm_a) 876 883 END DO 877 884 ! 878 zcorr = (tabres(ji,jj) - spgu(ji,jj)) * r1_hu _n(ji,jj)885 zcorr = (tabres(ji,jj) - spgu(ji,jj)) * r1_hu(ji,jj,Kmm_a) 879 886 DO jk=1,jpkm1 880 u n(ji,jj,jk) = un(ji,jj,jk) + zcorr * umask(ji,jj,jk)887 uu(ji,jj,jk,Kmm_a) = uu(ji,jj,jk,Kmm_a) + zcorr * umask(ji,jj,jk) 881 888 END DO 882 889 ! 883 890 ! Update barotropic velocities: 884 891 IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN 885 IF ( .NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN! Add asselin part886 zcorr = (tabres(ji,jj) - u n_b(ji,jj) * hu_a(ji,jj)) * r1_hu_b(ji,jj)887 u b_b(ji,jj) = ub_b(ji,jj) + atfp * zcorr * umask(ji,jj,1)892 IF ( .NOT.( lk_agrif_fstep .AND. (neuler==0) ) ) THEN ! Add asselin part 893 zcorr = (tabres(ji,jj) - uu_b(ji,jj,Kmm_a) * hu(ji,jj,Krhs_a)) * r1_hu(ji,jj,Kbb_a) 894 uu_b(ji,jj,Kbb_a) = uu_b(ji,jj,Kbb_a) + atfp * zcorr * umask(ji,jj,1) 888 895 END IF 889 896 ENDIF 890 u n_b(ji,jj) = tabres(ji,jj) * r1_hu_n(ji,jj) * umask(ji,jj,1)897 uu_b(ji,jj,Kmm_a) = tabres(ji,jj) * r1_hu(ji,jj,Kmm_a) * umask(ji,jj,1) 891 898 ! 892 899 ! Correct "before" velocities to hold correct bt component: 893 900 spgu(ji,jj) = 0.e0 894 901 DO jk=1,jpkm1 895 spgu(ji,jj) = spgu(ji,jj) + e3u _b(ji,jj,jk) * ub(ji,jj,jk)902 spgu(ji,jj) = spgu(ji,jj) + e3u(ji,jj,jk,Kbb_a) * uu(ji,jj,jk,Kbb_a) 896 903 END DO 897 904 ! 898 zcorr = u b_b(ji,jj) - spgu(ji,jj) * r1_hu_b(ji,jj)905 zcorr = uu_b(ji,jj,Kbb_a) - spgu(ji,jj) * r1_hu(ji,jj,Kbb_a) 899 906 DO jk=1,jpkm1 900 u b(ji,jj,jk) = ub(ji,jj,jk) + zcorr * umask(ji,jj,jk)907 uu(ji,jj,jk,Kbb_a) = uu(ji,jj,jk,Kbb_a) + zcorr * umask(ji,jj,jk) 901 908 END DO 902 909 ! … … 905 912 ! 906 913 IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 907 u b_b(i1:i2,j1:j2) = un_b(i1:i2,j1:j2)914 uu_b(i1:i2,j1:j2,Kbb_a) = uu_b(i1:i2,j1:j2,Kmm_a) 908 915 ENDIF 909 916 ENDIF … … 928 935 DO jj=j1,j2 929 936 DO ji=i1,i2 930 tabres(ji,jj) = zrhox * v n_b(ji,jj) * hv_n(ji,jj) * e1v(ji,jj)937 tabres(ji,jj) = zrhox * vv_b(ji,jj,Kmm_a) * hv(ji,jj,Kmm_a) * e1v(ji,jj) 931 938 END DO 932 939 END DO … … 939 946 spgv(ji,jj) = 0.e0 940 947 DO jk=1,jpkm1 941 spgv(ji,jj) = spgv(ji,jj) + e3v _n(ji,jj,jk) * vn(ji,jj,jk)948 spgv(ji,jj) = spgv(ji,jj) + e3v(ji,jj,jk,Kmm_a) * vv(ji,jj,jk,Kmm_a) 942 949 END DO 943 950 ! 944 zcorr = (tabres(ji,jj) - spgv(ji,jj)) * r1_hv _n(ji,jj)951 zcorr = (tabres(ji,jj) - spgv(ji,jj)) * r1_hv(ji,jj,Kmm_a) 945 952 DO jk=1,jpkm1 946 v n(ji,jj,jk) = vn(ji,jj,jk) + zcorr * vmask(ji,jj,jk)953 vv(ji,jj,jk,Kmm_a) = vv(ji,jj,jk,Kmm_a) + zcorr * vmask(ji,jj,jk) 947 954 END DO 948 955 ! 949 956 ! Update barotropic velocities: 950 IF ( .NOT.ln_dynspg_ts .OR. ( ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN951 IF ( .NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN! Add asselin part952 zcorr = (tabres(ji,jj) - v n_b(ji,jj) * hv_a(ji,jj)) * r1_hv_b(ji,jj)953 v b_b(ji,jj) = vb_b(ji,jj) + atfp * zcorr *vmask(ji,jj,1)957 IF ( .NOT.ln_dynspg_ts .OR. ( ln_dynspg_ts .AND. ( .NOT.ln_bt_fw ) ) ) THEN 958 IF ( .NOT. ( lk_agrif_fstep .AND. ( neuler==0 ) ) ) THEN ! Add asselin part 959 zcorr = (tabres(ji,jj) - vv_b(ji,jj,Kmm_a) * hv(ji,jj,Krhs_a)) * r1_hv(ji,jj,Kbb_a) 960 vv_b(ji,jj,Kbb_a) = vv_b(ji,jj,Kbb_a) + atfp * zcorr * vmask(ji,jj,1) 954 961 END IF 955 962 ENDIF 956 v n_b(ji,jj) = tabres(ji,jj) * r1_hv_n(ji,jj) * vmask(ji,jj,1)963 vv_b(ji,jj,Kmm_a) = tabres(ji,jj) * r1_hv(ji,jj,Kmm_a) * vmask(ji,jj,1) 957 964 ! 958 965 ! Correct "before" velocities to hold correct bt component: 959 966 spgv(ji,jj) = 0.e0 960 967 DO jk=1,jpkm1 961 spgv(ji,jj) = spgv(ji,jj) + e3v _b(ji,jj,jk) * vb(ji,jj,jk)968 spgv(ji,jj) = spgv(ji,jj) + e3v(ji,jj,jk,Kbb_a) * vv(ji,jj,jk,Kbb_a) 962 969 END DO 963 970 ! 964 zcorr = v b_b(ji,jj) - spgv(ji,jj) * r1_hv_b(ji,jj)971 zcorr = vv_b(ji,jj,Kbb_a) - spgv(ji,jj) * r1_hv(ji,jj,Kbb_a) 965 972 DO jk=1,jpkm1 966 v b(ji,jj,jk) = vb(ji,jj,jk) + zcorr * vmask(ji,jj,jk)973 vv(ji,jj,jk,Kbb_a) = vv(ji,jj,jk,Kbb_a) + zcorr * vmask(ji,jj,jk) 967 974 END DO 968 975 ! … … 971 978 ! 972 979 IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 973 v b_b(i1:i2,j1:j2) = vn_b(i1:i2,j1:j2)980 vv_b(i1:i2,j1:j2,Kbb_a) = vv_b(i1:i2,j1:j2,Kmm_a) 974 981 ENDIF 975 982 ! … … 993 1000 DO jj=j1,j2 994 1001 DO ji=i1,i2 995 tabres(ji,jj) = ssh n(ji,jj)1002 tabres(ji,jj) = ssh(ji,jj,Kmm_a) 996 1003 END DO 997 1004 END DO … … 1000 1007 DO jj=j1,j2 1001 1008 DO ji=i1,i2 1002 sshb(ji,jj) = sshb(ji,jj) & 1003 & + atfp * ( tabres(ji,jj) - sshn(ji,jj) ) * tmask(ji,jj,1) 1009 ssh(ji,jj,Kbb_a) = ssh(ji,jj,Kbb_a) & 1010 & + atfp * ( tabres(ji,jj) - ssh(ji,jj,Kmm_a) ) & 1011 & * tmask(ji,jj,1) 1004 1012 END DO 1005 1013 END DO … … 1008 1016 DO jj=j1,j2 1009 1017 DO ji=i1,i2 1010 ssh n(ji,jj) = tabres(ji,jj) * tmask(ji,jj,1)1018 ssh(ji,jj,Kmm_a) = tabres(ji,jj) * tmask(ji,jj,1) 1011 1019 END DO 1012 1020 END DO 1013 1021 ! 1014 1022 IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 1015 ssh b(i1:i2,j1:j2) = sshn(i1:i2,j1:j2)1023 ssh(i1:i2,j1:j2,Kbb_a) = ssh(i1:i2,j1:j2,Kmm_a) 1016 1024 ENDIF 1017 1025 ! … … 1094 1102 DO jj=j1,j2 1095 1103 zcor = rdt * r1_e1e2t(i1 ,jj) * e2u(i1,jj) * (ub2_b(i1,jj)-tabres(i1,jj)) 1096 sshn(i1 ,jj) = sshn(i1 ,jj) + zcor 1097 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(i1 ,jj) = sshb(i1 ,jj) + atfp * zcor 1104 ssh(i1 ,jj,Kmm_a) = ssh(i1 ,jj,Kmm_a) + zcor 1105 IF ( .NOT. ( lk_agrif_fstep .AND. ( neuler==0 ) ) ) & 1106 & ssh(i1 ,jj,Kbb_a) = ssh(i1 ,jj,Kbb_a) + atfp * zcor 1098 1107 END DO 1099 1108 ENDIF … … 1101 1110 DO jj=j1,j2 1102 1111 zcor = - rdt * r1_e1e2t(i2+1,jj) * e2u(i2,jj) * (ub2_b(i2,jj)-tabres(i2,jj)) 1103 sshn(i2+1,jj) = sshn(i2+1,jj) + zcor 1104 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(i2+1,jj) = sshb(i2+1,jj) + atfp * zcor 1112 ssh(i2+1,jj,Kmm_a) = ssh(i2+1,jj,Kmm_a) + zcor 1113 IF ( .NOT. ( lk_agrif_fstep .AND. ( neuler==0 ) ) ) & 1114 & ssh(i2+1,jj,Kbb_a) = ssh(i2+1,jj,Kbb_a) + atfp * zcor 1105 1115 END DO 1106 1116 ENDIF … … 1181 1191 IF (southern_side) THEN 1182 1192 DO ji=i1,i2 1183 zcor = rdt * r1_e1e2t(ji,j1 ) * e1v(ji,j1 ) * (vb2_b(ji,j1)-tabres(ji,j1)) 1184 sshn(ji,j1 ) = sshn(ji,j1 ) + zcor 1185 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(ji,j1 ) = sshb(ji,j1) + atfp * zcor 1193 zcor = rdt * r1_e1e2t(ji,j1 ) * e1v(ji,j1 ) * ( vb2_b(ji,j1)-tabres(ji,j1) ) 1194 ssh(ji,j1 ,Kmm_a) = ssh(ji,j1 ,Kmm_a) + zcor 1195 IF ( .NOT. ( lk_agrif_fstep .AND. ( neuler==0 ) ) ) & 1196 & ssh(ji,j1 ,Kbb_a) = ssh(ji,j1,Kbb_a) + atfp * zcor 1186 1197 END DO 1187 1198 ENDIF 1188 1199 IF (northern_side) THEN 1189 1200 DO ji=i1,i2 1190 zcor = - rdt * r1_e1e2t(ji,j2+1) * e1v(ji,j2 ) * (vb2_b(ji,j2)-tabres(ji,j2)) 1191 sshn(ji,j2+1) = sshn(ji,j2+1) + zcor 1192 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(ji,j2+1) = sshb(ji,j2+1) + atfp * zcor 1201 zcor = - rdt * r1_e1e2t(ji,j2+1) * e1v(ji,j2 ) * ( vb2_b(ji,j2)-tabres(ji,j2) ) 1202 ssh(ji,j2+1,Kmm_a) = ssh(ji,j2+1,Kmm_a) + zcor 1203 IF ( .NOT. ( lk_agrif_fstep .AND. ( neuler==0 ) ) ) & 1204 & ssh(ji,j2+1,Kbb_a) = ssh(ji,j2+1,Kbb_a) + atfp * zcor 1193 1205 END DO 1194 1206 ENDIF … … 1223 1235 END DO 1224 1236 END DO 1225 tabres(:,:,:,1) =tabres(:,:,:,1)*Agrif_Rhox()*Agrif_Rhoy()1226 tabres(:,:,:,2) =tabres(:,:,:,2)*Agrif_Rhox()1227 tabres(:,:,:,3) =tabres(:,:,:,3)*Agrif_Rhoy()1237 tabres(:,:,:,1) = tabres(:,:,:,1)*Agrif_Rhox()*Agrif_Rhoy() 1238 tabres(:,:,:,2) = tabres(:,:,:,2)*Agrif_Rhox() 1239 tabres(:,:,:,3) = tabres(:,:,:,3)*Agrif_Rhoy() 1228 1240 ELSE 1229 1241 DO jk=k1,k2 … … 1234 1246 print *,'VAL2 = ',ji,jj,jk,tabres(ji,jj,jk,2),e1t(ji,jj)*tmask(ji,jj,jk) 1235 1247 print *,'VAL3 = ',ji,jj,jk,tabres(ji,jj,jk,3),e2t(ji,jj)*tmask(ji,jj,jk) 1236 ztemp = sqrt(tabres(ji,jj,jk,1)/(tabres(ji,jj,jk,2)*tabres(ji,jj,jk,3)))1248 ztemp = SQRT( tabres(ji,jj,jk,1)/(tabres(ji,jj,jk,2)*tabres(ji,jj,jk,3)) ) 1237 1249 print *,'CORR = ',ztemp-1. 1238 1250 print *,'NEW VALS = ',tabres(ji,jj,jk,2)*ztemp,tabres(ji,jj,jk,3)*ztemp, & 1239 1251 tabres(ji,jj,jk,2)*ztemp*tabres(ji,jj,jk,3)*ztemp 1240 e1t(ji,jj) = tabres(ji,jj,jk,2) *ztemp1241 e2t(ji,jj) = tabres(ji,jj,jk,3) *ztemp1252 e1t(ji,jj) = tabres(ji,jj,jk,2) * ztemp 1253 e2t(ji,jj) = tabres(ji,jj,jk,3) * ztemp 1242 1254 END IF 1243 1255 END DO … … 1319 1331 DO jj=j1,j2 1320 1332 DO ji=i1,i2 1321 ptab(ji,jj,jk) = e3t_0(ji,jj,jk) * (1._wp + sshn(ji,jj) & 1322 & *ssmask(ji,jj)/(ht_0(ji,jj)-1._wp + ssmask(ji,jj))) 1333 ptab(ji,jj,jk) = e3t_0(ji,jj,jk) & 1334 & * ( 1._wp + ssh(ji,jj,Kmm_a) & 1335 & * ssmask(ji,jj) & 1336 & / ( ht_0(ji,jj)-1._wp + ssmask(ji,jj) ) ) 1323 1337 END DO 1324 1338 END DO … … 1330 1344 ! Save "old" scale factor (prior update) for subsequent asselin correction 1331 1345 ! of prognostic variables 1332 e3t _a(i1:i2,j1:j2,1:jpkm1) = e3t_n(i1:i2,j1:j2,1:jpkm1)1333 1334 ! One should also save e3t _b, but lacking of workspace...1335 ! hdiv n(i1:i2,j1:j2,1:jpkm1) = e3t_b(i1:i2,j1:j2,1:jpkm1)1346 e3t(i1:i2,j1:j2,1:jpkm1,Krhs_a) = e3t(i1:i2,j1:j2,1:jpkm1,Kmm_a) 1347 1348 ! One should also save e3t(:,:,:,Kbb_a), but lacking of workspace... 1349 ! hdiv(i1:i2,j1:j2,1:jpkm1) = e3t(i1:i2,j1:j2,1:jpkm1,Kbb_a) 1336 1350 1337 1351 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0) )) THEN … … 1339 1353 DO jj=j1,j2 1340 1354 DO ji=i1,i2 1341 e3t _b(ji,jj,jk) = e3t_b(ji,jj,jk)&1342 & + atfp * ( ptab(ji,jj,jk) - e3t_n(ji,jj,jk) )1355 e3t(ji,jj,jk,Kbb_a) = e3t(ji,jj,jk,Kbb_a) & 1356 & + atfp * ( ptab(ji,jj,jk) - e3t(ji,jj,jk,Kmm_a) ) 1343 1357 END DO 1344 1358 END DO 1345 1359 END DO 1346 1360 ! 1347 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)1348 gdepw _b(i1:i2,j1:j2,1) = 0.0_wp1349 gdept _b(i1:i2,j1:j2,1) = 0.5_wp * e3w_b(i1:i2,j1:j2,1)1361 e3w (i1:i2,j1:j2,1,Kbb_a) = e3w_0(i1:i2,j1:j2,1) + e3t(i1:i2,j1:j2,1,Kbb_a) - e3t_0(i1:i2,j1:j2,1) 1362 gdepw(i1:i2,j1:j2,1,Kbb_a) = 0.0_wp 1363 gdept(i1:i2,j1:j2,1,Kbb_a) = 0.5_wp * e3w(i1:i2,j1:j2,1,Kbb_a) 1350 1364 ! 1351 1365 DO jk = 2, jpk … … 1353 1367 DO ji = i1,i2 1354 1368 zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk)) 1355 e3w _b(ji,jj,jk) = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) *&1356 & ( e3t_b(ji,jj,jk-1) - e3t_0(ji,jj,jk-1) ) &1357 & + 0.5_wp * tmask(ji,jj,jk) *&1358 & ( e3t_b(ji,jj,jk) - e3t_0(ji,jj,jk ) )1359 gdepw _b(ji,jj,jk) = gdepw_b(ji,jj,jk-1) + e3t_b(ji,jj,jk-1)1360 gdept _b(ji,jj,jk) = zcoef * ( gdepw_b(ji,jj,jk ) + 0.5 * e3w_b(ji,jj,jk)) &1361 & + (1-zcoef) * ( gdept_b(ji,jj,jk-1) + e3w_b(ji,jj,jk))1369 e3w(ji,jj,jk,Kbb_a) = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) & 1370 & * ( e3t(ji,jj,jk-1,Kbb_a) - e3t_0(ji,jj,jk-1) ) & 1371 & + 0.5_wp * tmask(ji,jj,jk) & 1372 & * ( e3t(ji,jj,jk ,Kbb_a) - e3t_0(ji,jj,jk ) ) 1373 gdepw(ji,jj,jk,Kbb_a) = gdepw(ji,jj,jk-1,Kbb_a) + e3t(ji,jj,jk-1,Kbb_a) 1374 gdept(ji,jj,jk,Kbb_a) = zcoef * ( gdepw(ji,jj,jk ,Kbb_a) + 0.5 * e3w(ji,jj,jk ,Kbb_a) ) & 1375 & + (1._wp - zcoef) * ( gdept(ji,jj,jk-1,Kbb_a) + e3w(ji,jj,jk ,Kbb_a) ) 1362 1376 END DO 1363 1377 END DO … … 1370 1384 ! 1371 1385 ! Update vertical scale factor at T-points: 1372 e3t _n(i1:i2,j1:j2,1:jpkm1) = ptab(i1:i2,j1:j2,1:jpkm1)1386 e3t(i1:i2,j1:j2,1:jpkm1,Kmm_a) = ptab(i1:i2,j1:j2,1:jpkm1) 1373 1387 ! 1374 1388 ! Update total depth: 1375 ht _n(i1:i2,j1:j2) = 0._wp1389 ht(i1:i2,j1:j2) = 0._wp 1376 1390 DO jk = 1, jpkm1 1377 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)1391 ht(i1:i2,j1:j2) = ht(i1:i2,j1:j2) + e3t(i1:i2,j1:j2,jk,Kmm_a) * tmask(i1:i2,j1:j2,jk) 1378 1392 END DO 1379 1393 ! 1380 1394 ! Update vertical scale factor at W-points and depths: 1381 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)1382 gdept _n(i1:i2,j1:j2,1) = 0.5_wp * e3w_n(i1:i2,j1:j2,1)1383 gdepw _n(i1:i2,j1:j2,1) = 0.0_wp1384 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 ssh1395 e3w(i1:i2,j1:j2,1,Kmm_a) = e3w_0(i1:i2,j1:j2,1) + e3t(i1:i2,j1:j2,1,Kmm_a) - e3t_0(i1:i2,j1:j2,1) 1396 gdept(i1:i2,j1:j2,1,Kmm_a) = 0.5_wp * e3w(i1:i2,j1:j2,1,Kmm_a) 1397 gdepw(i1:i2,j1:j2,1,Kmm_a) = 0.0_wp 1398 gde3w(i1:i2,j1:j2,1 ) = gdept(i1:i2,j1:j2,1,Kmm_a) - ( ht(i1:i2,j1:j2) - ht_0(i1:i2,j1:j2) ) ! Last term in the rhs is ssh 1385 1399 ! 1386 1400 DO jk = 2, jpk … … 1388 1402 DO ji = i1,i2 1389 1403 zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk)) 1390 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) ) & 1391 & + 0.5_wp * tmask(ji,jj,jk) * ( e3t_n(ji,jj,jk ) - e3t_0(ji,jj,jk ) ) 1392 gdepw_n(ji,jj,jk) = gdepw_n(ji,jj,jk-1) + e3t_n(ji,jj,jk-1) 1393 gdept_n(ji,jj,jk) = zcoef * ( gdepw_n(ji,jj,jk ) + 0.5 * e3w_n(ji,jj,jk)) & 1394 & + (1-zcoef) * ( gdept_n(ji,jj,jk-1) + e3w_n(ji,jj,jk)) 1395 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 1404 e3w(ji,jj,jk,Kmm_a) = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) & 1405 & * ( e3t(ji,jj,jk-1,Kmm_a) - e3t_0(ji,jj,jk-1) ) & 1406 & + 0.5_wp * tmask(ji,jj,jk) & 1407 & * ( e3t(ji,jj,jk ,Kmm_a) - e3t_0(ji,jj,jk ) ) 1408 gdepw(ji,jj,jk,Kmm_a) = gdepw(ji,jj,jk-1,Kmm_a) + e3t(ji,jj,jk-1,Kmm_a) 1409 gdept(ji,jj,jk,Kmm_a) = zcoef * ( gdepw(ji,jj,jk ,Kmm_a) + 0.5 * e3w(ji,jj,jk ,Kmm_a) ) & 1410 & + ( 1._wp - zcoef ) * ( gdept(ji,jj,jk-1,Kmm_a) + e3w(ji,jj,jk ,Kmm_a) ) 1411 gde3w(ji,jj,jk ) = gdept(ji,jj,jk ,Kmm_a) - ( ht(ji,jj)-ht_0(ji,jj) ) ! Last term in the rhs is ssh 1396 1412 END DO 1397 1413 END DO … … 1399 1415 ! 1400 1416 IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 1401 e3t_b (i1:i2,j1:j2,1:jpk) = e3t_n (i1:i2,j1:j2,1:jpk)1402 e3w_b (i1:i2,j1:j2,1:jpk) = e3w_n (i1:i2,j1:j2,1:jpk)1403 gdepw _b(i1:i2,j1:j2,1:jpk) = gdepw_n(i1:i2,j1:j2,1:jpk)1404 gdept _b(i1:i2,j1:j2,1:jpk) = gdept_n(i1:i2,j1:j2,1:jpk)1417 e3t(i1:i2,j1:j2,1:jpk,Kbb_a) = e3t(i1:i2,j1:j2,1:jpk,Kmm_a) 1418 e3w(i1:i2,j1:j2,1:jpk,Kbb_a) = e3w(i1:i2,j1:j2,1:jpk,Kmm_a) 1419 gdepw(i1:i2,j1:j2,1:jpk,Kbb_a) = gdepw(i1:i2,j1:j2,1:jpk,Kmm_a) 1420 gdept(i1:i2,j1:j2,1:jpk,Kbb_a) = gdept(i1:i2,j1:j2,1:jpk,Kmm_a) 1405 1421 ENDIF 1406 1422 !
Note: See TracChangeset
for help on using the changeset viewer.
