- Timestamp:
- 2016-12-01T11:30:29+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90
r6152 r7412 33 33 USE dynvor ! vorticity term 34 34 USE wet_dry ! wetting/drying flux limter 35 USE bdy_ par ! for lk_bdy35 USE bdy_oce , ONLY: ln_bdy 36 36 USE bdytides ! open boundary condition data 37 37 USE bdydyn2d ! open boundary conditions on barotropic variables … … 156 156 REAL(wp), POINTER, DIMENSION(:,:) :: zhf 157 157 REAL(wp), POINTER, DIMENSION(:,:) :: zcpx, zcpy ! Wetting/Dying gravity filter coef. 158 REAL(wp), POINTER, DIMENSION(:,:) :: wduflt1, wdvflt1 ! Wetting/Dying velocity filter coef.159 158 !!---------------------------------------------------------------------- 160 159 ! … … 168 167 CALL wrk_alloc( jpi,jpj, zsshu_a, zsshv_a ) 169 168 CALL wrk_alloc( jpi,jpj, zhf ) 170 IF( ln_wd ) CALL wrk_alloc( jpi, jpj, zcpx, zcpy , wduflt1, wdvflt1)169 IF( ln_wd ) CALL wrk_alloc( jpi, jpj, zcpx, zcpy ) 171 170 ! 172 171 zmdi=1.e+20 ! missing data indicator for masking … … 374 373 IF( .NOT.ln_linssh ) THEN ! Variable volume : remove surface pressure gradient 375 374 IF( ln_wd ) THEN ! Calculating and applying W/D gravity filters 376 wduflt1(:,:) = 1.0_wp377 wdvflt1(:,:) = 1.0_wp378 DO jj = 2, jpjm1379 DO ji = 2, jpim1380 ll_tmp1 = MIN(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-bathy(ji,jj), -bathy(ji+1,jj))&381 & .and. MAX(sshn(ji,jj) + bathy(ji,jj), sshn(ji+1,jj) + bathy(ji+1,jj)) &382 & > rn_wdmin1 + rn_wdmin2383 ll_tmp2 = MAX(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-bathy(ji,jj), -bathy(ji+1,jj)) &384 & + rn_wdmin1 + rn_wdmin2375 DO jj = 2, jpjm1 376 DO ji = 2, jpim1 377 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji+1,jj) ) > & 378 & MAX( -bathy(ji,jj) , -bathy(ji+1,jj) ) .AND. & 379 & MAX( sshn(ji,jj) + bathy(ji,jj), sshn(ji+1,jj) + bathy(ji+1,jj) ) & 380 & > rn_wdmin1 + rn_wdmin2 381 ll_tmp2 = MAX( sshn(ji,jj) , sshn(ji+1,jj) ) > & 382 & MAX( -bathy(ji,jj) , -bathy(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 383 385 384 IF(ll_tmp1) THEN 386 zcpx(ji,jj) 387 ELSE IF(ll_tmp2) THEN388 ! no worries about sshn(ji+1,jj)-sshn(ji,jj) = 0, it won't happenhere389 zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + bathy(ji+1,jj) - sshn(ji,jj) - bathy(ji,jj)) &390 & /(sshn(ji+1,jj) - sshn(ji,jj)))385 zcpx(ji,jj) = 1.0_wp 386 ELSE IF(ll_tmp2) THEN 387 ! no worries about sshn(ji+1,jj) - sshn(ji ,jj) = 0, it won't happen ! here 388 zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + bathy(ji+1,jj) - sshn(ji,jj) - bathy(ji,jj)) & 389 & / (sshn(ji+1,jj) - sshn(ji ,jj)) ) 391 390 ELSE 392 zcpx(ji,jj) = 0._wp 393 wduflt1(ji,jj) = 0.0_wp 391 zcpx(ji,jj) = 0._wp 394 392 END IF 395 396 ll_tmp1 = MIN(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-bathy(ji,jj), -bathy(ji,jj+1)) & 397 & .and. MAX(sshn(ji,jj) + bathy(ji,jj), sshn(ji,jj+1) + bathy(ji,jj+1)) & 398 & > rn_wdmin1 + rn_wdmin2 399 ll_tmp2 = MAX(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-bathy(ji,jj), -bathy(ji,jj+1)) & 400 & + rn_wdmin1 + rn_wdmin2 393 394 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji,jj+1) ) > & 395 & MAX( -bathy(ji,jj) , -bathy(ji,jj+1) ) .AND. & 396 & MAX( sshn(ji,jj) + bathy(ji,jj), sshn(ji,jj+1) + bathy(ji,jj+1) ) & 397 & > rn_wdmin1 + rn_wdmin2 398 ll_tmp2 = MAX( sshn(ji,jj) , sshn(ji,jj+1) ) > & 399 & MAX( -bathy(ji,jj) , -bathy(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 400 401 401 IF(ll_tmp1) THEN 402 zcpy(ji,jj)= 1.0_wp403 ELSE IF(ll_tmp2) THEN404 ! no worries about sshn(ji,jj+1)-sshn(ji,jj) = 0, it won't happenhere405 zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + bathy(ji,jj+1) - sshn(ji,jj) - bathy(ji,jj)) &406 & /(sshn(ji,jj+1) - sshn(ji,jj)))402 zcpy(ji,jj) = 1.0_wp 403 ELSE IF(ll_tmp2) THEN 404 ! no worries about sshn(ji,jj+1) - sshn(ji,jj ) = 0, it won't happen ! here 405 zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + bathy(ji,jj+1) - sshn(ji,jj) - bathy(ji,jj)) & 406 & / (sshn(ji,jj+1) - sshn(ji,jj )) ) 407 407 ELSE 408 zcpy(ji,jj) = 0._wp 409 wdvflt1(ji,jj) = 0.0_wp 410 ENDIF 411 412 END DO 408 zcpy(ji,jj) = 0._wp 409 END IF 410 END DO 413 411 END DO 414 415 CALL lbc_lnk( zcpx, 'U', 1._wp ) ; CALL lbc_lnk( zcpy, 'V', 1._wp ) 416 412 417 413 DO jj = 2, jpjm1 418 414 DO ji = 2, jpim1 419 zu_trd(ji,jj) = (zu_trd(ji,jj) - grav * ( sshn(ji+1,jj ) - sshn(ji ,jj ) ) &420 & * r1_e1u(ji,jj) ) * zcpx(ji,jj) * wduflt1(ji,jj)421 zv_trd(ji,jj) = (zv_trd(ji,jj) - grav * ( sshn(ji ,jj+1) - sshn(ji ,jj ) ) &422 & * r1_e2v(ji,jj) ) * zcpy(ji,jj) * wdvflt1(ji,jj)415 zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( sshn(ji+1,jj ) - sshn(ji ,jj ) ) & 416 & * r1_e1u(ji,jj) * zcpx(ji,jj) 417 zv_trd(ji,jj) = zv_trd(ji,jj) - grav * ( sshn(ji ,jj+1) - sshn(ji ,jj ) ) & 418 & * r1_e2v(ji,jj) * zcpy(ji,jj) 423 419 END DO 424 420 END DO … … 567 563 ENDIF 568 564 569 IF( ln_wd ) THEN !preserve the positivity of water depth570 !ssh[b,n,a] should have already been processed for this571 sshbb_e(:,:) = MAX(sshbb_e(:,:), rn_wdmin1 - bathy(:,:))572 sshb_e(:,:) = MAX(sshb_e(:,:) , rn_wdmin1 - bathy(:,:))573 ENDIF574 565 ! 575 566 IF (ln_bt_fw) THEN ! FORWARD integration: start from NOW fields … … 607 598 ! ! ------------------ 608 599 ! Update only tidal forcing at open boundaries 609 #if defined key_tide 610 IF( lk_bdy .AND. lk_tide ) CALL bdy_dta_tides( kt, kit=jn, time_offset= noffset+1 ) 611 IF( ln_tide_pot .AND. lk_tide ) CALL upd_tide ( kt, kit=jn, time_offset= noffset ) 612 #endif 600 IF( ln_bdy .AND. ln_tide ) CALL bdy_dta_tides( kt, kit=jn, time_offset= noffset+1 ) 601 IF( ln_tide_pot .AND. ln_tide ) CALL upd_tide ( kt, kit=jn, time_offset= noffset ) 613 602 ! 614 603 ! Set extrapolation coefficients for predictor step: … … 646 635 zhup2_e (:,:) = hu_0(:,:) + zwx(:,:) ! Ocean depth at U- and V-points 647 636 zhvp2_e (:,:) = hv_0(:,:) + zwy(:,:) 648 IF( ln_wd ) THEN649 zhup2_e(:,:) = MAX(zhup2_e (:,:), rn_wdmin1)650 zhvp2_e(:,:) = MAX(zhvp2_e (:,:), rn_wdmin1)651 END IF652 637 ELSE 653 638 zhup2_e (:,:) = hu_n(:,:) … … 701 686 END DO 702 687 END DO 688 703 689 ssha_e(:,:) = ( sshn_e(:,:) - rdtbt * ( zssh_frc(:,:) + zhdiv(:,:) ) ) * ssmask(:,:) 704 IF( ln_wd ) ssha_e(:,:) = MAX(ssha_e(:,:), rn_wdmin1 - bathy(:,:))690 705 691 CALL lbc_lnk( ssha_e, 'T', 1._wp ) 706 692 707 #if defined key_bdy708 693 ! Duplicate sea level across open boundaries (this is only cosmetic if linssh=T) 709 IF( l k_bdy ) CALL bdy_ssh( ssha_e )710 #endif 694 IF( ln_bdy ) CALL bdy_ssh( ssha_e ) 695 711 696 #if defined key_agrif 712 697 IF( .NOT.Agrif_Root() ) CALL agrif_ssh_ts( jn ) … … 749 734 zsshp2_e(:,:) = za0 * ssha_e(:,:) + za1 * sshn_e (:,:) & 750 735 & + za2 * sshb_e(:,:) + za3 * sshbb_e(:,:) 736 751 737 IF( ln_wd ) THEN ! Calculating and applying W/D gravity filters 752 wduflt1(:,:) = 1._wp753 wdvflt1(:,:) = 1._wp754 738 DO jj = 2, jpjm1 755 DO ji = 2, jpim1 756 ll_tmp1 = MIN( zsshp2_e(ji,jj), zsshp2_e(ji+1,jj) ) > MAX( -bathy(ji,jj), -bathy(ji+1,jj) ) & 757 & .AND. MAX( zsshp2_e(ji,jj) + bathy(ji,jj), zsshp2_e(ji+1,jj) + bathy(ji+1,jj) ) & 758 & > rn_wdmin1 + rn_wdmin2 759 ll_tmp2 = MAX( zsshp2_e(ji,jj), zsshp2_e(ji+1,jj) ) > MAX( -bathy(ji,jj), -bathy(ji+1,jj) ) & 760 & + rn_wdmin1 + rn_wdmin2 761 IF(ll_tmp1) THEN 762 zcpx(ji,jj) = 1._wp 763 ELSE IF(ll_tmp2) THEN 764 ! no worries about zsshp2_e(ji+1,jj)-zsshp2_e(ji,jj) = 0, it won't happen here 765 zcpx(ji,jj) = ABS( (zsshp2_e(ji+1,jj) + bathy(ji+1,jj) - zsshp2_e(ji,jj) - bathy(ji,jj)) & 766 & / (zsshp2_e(ji+1,jj) - zsshp2_e(ji,jj)) ) 767 ELSE 768 zcpx(ji,jj) = 0._wp 769 wduflt1(ji,jj) = 0._wp 770 END IF 771 772 ll_tmp1 = MIN( zsshp2_e(ji,jj), zsshp2_e(ji,jj+1) ) > MAX( -bathy(ji,jj), -bathy(ji,jj+1) ) & 773 & .AND. MAX( zsshp2_e(ji,jj) + bathy(ji,jj), zsshp2_e(ji,jj+1) + bathy(ji,jj+1) ) & 774 & > rn_wdmin1 + rn_wdmin2 775 ll_tmp2 = MAX( zsshp2_e(ji,jj), zsshp2_e(ji,jj+1) ) > MAX( -bathy(ji,jj), -bathy(ji,jj+1) ) & 776 & + rn_wdmin1 + rn_wdmin2 777 IF(ll_tmp1) THEN 778 zcpy(ji,jj) = 1._wp 779 ELSE IF(ll_tmp2) THEN 780 ! no worries about zsshp2_e(ji,jj+1)-zsshp2_e(ji,jj) = 0, it won't happen here 781 zcpy(ji,jj) = ABS( (zsshp2_e(ji,jj+1) + bathy(ji,jj+1) - zsshp2_e(ji,jj) - bathy(ji,jj)) & 782 & / (zsshp2_e(ji,jj+1) - zsshp2_e(ji,jj)) ) 783 ELSE 784 zcpy(ji,jj) = 0._wp 785 wdvflt1(ji,jj) = 0._wp 786 END IF 739 DO ji = 2, jpim1 740 ll_tmp1 = MIN( zsshp2_e(ji,jj) , zsshp2_e(ji+1,jj) ) > & 741 & MAX( -bathy(ji,jj) , -bathy(ji+1,jj) ) .AND. & 742 & MAX( zsshp2_e(ji,jj) + bathy(ji,jj), zsshp2_e(ji+1,jj) + bathy(ji+1,jj) ) & 743 & > rn_wdmin1 + rn_wdmin2 744 ll_tmp2 = MAX( zsshp2_e(ji,jj) , zsshp2_e(ji+1,jj) ) > & 745 & MAX( -bathy(ji,jj) , -bathy(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 746 747 IF(ll_tmp1) THEN 748 zcpx(ji,jj) = 1.0_wp 749 ELSE IF(ll_tmp2) THEN 750 ! no worries about zsshp2_e(ji+1,jj) - zsshp2_e(ji ,jj) = 0, it won't happen ! here 751 zcpx(ji,jj) = ABS( (zsshp2_e(ji+1,jj) + bathy(ji+1,jj) - zsshp2_e(ji,jj) - bathy(ji,jj)) & 752 & / (zsshp2_e(ji+1,jj) - zsshp2_e(ji ,jj)) ) 753 ELSE 754 zcpx(ji,jj) = 0._wp 755 END IF 756 757 ll_tmp1 = MIN( zsshp2_e(ji,jj) , zsshp2_e(ji,jj+1) ) > & 758 & MAX( -bathy(ji,jj) , -bathy(ji,jj+1) ) .AND. & 759 & MAX( zsshp2_e(ji,jj) + bathy(ji,jj), zsshp2_e(ji,jj+1) + bathy(ji,jj+1) ) & 760 & > rn_wdmin1 + rn_wdmin2 761 ll_tmp2 = MAX( zsshp2_e(ji,jj) , zsshp2_e(ji,jj+1) ) > & 762 & MAX( -bathy(ji,jj) , -bathy(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 763 764 IF(ll_tmp1) THEN 765 zcpy(ji,jj) = 1.0_wp 766 ELSE IF(ll_tmp2) THEN 767 ! no worries about zsshp2_e(ji,jj+1) - zsshp2_e(ji,jj ) = 0, it won't happen ! here 768 zcpy(ji,jj) = ABS( (zsshp2_e(ji,jj+1) + bathy(ji,jj+1) - zsshp2_e(ji,jj) - bathy(ji,jj)) & 769 & / (zsshp2_e(ji,jj+1) - zsshp2_e(ji,jj )) ) 770 ELSE 771 zcpy(ji,jj) = 0._wp 772 END IF 787 773 END DO 788 END DO 789 CALL lbc_lnk( zcpx, 'U', 1._wp ) ; CALL lbc_lnk( zcpy, 'V', 1._wp ) 790 ENDIF 774 END DO 775 END IF 791 776 ! 792 777 ! Compute associated depths at U and V points: … … 806 791 END DO 807 792 808 IF( ln_wd ) THEN809 zhust_e(:,:) = MAX(zhust_e (:,:), rn_wdmin1 )810 zhvst_e(:,:) = MAX(zhvst_e (:,:), rn_wdmin1 )811 END IF812 813 793 ENDIF 814 794 ! … … 861 841 ! 862 842 ! Add tidal astronomical forcing if defined 863 IF ( l k_tide.AND.ln_tide_pot ) THEN843 IF ( ln_tide.AND.ln_tide_pot ) THEN 864 844 DO jj = 2, jpjm1 865 845 DO ji = fs_2, fs_jpim1 ! vector opt. … … 888 868 zu_spg = - grav * ( zsshp2_e(ji+1,jj) - zsshp2_e(ji,jj) ) * r1_e1u(ji,jj) 889 869 zv_spg = - grav * ( zsshp2_e(ji,jj+1) - zsshp2_e(ji,jj) ) * r1_e2v(ji,jj) 890 zwx(ji,jj) = zu_spg * zcpx(ji,jj) 891 zwy(ji,jj) = zv_spg * zcpy(ji,jj) 870 zwx(ji,jj) = zu_spg * zcpx(ji,jj) * wdmask(ji,jj) * wdmask(ji+1, jj) 871 zwy(ji,jj) = zv_spg * zcpy(ji,jj) * wdmask(ji,jj) * wdmask(ji, jj+1) 892 872 END DO 893 873 END DO … … 927 907 DO ji = fs_2, fs_jpim1 ! vector opt. 928 908 929 IF( ln_wd ) THEN 930 zhura = MAX(hu_0(ji,jj) + zsshu_a(ji,jj), rn_wdmin1) 931 zhvra = MAX(hv_0(ji,jj) + zsshv_a(ji,jj), rn_wdmin1) 932 ELSE 933 zhura = hu_0(ji,jj) + zsshu_a(ji,jj) 934 zhvra = hv_0(ji,jj) + zsshv_a(ji,jj) 935 END IF 909 zhura = hu_0(ji,jj) + zsshu_a(ji,jj) 910 zhvra = hv_0(ji,jj) + zsshv_a(ji,jj) 936 911 zhura = ssumask(ji,jj)/(zhura + 1._wp - ssumask(ji,jj)) 937 912 zhvra = ssvmask(ji,jj)/(zhvra + 1._wp - ssvmask(ji,jj)) … … 953 928 ! 954 929 IF( .NOT.ln_linssh ) THEN !* Update ocean depth (variable volume case only) 955 IF( ln_wd ) THEN 956 hu_e (:,:) = MAX(hu_0(:,:) + zsshu_a(:,:), rn_wdmin1) 957 hv_e (:,:) = MAX(hv_0(:,:) + zsshv_a(:,:), rn_wdmin1) 958 ELSE 959 hu_e (:,:) = hu_0(:,:) + zsshu_a(:,:) 960 hv_e (:,:) = hv_0(:,:) + zsshv_a(:,:) 961 END IF 930 hu_e (:,:) = hu_0(:,:) + zsshu_a(:,:) 931 hv_e (:,:) = hv_0(:,:) + zsshv_a(:,:) 962 932 hur_e(:,:) = ssumask(:,:) / ( hu_e(:,:) + 1._wp - ssumask(:,:) ) 963 933 hvr_e(:,:) = ssvmask(:,:) / ( hv_e(:,:) + 1._wp - ssvmask(:,:) ) … … 967 937 CALL lbc_lnk_multi( ua_e, 'U', -1._wp, va_e , 'V', -1._wp ) 968 938 ! 969 #if defined key_bdy970 939 ! ! open boundaries 971 IF( l k_bdy ) CALL bdy_dyn2d( jn, ua_e, va_e, un_e, vn_e, hur_e, hvr_e, ssha_e )972 #endif 940 IF( ln_bdy ) CALL bdy_dyn2d( jn, ua_e, va_e, un_e, vn_e, hur_e, hvr_e, ssha_e ) 941 973 942 #if defined key_agrif 974 943 IF( .NOT.Agrif_Root() ) CALL agrif_dyn_ts( jn ) ! Agrif … … 1024 993 ! 1025 994 ! Update barotropic trend: 1026 IF( ln_dynadv_vec .OR. ln_linssh ) THEN 1027 DO jk=1,jpkm1 1028 ua(:,:,jk) = ua(:,:,jk) + ( ua_b(:,:) - ub_b(:,:) ) * z1_2dt_b 1029 va(:,:,jk) = va(:,:,jk) + ( va_b(:,:) - vb_b(:,:) ) * z1_2dt_b 1030 END DO 995 IF(ln_wd) THEN 996 IF( ln_dynadv_vec .OR. ln_linssh ) THEN 997 DO jk=1,jpkm1 998 ua(:,:,jk) = ua(:,:,jk) + ( ua_b(:,:) - ub_b(:,:) ) * z1_2dt_b * wdmask(:,:) 999 va(:,:,jk) = va(:,:,jk) + ( va_b(:,:) - vb_b(:,:) ) * z1_2dt_b * wdmask(:,:) 1000 END DO 1001 ELSE 1002 ! At this stage, ssha has been corrected: compute new depths at velocity points 1003 DO jj = 1, jpjm1 1004 DO ji = 1, jpim1 ! NO Vector Opt. 1005 zsshu_a(ji,jj) = z1_2 * umask(ji,jj,1) * r1_e1e2u(ji,jj) & 1006 & * ( e1e2t(ji ,jj) * ssha(ji ,jj) & 1007 & + e1e2t(ji+1,jj) * ssha(ji+1,jj) ) 1008 zsshv_a(ji,jj) = z1_2 * vmask(ji,jj,1) * r1_e1e2v(ji,jj) & 1009 & * ( e1e2t(ji,jj ) * ssha(ji,jj ) & 1010 & + e1e2t(ji,jj+1) * ssha(ji,jj+1) ) 1011 END DO 1012 END DO 1013 CALL lbc_lnk_multi( zsshu_a, 'U', 1._wp, zsshv_a, 'V', 1._wp ) ! Boundary conditions 1014 ! 1015 DO jk=1,jpkm1 1016 ua(:,:,jk) = ua(:,:,jk) + r1_hu_n(:,:) * ( ua_b(:,:) - ub_b(:,:) * hu_b(:,:) ) * z1_2dt_b * wdmask(:,:) 1017 va(:,:,jk) = va(:,:,jk) + r1_hv_n(:,:) * ( va_b(:,:) - vb_b(:,:) * hv_b(:,:) ) * z1_2dt_b * wdmask(:,:) 1018 END DO 1019 ! Save barotropic velocities not transport: 1020 ua_b(:,:) = ua_b(:,:) / ( hu_0(:,:) + zsshu_a(:,:) + 1._wp - ssumask(:,:) ) 1021 va_b(:,:) = va_b(:,:) / ( hv_0(:,:) + zsshv_a(:,:) + 1._wp - ssvmask(:,:) ) 1022 ENDIF 1031 1023 ELSE 1032 ! At this stage, ssha has been corrected: compute new depths at velocity points 1033 DO jj = 1, jpjm1 1034 DO ji = 1, jpim1 ! NO Vector Opt. 1035 zsshu_a(ji,jj) = z1_2 * umask(ji,jj,1) * r1_e1e2u(ji,jj) & 1036 & * ( e1e2t(ji ,jj) * ssha(ji ,jj) & 1037 & + e1e2t(ji+1,jj) * ssha(ji+1,jj) ) 1038 zsshv_a(ji,jj) = z1_2 * vmask(ji,jj,1) * r1_e1e2v(ji,jj) & 1039 & * ( e1e2t(ji,jj ) * ssha(ji,jj ) & 1040 & + e1e2t(ji,jj+1) * ssha(ji,jj+1) ) 1041 END DO 1042 END DO 1043 CALL lbc_lnk_multi( zsshu_a, 'U', 1._wp, zsshv_a, 'V', 1._wp ) ! Boundary conditions 1044 ! 1045 DO jk=1,jpkm1 1046 ua(:,:,jk) = ua(:,:,jk) + r1_hu_n(:,:) * ( ua_b(:,:) - ub_b(:,:) * hu_b(:,:) ) * z1_2dt_b 1047 va(:,:,jk) = va(:,:,jk) + r1_hv_n(:,:) * ( va_b(:,:) - vb_b(:,:) * hv_b(:,:) ) * z1_2dt_b 1048 END DO 1049 ! Save barotropic velocities not transport: 1050 ua_b(:,:) = ua_b(:,:) / ( hu_0(:,:) + zsshu_a(:,:) + 1._wp - ssumask(:,:) ) 1051 va_b(:,:) = va_b(:,:) / ( hv_0(:,:) + zsshv_a(:,:) + 1._wp - ssvmask(:,:) ) 1052 ENDIF 1024 IF( ln_dynadv_vec .OR. ln_linssh ) THEN 1025 DO jk=1,jpkm1 1026 ua(:,:,jk) = ua(:,:,jk) + ( ua_b(:,:) - ub_b(:,:) ) * z1_2dt_b 1027 va(:,:,jk) = va(:,:,jk) + ( va_b(:,:) - vb_b(:,:) ) * z1_2dt_b 1028 END DO 1029 ELSE 1030 ! At this stage, ssha has been corrected: compute new depths at velocity points 1031 DO jj = 1, jpjm1 1032 DO ji = 1, jpim1 ! NO Vector Opt. 1033 zsshu_a(ji,jj) = z1_2 * umask(ji,jj,1) * r1_e1e2u(ji,jj) & 1034 & * ( e1e2t(ji ,jj) * ssha(ji ,jj) & 1035 & + e1e2t(ji+1,jj) * ssha(ji+1,jj) ) 1036 zsshv_a(ji,jj) = z1_2 * vmask(ji,jj,1) * r1_e1e2v(ji,jj) & 1037 & * ( e1e2t(ji,jj ) * ssha(ji,jj ) & 1038 & + e1e2t(ji,jj+1) * ssha(ji,jj+1) ) 1039 END DO 1040 END DO 1041 CALL lbc_lnk_multi( zsshu_a, 'U', 1._wp, zsshv_a, 'V', 1._wp ) ! Boundary conditions 1042 ! 1043 DO jk=1,jpkm1 1044 ua(:,:,jk) = ua(:,:,jk) + r1_hu_n(:,:) * ( ua_b(:,:) - ub_b(:,:) * hu_b(:,:) ) * z1_2dt_b 1045 va(:,:,jk) = va(:,:,jk) + r1_hv_n(:,:) * ( va_b(:,:) - vb_b(:,:) * hv_b(:,:) ) * z1_2dt_b 1046 END DO 1047 ! Save barotropic velocities not transport: 1048 ua_b(:,:) = ua_b(:,:) / ( hu_0(:,:) + zsshu_a(:,:) + 1._wp - ssumask(:,:) ) 1049 va_b(:,:) = va_b(:,:) / ( hv_0(:,:) + zsshv_a(:,:) + 1._wp - ssvmask(:,:) ) 1050 ENDIF 1051 1052 END IF 1053 1053 ! 1054 1054 DO jk = 1, jpkm1 … … 1086 1086 CALL wrk_dealloc( jpi,jpj, zsshu_a, zsshv_a ) 1087 1087 CALL wrk_dealloc( jpi,jpj, zhf ) 1088 IF( ln_wd ) CALL wrk_dealloc( jpi, jpj, zcpx, zcpy , wduflt1, wdvflt1)1088 IF( ln_wd ) CALL wrk_dealloc( jpi, jpj, zcpx, zcpy ) 1089 1089 ! 1090 1090 IF ( ln_diatmb ) THEN
Note: See TracChangeset
for help on using the changeset viewer.
