Changeset 10989 for NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/NST/agrif_oce_update.F90
- Timestamp:
- 2019-05-16T17:45:46+02:00 (5 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/NST/agrif_oce_update.F90
r10068 r10989 230 230 ! ----------------------- 231 231 ! 232 e3u _a(:,:,:) = e3u_n(:,:,:)233 e3v _a(:,:,:) = e3v_n(:,:,:)234 ! u a(:,:,:) = e3u_b(:,:,:)235 ! v a(:,:,:) = e3v_b(:,:,:)232 e3u(:,:,:,Krhs) = e3u(:,:,:,Kmm) 233 e3v(:,:,:,Krhs) = e3v(:,:,:,Kmm) 234 ! uu(:,:,:,Krhs) = e3u(:,:,:,Kbb) 235 ! vv(:,:,:,Krhs) = e3v(:,:,:,Kbb) 236 236 hu_a(:,:) = hu_n(:,:) 237 237 hv_a(:,:) = hv_n(:,:) … … 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), e3u(:,:,:,Kmm) , 'U' ) 245 CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3v(:,:,:,Kmm) , 'V' ) 246 CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3f(:,:,:) , 'F' ) 247 248 CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3uw(:,:,:,Kmm), 'UW' ) 249 CALL dom_vvl_interpol( e3v(:,:,:,Kmm), e3vw(:,:,:,Kmm), 'VW' ) 250 250 251 251 ! Update total depths: … … 254 254 hv_n(:,:) = 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_n(:,:) = hu_n(:,:) + e3u(:,:,jk,Kmm) * umask(:,:,jk) 257 hv_n(:,:) = hv_n(:,:) + e3v(:,:,jk,Kmm) * vmask(:,:,jk) 258 258 END DO 259 259 ! ! Inverse of the local depth … … 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), e3u(:,:,:,Kbb), 'U' ) 271 CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3v(:,:,:,Kbb), 'V' ) 272 273 CALL dom_vvl_interpol( e3u(:,:,:,Kbb), e3uw(:,:,:,Kbb), 'UW' ) 274 CALL dom_vvl_interpol( e3v(:,:,:,Kbb), e3vw(:,:,:,Kbb), 'VW' ) 275 275 276 276 ! Update total depths: … … 279 279 hv_b(:,:) = 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_b(:,:) = hu_b(:,:) + e3u(:,:,jk,Kbb) * umask(:,:,jk) 282 hv_b(:,:) = hv_b(:,:) + e3v(:,:,jk,Kbb) * vmask(:,:,jk) 283 283 END DO 284 284 ! ! Inverse of the local depth … … 315 315 DO jj=j1,j2 316 316 DO ji=i1,i2 317 tabres(ji,jj,jk,jn) = (ts n(ji,jj,jk,jn) * e3t_n(ji,jj,jk) ) &317 tabres(ji,jj,jk,jn) = (ts(ji,jj,jk,jn,Kmm) * e3t(ji,jj,jk,Kmm) ) & 318 318 * tmask(ji,jj,jk) + (tmask(ji,jj,jk)-1)*999._wp 319 319 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) &326 tabres(ji,jj,jk,n2) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm) & 327 327 + (tmask(ji,jj,jk)-1)*999._wp 328 328 END DO … … 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) 348 348 ENDDO 349 349 IF (N_in > 0) THEN !Remove this? … … 369 369 DO ji=i1,i2 370 370 IF( tabres_child(ji,jj,jk,jn) .NE. 0. ) THEN 371 ts b(ji,jj,jk,jn) = tsb(ji,jj,jk,jn) &371 ts(ji,jj,jk,jn,Kbb) = ts(ji,jj,jk,jn,Kbb) & 372 372 & + atfp * ( tabres_child(ji,jj,jk,jn) & 373 & - ts n(ji,jj,jk,jn) ) * tmask(ji,jj,jk)373 & - ts(ji,jj,jk,jn,Kmm) ) * tmask(ji,jj,jk) 374 374 ENDIF 375 375 ENDDO … … 383 383 DO ji=i1,i2 384 384 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)385 ts(ji,jj,jk,jn,Kmm) = tabres_child(ji,jj,jk,jn) * tmask(ji,jj,jk) 386 386 END IF 387 387 END DO … … 413 413 DO ji=i1,i2 414 414 !> 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) = ts n(ji,jj,jk,jn) * e3t_n(ji,jj,jk)415 tabres(ji,jj,jk,jn) = ts(ji,jj,jk,jn,Kmm) * e3t(ji,jj,jk,Kmm) / e3t_0(ji,jj,jk) 416 ! tabres(ji,jj,jk,jn) = ts(ji,jj,jk,jn,Kmm) * e3t(ji,jj,jk,Kmm) 417 417 !< jc tmp 418 418 END DO … … 434 434 DO ji = i1, i2 435 435 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 used436 ztb = ts(ji,jj,jk,jn,Kbb) * e3t(ji,jj,jk,Kbb) ! fse3t_b prior update should be used 437 437 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)438 ztno = ts(ji,jj,jk,jn,Kmm) * e3t(ji,jj,jk,Krhs) 439 ts(ji,jj,jk,jn,Kbb) = ( ztb + atfp * ( ztnu - ztno) ) & 440 & * tmask(ji,jj,jk) / e3t(ji,jj,jk,Kbb) 441 441 ENDIF 442 442 END DO … … 450 450 DO ji=i1,i2 451 451 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)452 ts(ji,jj,jk,jn,Kmm) = tabres(ji,jj,jk,jn) / e3t(ji,jj,jk,Kmm) 453 453 END IF 454 454 END DO … … 458 458 ! 459 459 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)460 ts(i1:i2,j1:j2,k1:k2,1:jpts,Kbb) = ts(i1:i2,j1:j2,k1:k2,1:jpts,Kmm) 461 461 ENDIF 462 462 ! … … 495 495 DO jj=j1,j2 496 496 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) &497 tabres(ji,jj,jk,1) = zrhoy * e2u(ji,jj) * e3u(ji,jj,jk,Kmm) * umask(ji,jj,jk) * uu(ji,jj,jk,Kmm) & 498 498 + (umask(ji,jj,jk)-1)*999._wp 499 tabres(ji,jj,jk,2) = zrhoy * umask(ji,jj,jk) * e2u(ji,jj) * e3u _n(ji,jj,jk) &499 tabres(ji,jj,jk,2) = zrhoy * umask(ji,jj,jk) * e2u(ji,jj) * e3u(ji,jj,jk,Kmm) & 500 500 + (umask(ji,jj,jk)-1)*999._wp 501 501 END DO … … 520 520 IF (umask(ji,jj,jk) == 0) EXIT 521 521 N_out = N_out + 1 522 h_out(N_out) = e3u _n(ji,jj,jk)522 h_out(N_out) = e3u(ji,jj,jk,Kmm) 523 523 ENDDO 524 524 IF (N_in * N_out > 0) THEN … … 550 550 DO ji=i1,i2 551 551 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) - u n(ji,jj,jk) ) * umask(ji,jj,jk)552 uu(ji,jj,jk,Kbb) = uu(ji,jj,jk,Kbb) & 553 & + atfp * ( tabres_child(ji,jj,jk) - uu(ji,jj,jk,Kmm) ) * umask(ji,jj,jk) 554 554 ENDIF 555 555 ! 556 u n(ji,jj,jk) = tabres_child(ji,jj,jk) * umask(ji,jj,jk)556 uu(ji,jj,jk,Kmm) = tabres_child(ji,jj,jk) * umask(ji,jj,jk) 557 557 END DO 558 558 END DO … … 579 579 zrhoy = Agrif_Rhoy() 580 580 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)581 tabres(i1:i2,j1:j2,jk,1) = zrhoy * e2u(i1:i2,j1:j2) * e3u(i1:i2,j1:j2,jk,Kmm) * uu(i1:i2,j1:j2,jk,Kmm) 582 582 END DO 583 583 ELSE … … 588 588 ! 589 589 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)590 zub = uu(ji,jj,jk,Kbb) * e3u(ji,jj,jk,Kbb) ! fse3t_b prior update should be used 591 zuno = uu(ji,jj,jk,Kmm) * e3u(ji,jj,jk,Krhs) 592 592 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)593 uu(ji,jj,jk,Kbb) = ( zub + atfp * ( zunu - zuno) ) & 594 & * umask(ji,jj,jk) / e3u(ji,jj,jk,Kbb) 595 595 ENDIF 596 596 ! 597 u n(ji,jj,jk) = tabres(ji,jj,jk,1) * umask(ji,jj,jk) / e3u_n(ji,jj,jk)597 uu(ji,jj,jk,Kmm) = tabres(ji,jj,jk,1) * umask(ji,jj,jk) / e3u(ji,jj,jk,Kmm) 598 598 END DO 599 599 END DO … … 601 601 ! 602 602 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)603 uu(i1:i2,j1:j2,k1:k2,Kbb) = uu(i1:i2,j1:j2,k1:k2,Kmm) 604 604 ENDIF 605 605 ! … … 632 632 IF (western_side) THEN 633 633 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) + zcor634 zcor = uu_b(i1-1,jj,Kmm) * hu_a(i1-1,jj) * r1_hu_n(i1-1,jj) - uu_b(i1-1,jj,Kmm) 635 uu_b(i1-1,jj,Kmm) = uu_b(i1-1,jj,Kmm) + zcor 636 636 DO jk=1,jpkm1 637 u n(i1-1,jj,jk) = un(i1-1,jj,jk) + zcor * umask(i1-1,jj,jk)637 uu(i1-1,jj,jk,Kmm) = uu(i1-1,jj,jk,Kmm) + zcor * umask(i1-1,jj,jk) 638 638 END DO 639 639 END DO … … 642 642 IF (eastern_side) THEN 643 643 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) + zcor644 zcor = uu_b(i2+1,jj,Kmm) * hu_a(i2+1,jj) * r1_hu_n(i2+1,jj) - uu_b(i2+1,jj,Kmm) 645 uu_b(i2+1,jj,Kmm) = uu_b(i2+1,jj,Kmm) + zcor 646 646 DO jk=1,jpkm1 647 u n(i2+1,jj,jk) = un(i2+1,jj,jk) + zcor * umask(i2+1,jj,jk)647 uu(i2+1,jj,jk,Kmm) = uu(i2+1,jj,jk,Kmm) + zcor * umask(i2+1,jj,jk) 648 648 END DO 649 649 END DO … … 682 682 DO jj=j1,j2 683 683 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) &684 tabres(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v(ji,jj,jk,Kmm) * vmask(ji,jj,jk) * vv(ji,jj,jk,Kmm) & 685 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) &686 tabres(ji,jj,jk,2) = vmask(ji,jj,jk) * zrhox * e1v(ji,jj) * e3v(ji,jj,jk,Kmm) & 687 687 + (vmask(ji,jj,jk)-1)*999._wp 688 688 END DO … … 705 705 IF (vmask(ji,jj,jk) == 0) EXIT 706 706 N_out = N_out + 1 707 h_out(N_out) = e3v _n(ji,jj,jk)707 h_out(N_out) = e3v(ji,jj,jk,Kmm) 708 708 ENDDO 709 709 IF (N_in * N_out > 0) THEN … … 736 736 ! 737 737 IF( .NOT.(lk_agrif_fstep.AND.(neuler==0)) ) THEN ! Add asselin part 738 v b(ji,jj,jk) = vb(ji,jj,jk) &739 & + atfp * ( tabres_child(ji,jj,jk) - v n(ji,jj,jk) ) * vmask(ji,jj,jk)738 vv(ji,jj,jk,Kbb) = vv(ji,jj,jk,Kbb) & 739 & + atfp * ( tabres_child(ji,jj,jk) - vv(ji,jj,jk,Kmm) ) * vmask(ji,jj,jk) 740 740 ENDIF 741 741 ! 742 v n(ji,jj,jk) = tabres_child(ji,jj,jk) * vmask(ji,jj,jk)742 vv(ji,jj,jk,Kmm) = tabres_child(ji,jj,jk) * vmask(ji,jj,jk) 743 743 END DO 744 744 END DO … … 767 767 DO jj=j1,j2 768 768 DO ji=i1,i2 769 tabres(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v _n(ji,jj,jk) * vn(ji,jj,jk)769 tabres(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v(ji,jj,jk,Kmm) * vv(ji,jj,jk,Kmm) 770 770 END DO 771 771 END DO … … 778 778 ! 779 779 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)780 zvb = vv(ji,jj,jk,Kbb) * e3v(ji,jj,jk,Kbb) ! fse3t_b prior update should be used 781 zvno = vv(ji,jj,jk,Kmm) * e3v(ji,jj,jk,Krhs) 782 782 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)783 vv(ji,jj,jk,Kbb) = ( zvb + atfp * ( zvnu - zvno) ) & 784 & * vmask(ji,jj,jk) / e3v(ji,jj,jk,Kbb) 785 785 ENDIF 786 786 ! 787 v n(ji,jj,jk) = tabres(ji,jj,jk,1) * vmask(ji,jj,jk) / e3v_n(ji,jj,jk)787 vv(ji,jj,jk,Kmm) = tabres(ji,jj,jk,1) * vmask(ji,jj,jk) / e3v(ji,jj,jk,Kmm) 788 788 END DO 789 789 END DO … … 791 791 ! 792 792 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)793 vv(i1:i2,j1:j2,k1:k2,Kbb) = vv(i1:i2,j1:j2,k1:k2,Kmm) 794 794 ENDIF 795 795 ! … … 822 822 IF (southern_side) THEN 823 823 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) + zcor824 zcor = vv_b(ji,j1-1,Kmm) * hv_a(ji,j1-1) * r1_hv_n(ji,j1-1) - vv_b(ji,j1-1,Kmm) 825 vv_b(ji,j1-1,Kmm) = vv_b(ji,j1-1,Kmm) + zcor 826 826 DO jk=1,jpkm1 827 v n(ji,j1-1,jk) = vn(ji,j1-1,jk) + zcor * vmask(ji,j1-1,jk)827 vv(ji,j1-1,jk,Kmm) = vv(ji,j1-1,jk,Kmm) + zcor * vmask(ji,j1-1,jk) 828 828 END DO 829 829 END DO … … 832 832 IF (northern_side) THEN 833 833 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) + zcor834 zcor = vv_b(ji,j2+1,Kmm) * hv_a(ji,j2+1) * r1_hv_n(ji,j2+1) - vv_b(ji,j2+1,Kmm) 835 vv_b(ji,j2+1,Kmm) = vv_b(ji,j2+1,Kmm) + zcor 836 836 DO jk=1,jpkm1 837 v n(ji,j2+1,jk) = vn(ji,j2+1,jk) + zcor * vmask(ji,j2+1,jk)837 vv(ji,j2+1,jk,Kmm) = vv(ji,j2+1,jk,Kmm) + zcor * vmask(ji,j2+1,jk) 838 838 END DO 839 839 END DO … … 862 862 DO jj=j1,j2 863 863 DO ji=i1,i2 864 tabres(ji,jj) = zrhoy * u n_b(ji,jj) * hu_n(ji,jj) * e2u(ji,jj)864 tabres(ji,jj) = zrhoy * uu_b(ji,jj,Kmm) * hu_n(ji,jj) * e2u(ji,jj) 865 865 END DO 866 866 END DO … … 873 873 spgu(ji,jj) = 0._wp 874 874 DO jk=1,jpkm1 875 spgu(ji,jj) = spgu(ji,jj) + e3u _n(ji,jj,jk) * un(ji,jj,jk)875 spgu(ji,jj) = spgu(ji,jj) + e3u(ji,jj,jk,Kmm) * uu(ji,jj,jk,Kmm) 876 876 END DO 877 877 ! 878 878 zcorr = (tabres(ji,jj) - spgu(ji,jj)) * r1_hu_n(ji,jj) 879 879 DO jk=1,jpkm1 880 u n(ji,jj,jk) = un(ji,jj,jk) + zcorr * umask(ji,jj,jk)880 uu(ji,jj,jk,Kmm) = uu(ji,jj,jk,Kmm) + zcorr * umask(ji,jj,jk) 881 881 END DO 882 882 ! … … 884 884 IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN 885 885 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 886 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)886 zcorr = (tabres(ji,jj) - uu_b(ji,jj,Kmm) * hu_a(ji,jj)) * r1_hu_b(ji,jj) 887 uu_b(ji,jj,Kbb) = uu_b(ji,jj,Kbb) + atfp * zcorr * umask(ji,jj,1) 888 888 END IF 889 889 ENDIF 890 u n_b(ji,jj) = tabres(ji,jj) * r1_hu_n(ji,jj) * umask(ji,jj,1)890 uu_b(ji,jj,Kmm) = tabres(ji,jj) * r1_hu_n(ji,jj) * umask(ji,jj,1) 891 891 ! 892 892 ! Correct "before" velocities to hold correct bt component: 893 893 spgu(ji,jj) = 0.e0 894 894 DO jk=1,jpkm1 895 spgu(ji,jj) = spgu(ji,jj) + e3u _b(ji,jj,jk) * ub(ji,jj,jk)895 spgu(ji,jj) = spgu(ji,jj) + e3u(ji,jj,jk,Kbb) * uu(ji,jj,jk,Kbb) 896 896 END DO 897 897 ! 898 zcorr = u b_b(ji,jj) - spgu(ji,jj) * r1_hu_b(ji,jj)898 zcorr = uu_b(ji,jj,Kbb) - spgu(ji,jj) * r1_hu_b(ji,jj) 899 899 DO jk=1,jpkm1 900 u b(ji,jj,jk) = ub(ji,jj,jk) + zcorr * umask(ji,jj,jk)900 uu(ji,jj,jk,Kbb) = uu(ji,jj,jk,Kbb) + zcorr * umask(ji,jj,jk) 901 901 END DO 902 902 ! … … 905 905 ! 906 906 IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 907 u b_b(i1:i2,j1:j2) = un_b(i1:i2,j1:j2)907 uu_b(i1:i2,j1:j2,Kbb) = uu_b(i1:i2,j1:j2,Kmm) 908 908 ENDIF 909 909 ENDIF … … 928 928 DO jj=j1,j2 929 929 DO ji=i1,i2 930 tabres(ji,jj) = zrhox * v n_b(ji,jj) * hv_n(ji,jj) * e1v(ji,jj)930 tabres(ji,jj) = zrhox * vv_b(ji,jj,Kmm) * hv_n(ji,jj) * e1v(ji,jj) 931 931 END DO 932 932 END DO … … 939 939 spgv(ji,jj) = 0.e0 940 940 DO jk=1,jpkm1 941 spgv(ji,jj) = spgv(ji,jj) + e3v _n(ji,jj,jk) * vn(ji,jj,jk)941 spgv(ji,jj) = spgv(ji,jj) + e3v(ji,jj,jk,Kmm) * vv(ji,jj,jk,Kmm) 942 942 END DO 943 943 ! 944 944 zcorr = (tabres(ji,jj) - spgv(ji,jj)) * r1_hv_n(ji,jj) 945 945 DO jk=1,jpkm1 946 v n(ji,jj,jk) = vn(ji,jj,jk) + zcorr * vmask(ji,jj,jk)946 vv(ji,jj,jk,Kmm) = vv(ji,jj,jk,Kmm) + zcorr * vmask(ji,jj,jk) 947 947 END DO 948 948 ! … … 950 950 IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN 951 951 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 952 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)952 zcorr = (tabres(ji,jj) - vv_b(ji,jj,Kmm) * hv_a(ji,jj)) * r1_hv_b(ji,jj) 953 vv_b(ji,jj,Kbb) = vv_b(ji,jj,Kbb) + atfp * zcorr * vmask(ji,jj,1) 954 954 END IF 955 955 ENDIF 956 v n_b(ji,jj) = tabres(ji,jj) * r1_hv_n(ji,jj) * vmask(ji,jj,1)956 vv_b(ji,jj,Kmm) = tabres(ji,jj) * r1_hv_n(ji,jj) * vmask(ji,jj,1) 957 957 ! 958 958 ! Correct "before" velocities to hold correct bt component: 959 959 spgv(ji,jj) = 0.e0 960 960 DO jk=1,jpkm1 961 spgv(ji,jj) = spgv(ji,jj) + e3v _b(ji,jj,jk) * vb(ji,jj,jk)961 spgv(ji,jj) = spgv(ji,jj) + e3v(ji,jj,jk,Kbb) * vv(ji,jj,jk,Kbb) 962 962 END DO 963 963 ! 964 zcorr = v b_b(ji,jj) - spgv(ji,jj) * r1_hv_b(ji,jj)964 zcorr = vv_b(ji,jj,Kbb) - spgv(ji,jj) * r1_hv_b(ji,jj) 965 965 DO jk=1,jpkm1 966 v b(ji,jj,jk) = vb(ji,jj,jk) + zcorr * vmask(ji,jj,jk)966 vv(ji,jj,jk,Kbb) = vv(ji,jj,jk,Kbb) + zcorr * vmask(ji,jj,jk) 967 967 END DO 968 968 ! … … 971 971 ! 972 972 IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 973 v b_b(i1:i2,j1:j2) = vn_b(i1:i2,j1:j2)973 vv_b(i1:i2,j1:j2,Kbb) = vv_b(i1:i2,j1:j2,Kmm) 974 974 ENDIF 975 975 ! … … 993 993 DO jj=j1,j2 994 994 DO ji=i1,i2 995 tabres(ji,jj) = ssh n(ji,jj)995 tabres(ji,jj) = ssh(ji,jj,Kmm) 996 996 END DO 997 997 END DO … … 1000 1000 DO jj=j1,j2 1001 1001 DO ji=i1,i2 1002 ssh b(ji,jj) = sshb(ji,jj) &1003 & + atfp * ( tabres(ji,jj) - ssh n(ji,jj) ) * tmask(ji,jj,1)1002 ssh(ji,jj,Kbb) = ssh(ji,jj,Kbb) & 1003 & + atfp * ( tabres(ji,jj) - ssh(ji,jj,Kmm) ) * tmask(ji,jj,1) 1004 1004 END DO 1005 1005 END DO … … 1008 1008 DO jj=j1,j2 1009 1009 DO ji=i1,i2 1010 ssh n(ji,jj) = tabres(ji,jj) * tmask(ji,jj,1)1010 ssh(ji,jj,Kmm) = tabres(ji,jj) * tmask(ji,jj,1) 1011 1011 END DO 1012 1012 END DO 1013 1013 ! 1014 1014 IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 1015 ssh b(i1:i2,j1:j2) = sshn(i1:i2,j1:j2)1015 ssh(i1:i2,j1:j2,Kbb) = ssh(i1:i2,j1:j2,Kmm) 1016 1016 ENDIF 1017 1017 ! … … 1094 1094 DO jj=j1,j2 1095 1095 zcor = rdt * r1_e1e2t(i1 ,jj) * e2u(i1,jj) * (ub2_b(i1,jj)-tabres(i1,jj)) 1096 ssh n(i1 ,jj) = sshn(i1 ,jj) + zcor1097 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) ssh b(i1 ,jj) = sshb(i1 ,jj) + atfp * zcor1096 ssh(i1 ,jj,Kmm) = ssh(i1 ,jj,Kmm) + zcor 1097 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) ssh(i1 ,jj,Kbb) = ssh(i1 ,jj,Kbb) + atfp * zcor 1098 1098 END DO 1099 1099 ENDIF … … 1101 1101 DO jj=j1,j2 1102 1102 zcor = - rdt * r1_e1e2t(i2+1,jj) * e2u(i2,jj) * (ub2_b(i2,jj)-tabres(i2,jj)) 1103 ssh n(i2+1,jj) = sshn(i2+1,jj) + zcor1104 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) ssh b(i2+1,jj) = sshb(i2+1,jj) + atfp * zcor1103 ssh(i2+1,jj,Kmm) = ssh(i2+1,jj,Kmm) + zcor 1104 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) ssh(i2+1,jj,Kbb) = ssh(i2+1,jj,Kbb) + atfp * zcor 1105 1105 END DO 1106 1106 ENDIF … … 1182 1182 DO ji=i1,i2 1183 1183 zcor = rdt * r1_e1e2t(ji,j1 ) * e1v(ji,j1 ) * (vb2_b(ji,j1)-tabres(ji,j1)) 1184 ssh n(ji,j1 ) = sshn(ji,j1) + zcor1185 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) ssh b(ji,j1 ) = sshb(ji,j1) + atfp * zcor1184 ssh(ji,j1 ,Kmm) = ssh(ji,j1 ,Kmm) + zcor 1185 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) ssh(ji,j1 ,Kbb) = ssh(ji,j1,Kbb) + atfp * zcor 1186 1186 END DO 1187 1187 ENDIF … … 1189 1189 DO ji=i1,i2 1190 1190 zcor = - rdt * r1_e1e2t(ji,j2+1) * e1v(ji,j2 ) * (vb2_b(ji,j2)-tabres(ji,j2)) 1191 ssh n(ji,j2+1) = sshn(ji,j2+1) + zcor1192 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) ssh b(ji,j2+1) = sshb(ji,j2+1) + atfp * zcor1191 ssh(ji,j2+1,Kmm) = ssh(ji,j2+1,Kmm) + zcor 1192 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) ssh(ji,j2+1,Kbb) = ssh(ji,j2+1,Kbb) + atfp * zcor 1193 1193 END DO 1194 1194 ENDIF … … 1319 1319 DO jj=j1,j2 1320 1320 DO ji=i1,i2 1321 ptab(ji,jj,jk) = e3t_0(ji,jj,jk) * (1._wp + ssh n(ji,jj) &1321 ptab(ji,jj,jk) = e3t_0(ji,jj,jk) * (1._wp + ssh(ji,jj,Kmm) & 1322 1322 & *ssmask(ji,jj)/(ht_0(ji,jj)-1._wp + ssmask(ji,jj))) 1323 1323 END DO … … 1330 1330 ! Save "old" scale factor (prior update) for subsequent asselin correction 1331 1331 ! 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)1332 e3t(i1:i2,j1:j2,1:jpkm1,Krhs) = e3t(i1:i2,j1:j2,1:jpkm1,Kmm) 1333 1334 ! One should also save e3t(:,:,:,Kbb), but lacking of workspace... 1335 ! hdiv(i1:i2,j1:j2,1:jpkm1) = e3t(i1:i2,j1:j2,1:jpkm1,Kbb) 1336 1336 1337 1337 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0) )) THEN … … 1339 1339 DO jj=j1,j2 1340 1340 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) )1341 e3t(ji,jj,jk,Kbb) = e3t(ji,jj,jk,Kbb) & 1342 & + atfp * ( ptab(ji,jj,jk) - e3t(ji,jj,jk,Kmm) ) 1343 1343 END DO 1344 1344 END DO 1345 1345 END DO 1346 1346 ! 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)1347 e3w (i1:i2,j1:j2,1,Kbb) = e3w_0(i1:i2,j1:j2,1) + e3t(i1:i2,j1:j2,1,Kbb) - e3t_0(i1:i2,j1:j2,1) 1348 gdepw(i1:i2,j1:j2,1,Kbb) = 0.0_wp 1349 gdept(i1:i2,j1:j2,1,Kbb) = 0.5_wp * e3w(i1:i2,j1:j2,1,Kbb) 1350 1350 ! 1351 1351 DO jk = 2, jpk … … 1353 1353 DO ji = i1,i2 1354 1354 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) ) &1355 e3w(ji,jj,jk,Kbb) = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) * & 1356 & ( e3t(ji,jj,jk-1,Kbb) - e3t_0(ji,jj,jk-1) ) & 1357 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))1358 & ( e3t(ji,jj,jk ,Kbb) - e3t_0(ji,jj,jk ) ) 1359 gdepw(ji,jj,jk,Kbb) = gdepw(ji,jj,jk-1,Kbb) + e3t(ji,jj,jk-1,Kbb) 1360 gdept(ji,jj,jk,Kbb) = zcoef * ( gdepw(ji,jj,jk ,Kbb) + 0.5 * e3w(ji,jj,jk,Kbb)) & 1361 & + (1-zcoef) * ( gdept(ji,jj,jk-1,Kbb) + e3w(ji,jj,jk,Kbb)) 1362 1362 END DO 1363 1363 END DO … … 1370 1370 ! 1371 1371 ! Update vertical scale factor at T-points: 1372 e3t _n(i1:i2,j1:j2,1:jpkm1) = ptab(i1:i2,j1:j2,1:jpkm1)1372 e3t(i1:i2,j1:j2,1:jpkm1,Kmm) = ptab(i1:i2,j1:j2,1:jpkm1) 1373 1373 ! 1374 1374 ! Update total depth: 1375 1375 ht_n(i1:i2,j1:j2) = 0._wp 1376 1376 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)1377 ht_n(i1:i2,j1:j2) = ht_n(i1:i2,j1:j2) + e3t(i1:i2,j1:j2,jk,Kmm) * tmask(i1:i2,j1:j2,jk) 1378 1378 END DO 1379 1379 ! 1380 1380 ! 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 ssh1381 e3w (i1:i2,j1:j2,1,Kmm) = e3w_0(i1:i2,j1:j2,1) + e3t(i1:i2,j1:j2,1,Kmm) - e3t_0(i1:i2,j1:j2,1) 1382 gdept(i1:i2,j1:j2,1,Kmm) = 0.5_wp * e3w(i1:i2,j1:j2,1,Kmm) 1383 gdepw(i1:i2,j1:j2,1,Kmm) = 0.0_wp 1384 gde3w(i1:i2,j1:j2,1) = gdept(i1:i2,j1:j2,1,Kmm) - (ht_n(i1:i2,j1:j2)-ht_0(i1:i2,j1:j2)) ! Last term in the rhs is ssh 1385 1385 ! 1386 1386 DO jk = 2, jpk … … 1388 1388 DO ji = i1,i2 1389 1389 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 ssh1390 e3w(ji,jj,jk,Kmm) = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) * ( e3t(ji,jj,jk-1,Kmm) - e3t_0(ji,jj,jk-1) ) & 1391 & + 0.5_wp * tmask(ji,jj,jk) * ( e3t(ji,jj,jk ,Kmm) - e3t_0(ji,jj,jk ) ) 1392 gdepw(ji,jj,jk,Kmm) = gdepw(ji,jj,jk-1,Kmm) + e3t(ji,jj,jk-1,Kmm) 1393 gdept(ji,jj,jk,Kmm) = zcoef * ( gdepw(ji,jj,jk ,Kmm) + 0.5 * e3w(ji,jj,jk,Kmm)) & 1394 & + (1-zcoef) * ( gdept(ji,jj,jk-1,Kmm) + e3w(ji,jj,jk,Kmm)) 1395 gde3w(ji,jj,jk) = gdept(ji,jj,jk,Kmm) - (ht_n(ji,jj)-ht_0(ji,jj)) ! Last term in the rhs is ssh 1396 1396 END DO 1397 1397 END DO … … 1399 1399 ! 1400 1400 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)1401 e3t (i1:i2,j1:j2,1:jpk,Kbb) = e3t (i1:i2,j1:j2,1:jpk,Kmm) 1402 e3w (i1:i2,j1:j2,1:jpk,Kbb) = e3w (i1:i2,j1:j2,1:jpk,Kmm) 1403 gdepw(i1:i2,j1:j2,1:jpk,Kbb) = gdepw(i1:i2,j1:j2,1:jpk,Kmm) 1404 gdept(i1:i2,j1:j2,1:jpk,Kbb) = gdept(i1:i2,j1:j2,1:jpk,Kmm) 1405 1405 ENDIF 1406 1406 !
Note: See TracChangeset
for help on using the changeset viewer.
