Changeset 15173
- Timestamp:
- 2021-08-04T16:43:12+02:00 (3 years ago)
- Location:
- NEMO/branches/NERC/NEMO_4.0.4_CO9_package_tides/src/OCE
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/NERC/NEMO_4.0.4_CO9_package_tides/src/OCE/DOM/dom_oce.F90
r14075 r15173 152 152 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: r1_hv_b , r1_hv_n , r1_hv_a !: inverse of v-depth [1/m] 153 153 154 !CEOD Scale of water column down to shallowest of neighbourinbg points over total 155 !water depth 156 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: scaled_e3t_0_ik , scaled_e3t_0_jk 157 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: scaled_e3t_0_ip1k , scaled_e3t_0_jp1k 158 !CEOD 154 159 INTEGER, PUBLIC :: nla10 !: deepest W level Above ~10m (nlb10 - 1) 155 160 INTEGER, PUBLIC :: nlb10 !: shallowest W level Bellow ~10m (nla10 + 1) … … 294 299 ! 295 300 ALLOCATE( wmask(jpi,jpj,jpk) , wumask(jpi,jpj,jpk), wvmask(jpi,jpj,jpk) , STAT=ierr(12) ) 301 !CEOD 302 ALLOCATE( scaled_e3t_0_ik (jpi,jpj) , scaled_e3t_0_jk (jpi,jpj) , & 303 & scaled_e3t_0_ip1k(jpi,jpj) , scaled_e3t_0_jp1k(jpi,jpj) , STAT=ierr(13) ) 304 ! 296 305 ! 297 306 dom_oce_alloc = MAXVAL(ierr) -
NEMO/branches/NERC/NEMO_4.0.4_CO9_package_tides/src/OCE/DOM/domain.F90
r14078 r15173 79 79 INTEGER , DIMENSION(jpi,jpj) :: ik_top , ik_bot ! top and bottom ocean level 80 80 REAL(wp), DIMENSION(jpi,jpj) :: z1_hu_0, z1_hv_0 81 !CEOD 82 REAL(wp), DIMENSION(jpi,jpj) :: k_bot_i_min, k_bot_j_min 83 REAL(wp), DIMENSION(jpi,jpj) :: sum_e3t_0_min_ik 84 REAL(wp), DIMENSION(jpi,jpj) :: sum_e3t_0_min_jk 85 REAL(wp), DIMENSION(jpi,jpj) :: sum_e3t_0_min_ip1k 86 REAL(wp), DIMENSION(jpi,jpj) :: sum_e3t_0_min_jp1k 81 87 !!---------------------------------------------------------------------- 82 88 ! … … 155 161 hv_0(:,:) = hv_0(:,:) + e3v_0(:,:,jk) * vmask(:,:,jk) 156 162 END DO 163 !CEOD Here we need to work out which is the shallowest point in terms of k 164 ! levels at neighbouring points in i and j directions. 165 ! Then we need to work out the proportion of the water column down to ht_0 that 166 ! that level is for both i points. 167 ! e.g. point i could say go down to level 10 but point i+1 only level 5 168 ! then work out depth of bottom of level 5 at point i over depth at point i to 169 ! level 10 170 ! This is needed later to work out where a u point depth should be when using 171 ! enveloping coordinates, as it will be 1/2 depth at t point at level 5 for i,i+1 172 ! points Because level 5 at i can change in time, thus u bed can change in time! 173 ! What about under ice shelves? Need rto return to that later 174 175 ! 1) get the shallowest k level at neighbouring i and j points store in ! k_bot_i_min 176 177 DO jj = 1, jpjm1 178 DO ji = 1, jpim1 ! SPG with the application of W/D gravity filters 179 k_bot_i_min(ji,jj) = MIN( ik_bot(ji,jj), ik_bot(ji+1,jj )) 180 k_bot_j_min(ji,jj) = MIN( ik_bot(ji,jj), ik_bot(ji, jj+1)) 181 ENDDO 182 ENDDO 183 ! 2) Work out the depth down to the shallowest k level both at point i and i+1 184 ! (likewise for j) 185 186 sum_e3t_0_min_ik(:,:) = 0 187 sum_e3t_0_min_jk(:,:) = 0 188 sum_e3t_0_min_ip1k(:,:) = 0 189 sum_e3t_0_min_jp1k(:,:) = 0 190 191 !Get sum of e3t_0s down to local min 192 DO jj = 1, jpjm1 193 DO ji = 1, jpim1 194 DO jk = 1, k_bot_i_min(ji,jj) 195 196 sum_e3t_0_min_ik (ji,jj) = sum_e3t_0_min_ik (ji,jj) + e3t_0(ji ,jj,jk)*tmask(ji,jj,jk) 197 sum_e3t_0_min_ip1k(ji,jj) = sum_e3t_0_min_ip1k(ji,jj) + e3t_0(ji+1,jj,jk)*tmask(ji+1,jj,jk) 198 ENDDO 199 200 DO jk = 1, k_bot_j_min(ji,jj) 201 sum_e3t_0_min_jk (ji,jj) = sum_e3t_0_min_jk (ji,jj) + e3t_0(ji,jj ,jk)*tmask(ji,jj,jk) 202 sum_e3t_0_min_jp1k(ji,jj) = sum_e3t_0_min_jp1k(ji,jj) + e3t_0(ji,jj+1,jk)*tmask(ji,jj+1,jk) 203 ENDDO 204 ! 3) Then work out what fraction of the at rest water column that is, we later 205 ! multiply the now water depth by this scale to work out the bottom of the kth 206 ! level at time now in dynspg_ts 207 scaled_e3t_0_ik (ji,jj) = sum_e3t_0_min_ik (ji,jj)/( ht_0(ji ,jj ) + 1._wp - ssmask(ji,jj)) 208 scaled_e3t_0_jk (ji,jj) = sum_e3t_0_min_jk (ji,jj)/( ht_0(ji ,jj ) + 1._wp - ssmask(ji,jj)) 209 scaled_e3t_0_ip1k(ji,jj) = sum_e3t_0_min_ip1k(ji,jj)/( ht_0(ji+1,jj ) + 1._wp - ssmask(ji+1,jj)) 210 scaled_e3t_0_jp1k(ji,jj) = sum_e3t_0_min_jp1k(ji,jj)/( ht_0(ji ,jj+1) + 1._wp - ssmask(ji,jj+1)) 211 ENDDO 212 ENDDO 157 213 ! 158 214 ! !== time varying part of coordinate system ==! -
NEMO/branches/NERC/NEMO_4.0.4_CO9_package_tides/src/OCE/DOM/dtatsd.F90
r14075 r15173 179 179 ptsd(:,:,:,jp_sal) = sf_tsd(jp_sal)%fnow(:,:,:) 180 180 ! 181 IF( ln_sco ) THEN !== s- or mixed s-zps-coordinate ==! 182 ! 183 IF( kt == nit000 .AND. lwp )THEN 184 WRITE(numout,*) 185 WRITE(numout,*) 'dta_tsd: interpolates T & S data onto the s- or mixed s-z-coordinate mesh' 186 ENDIF 187 ! 188 DO jj = 1, jpj ! vertical interpolation of T & S 189 DO ji = 1, jpi 190 DO jk = 1, jpk ! determines the intepolated T-S profiles at each (i,j) points 191 zl = gdept_0(ji,jj,jk) 192 IF( zl < gdept_1d(1 ) ) THEN ! above the first level of data 193 ztp(jk) = ptsd(ji,jj,1 ,jp_tem) 194 zsp(jk) = ptsd(ji,jj,1 ,jp_sal) 195 ELSEIF( zl > gdept_1d(jpk) ) THEN ! below the last level of data 196 ztp(jk) = ptsd(ji,jj,jpkm1,jp_tem) 197 zsp(jk) = ptsd(ji,jj,jpkm1,jp_sal) 198 ELSE ! inbetween : vertical interpolation between jkk & jkk+1 199 DO jkk = 1, jpkm1 ! when gdept(jkk) < zl < gdept(jkk+1) 200 IF( (zl-gdept_1d(jkk)) * (zl-gdept_1d(jkk+1)) <= 0._wp ) THEN 201 zi = ( zl - gdept_1d(jkk) ) / (gdept_1d(jkk+1)-gdept_1d(jkk)) 202 ztp(jk) = ptsd(ji,jj,jkk,jp_tem) + ( ptsd(ji,jj,jkk+1,jp_tem) - ptsd(ji,jj,jkk,jp_tem) ) * zi 203 zsp(jk) = ptsd(ji,jj,jkk,jp_sal) + ( ptsd(ji,jj,jkk+1,jp_sal) - ptsd(ji,jj,jkk,jp_sal) ) * zi 204 ENDIF 205 END DO 206 ENDIF 207 END DO 208 DO jk = 1, jpkm1 209 ptsd(ji,jj,jk,jp_tem) = ztp(jk) * tmask(ji,jj,jk) ! mask required for mixed zps-s-coord 210 ptsd(ji,jj,jk,jp_sal) = zsp(jk) * tmask(ji,jj,jk) 211 END DO 212 ptsd(ji,jj,jpk,jp_tem) = 0._wp 213 ptsd(ji,jj,jpk,jp_sal) = 0._wp 214 END DO 215 END DO 216 ! 217 ELSE !== z- or zps- coordinate ==! 181 !CEOD, We think this is incorrect we dont want to do this. 218 182 ! 219 183 ptsd(:,:,:,jp_tem) = ptsd(:,:,:,jp_tem) * tmask(:,:,:) ! Mask … … 239 203 ENDIF 240 204 ! 241 ENDIF 205 !CEOD 242 206 ! 243 207 IF( .NOT.ln_tsd_dmp ) THEN !== deallocate T & S structure ==! -
NEMO/branches/NERC/NEMO_4.0.4_CO9_package_tides/src/OCE/DYN/dynspg_ts.F90
r14075 r15173 150 150 REAL(wp) :: r1_2dt_b, z1_hu, z1_hv ! local scalars 151 151 REAL(wp) :: za0, za1, za2, za3 ! - - 152 REAL(wp), DIMENSION(jpi,jpj) :: zdep_u, zdep_v 152 153 REAL(wp) :: zztmp, zldg ! - - 153 154 REAL(wp) :: zhu_bck, zhv_bck, zhdiv ! - - … … 459 460 DO jj = 1, jpj 460 461 DO ji = 1, jpim1 ! not jpi-column 461 zhup2_e(ji,jj) = hu_0(ji,jj) + r1_2 * r1_e1e2u(ji,jj) & 462 & * ( e1e2t(ji ,jj) * zsshp2_e(ji ,jj) & 463 & + e1e2t(ji+1,jj) * zsshp2_e(ji+1,jj) ) * ssumask(ji,jj) 462 zhup2_e(ji,jj) = r1_2 * r1_e1e2u(ji,jj) * ( e1e2t(ji ,jj) * zhtp2_e(ji ,jj)*scaled_e3t_0_ik (ji,jj) & 463 & + e1e2t(ji+1,jj) * zhtp2_e(ji+1,jj)*scaled_e3t_0_ip1k(ji,jj) ) *ssumask(ji,jj) 464 464 END DO 465 465 END DO 466 466 DO jj = 1, jpjm1 ! not jpj-row 467 467 DO ji = 1, jpi 468 zhvp2_e(ji,jj) = hv_0(ji,jj) + r1_2 * r1_e1e2v(ji,jj) & 469 & * ( e1e2t(ji,jj ) * zsshp2_e(ji,jj ) & 470 & + e1e2t(ji,jj+1) * zsshp2_e(ji,jj+1) ) * ssvmask(ji,jj) 468 zhvp2_e(ji,jj) = r1_2 * r1_e1e2v(ji,jj) * ( e1e2t(ji,jj ) * zhtp2_e(ji, jj)*scaled_e3t_0_jk (ji,jj) & 469 & + e1e2t(ji,jj+1) * zhtp2_e(ji,jj+1)*scaled_e3t_0_jp1k(ji,jj) ) *ssvmask(ji,jj) 471 470 END DO 472 471 END DO … … 647 646 ! ! hu_e, hv_e hold depth at jn, zhup2_e, zhvp2_e hold extrapolated depth at jn+1/2 648 647 ! ! backward interpolated depth used in spg terms at jn+1/2 649 zhu_bck = hu_0(ji,jj) + r1_2*r1_e1e2u(ji,jj) * ( e1e2t(ji ,jj) * zsshp2_e(ji ,jj) & 650 & + e1e2t(ji+1,jj) * zsshp2_e(ji+1,jj) ) * ssumask(ji,jj) 651 zhv_bck = hv_0(ji,jj) + r1_2*r1_e1e2v(ji,jj) * ( e1e2t(ji,jj ) * zsshp2_e(ji,jj ) & 652 & + e1e2t(ji,jj+1) * zsshp2_e(ji,jj+1) ) * ssvmask(ji,jj) 653 ! ! inverse depth at jn+1 654 z1_hu = ssumask(ji,jj) / ( hu_0(ji,jj) + zsshu_a(ji,jj) + 1._wp - ssumask(ji,jj) ) 655 z1_hv = ssvmask(ji,jj) / ( hv_0(ji,jj) + zsshv_a(ji,jj) + 1._wp - ssvmask(ji,jj) ) 656 ! 648 zhu_bck = r1_2 * r1_e1e2u(ji,jj) * ( e1e2t(ji ,jj)*(ht_0(ji ,jj) + zsshp2_e(ji ,jj))*scaled_e3t_0_ik (ji,jj) & 649 & + e1e2t(ji+1,jj)*(ht_0(ji+1,jj) + zsshp2_e(ji+1,jj))*scaled_e3t_0_ip1k(ji,jj) & 650 & ) * ssumask(ji,jj) 651 652 zhv_bck = r1_2 * r1_e1e2v(ji,jj) * ( e1e2t(ji,jj )*(ht_0(ji,jj ) + zsshp2_e(ji,jj ))*scaled_e3t_0_jk (ji,jj) & 653 & + e1e2t(ji,jj+1)*(ht_0(ji,jj+1) + zsshp2_e(ji,jj+1))*scaled_e3t_0_jp1k(ji,jj) & 654 & ) * ssvmask(ji,jj) 655 656 z1_hu = ssumask(ji,jj) / ( r1_2 * r1_e1e2u(ji,jj) *(e1e2t(ji ,jj)*(ht_0(ji,jj) + ssha_e(ji,jj) )*scaled_e3t_0_ik(ji ,jj) & 657 & +e1e2t(ji+1,jj)*(ht_0(ji+1,jj) + ssha_e(ji+1,jj))*scaled_e3t_0_ip1k(ji,jj)) + 1._wp - ssumask(ji,jj) ) 658 659 z1_hv = ssvmask(ji,jj) / ( r1_2 * r1_e1e2v(ji,jj)* ( e1e2t(ji,jj )*(ht_0(ji,jj) + ssha_e(ji,jj) )*scaled_e3t_0_jk(ji ,jj) & 660 & + e1e2t(ji,jj+1)* (ht_0(ji,jj+1) + ssha_e(ji,jj+1))*scaled_e3t_0_jp1k(ji,jj)) + 1._wp - ssvmask(ji,jj) ) 661 662 657 663 ua_e(ji,jj) = ( hu_e (ji,jj) * un_e (ji,jj) & 658 664 & + rdtbt * ( zhu_bck * zu_spg (ji,jj) & ! … … 678 684 679 685 IF( .NOT.ln_linssh ) THEN !* Update ocean depth (variable volume case only) 680 hu_e (2:jpim1,2:jpjm1) = hu_0(2:jpim1,2:jpjm1) + zsshu_a(2:jpim1,2:jpjm1) 686 DO jj = 1, jpjm1 687 DO ji = 1, jpim1 ! NO Vector Opt. 688 zdep_u(ji,jj) = r1_2 * r1_e1e2u(ji,jj) * ( e1e2t(ji,jj ) * (ssha_e(ji, jj)+ht_0(ji, jj))*scaled_e3t_0_ik (ji,jj) & 689 & + e1e2t(ji+1,jj) * (ssha_e(ji+1,jj)+ht_0(ji+1,jj))*scaled_e3t_0_ip1k(ji,jj) )*ssumask(ji,jj) 690 zdep_v(ji,jj) = r1_2 * r1_e1e2v(ji,jj) * ( e1e2t(ji,jj ) * (ssha_e(ji, jj)+ht_0(ji, jj))*scaled_e3t_0_jk (ji,jj ) & 691 & + e1e2t(ji,jj+1) * (ssha_e(ji,jj+1)+ht_0(ji,jj+1))*scaled_e3t_0_jp1k(ji,jj) )*ssvmask(ji,jj) 692 ENDDO 693 ENDDO 694 CALL lbc_lnk_multi( 'dynspg_ts', zdep_u, 'U', -1._wp ) 695 CALL lbc_lnk_multi( 'dynspg_ts', zdep_v, 'V', -1._wp ) 696 697 hu_e (2:jpim1,2:jpjm1) = zdep_u(2:jpim1,2:jpjm1) 681 698 hur_e(2:jpim1,2:jpjm1) = ssumask(2:jpim1,2:jpjm1) / ( hu_e(2:jpim1,2:jpjm1) + 1._wp - ssumask(2:jpim1,2:jpjm1) ) 682 hv_e (2:jpim1,2:jpjm1) = hv_0(2:jpim1,2:jpjm1) + zsshv_a(2:jpim1,2:jpjm1)699 hv_e (2:jpim1,2:jpjm1) = zdep_v(2:jpim1,2:jpjm1) 683 700 hvr_e(2:jpim1,2:jpjm1) = ssvmask(2:jpim1,2:jpjm1) / ( hv_e(2:jpim1,2:jpjm1) + 1._wp - ssvmask(2:jpim1,2:jpjm1) ) 684 701 CALL lbc_lnk_multi( 'dynspg_ts', ua_e , 'U', -1._wp, va_e , 'V', -1._wp & … … 778 795 & * ( e1e2t(ji,jj ) * ssha(ji,jj ) & 779 796 & + e1e2t(ji,jj+1) * ssha(ji,jj+1) ) 780 END DO 781 END DO 797 zdep_u(ji,jj) = r1_2 * r1_e1e2u(ji,jj) * ( e1e2t(ji,jj ) * (ssha(ji, jj)+ht_0(ji, jj))*scaled_e3t_0_ik(ji ,jj) & 798 & + e1e2t(ji+1,jj) * (ssha(ji+1,jj)+ht_0(ji+1,jj))*scaled_e3t_0_ip1k(ji,jj) ) 799 zdep_v(ji,jj) = r1_2 * r1_e1e2v(ji,jj) * ( e1e2t(ji,jj ) * (ssha(ji, jj)+ht_0(ji, jj))*scaled_e3t_0_jk(ji,jj ) & 800 & + e1e2t(ji,jj+1) * (ssha(ji,jj+1)+ht_0(ji,jj+1))*scaled_e3t_0_jp1k(ji,jj) ) 801 END DO 802 END DO 803 CALL lbc_lnk_multi( 'dynspg_ts', zdep_u, 'U', -1._wp ) 804 CALL lbc_lnk_multi( 'dynspg_ts', zdep_v, 'V', -1._wp ) 805 782 806 CALL lbc_lnk_multi( 'dynspg_ts', zsshu_a, 'U', 1._wp, zsshv_a, 'V', 1._wp ) ! Boundary conditions 783 807 ! … … 787 811 END DO 788 812 ! Save barotropic velocities not transport: 789 ua_b(:,:) = ua_b(:,:) / ( hu_0(:,:) + zsshu_a(:,:) + 1._wp - ssumask(:,:) )790 va_b(:,:) = va_b(:,:) / ( hv_0(:,:) + zsshv_a(:,:) + 1._wp - ssvmask(:,:) )813 ua_b(:,:) = ua_b(:,:) / ( zdep_u(:,:) + 1._wp - ssumask(:,:) ) 814 va_b(:,:) = va_b(:,:) / ( zdep_v(:,:) + 1._wp - ssvmask(:,:) ) 791 815 ENDIF 792 816
Note: See TracChangeset
for help on using the changeset viewer.