Changeset 13696
- Timestamp:
- 2020-10-28T19:09:27+01:00 (3 years ago)
- Location:
- NEMO/branches/2020/dev_r13327_KERNEL-06_2_techene_e3/src/OCE
- Files:
-
- 3 edited
-
DOM/domzgr_substitute.h90 (modified) (1 diff)
-
DYN/dynspg_ts.F90 (modified) (7 diffs)
-
DYN/dynvor.F90 (modified) (11 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/dev_r13327_KERNEL-06_2_techene_e3/src/OCE/DOM/domzgr_substitute.h90
r13427 r13696 16 16 # define e3v(i,j,k,t) (e3v_0(i,j,k)*(1._wp+r3v(i,j,t)*vmask(i,j,k))) 17 17 # define e3f(i,j,k) (e3f_0(i,j,k)*(1._wp+r3f(i,j)*fmask(i,j,k))) 18 # define e3f_vor(i,j,k) (e3f_0vor(i,j,k)*(1._wp+r3f(i,j)*fmask(i,j,k))) 18 19 # define e3w(i,j,k,t) (e3w_0(i,j,k)*(1._wp+r3t(i,j,t))) 19 20 # define e3uw(i,j,k,t) (e3uw_0(i,j,k)*(1._wp+r3u(i,j,t))) -
NEMO/branches/2020/dev_r13327_KERNEL-06_2_techene_e3/src/OCE/DYN/dynspg_ts.F90
r13427 r13696 1119 1119 !! although they should be updated in the variable volume case. Not a big approximation. 1120 1120 !! To remove this approximation, copy lines below inside barotropic loop 1121 !! and update depths at T- F points (ht and zhf resp.) at each barotropic time step1121 !! and update depths at T- points (ht) at each barotropic time step 1122 1122 !! 1123 1123 !! Compute zwz = f / ( height of the water colomn ) … … 1126 1126 INTEGER :: ji ,jj, jk ! dummy loop indices 1127 1127 REAL(wp) :: z1_ht 1128 REAL(wp), DIMENSION(jpi,jpj) :: zhf1129 1128 !!---------------------------------------------------------------------- 1130 1129 ! 1131 1130 SELECT CASE( nvor_scheme ) 1132 CASE( np_EEN ) != EEN scheme using e3f energy & enstrophy scheme1133 SELECT CASE( nn_e en_e3f )!* ff_f/e3 at F-point1131 CASE( np_EEN, np_ENE, np_ENS , np_MIX ) != schemes using the same e3f definition 1132 SELECT CASE( nn_e3f_typ ) !* ff_f/e3 at F-point 1134 1133 CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4) 1135 DO_2D( 1, 0, 1, 0 )1134 DO_2D( 0, 0, 0, 0 ) 1136 1135 zwz(ji,jj) = ( ht(ji,jj+1) + ht(ji+1,jj+1) & 1137 1136 & + ht(ji,jj ) + ht(ji+1,jj ) ) * 0.25_wp … … 1139 1138 END_2D 1140 1139 CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask) 1141 DO_2D( 1, 0, 1, 0 )1140 DO_2D( 0, 0, 0, 0 ) 1142 1141 zwz(ji,jj) = ( ht(ji,jj+1) + ht(ji+1,jj+1) & 1143 1142 & + ht(ji,jj ) + ht(ji+1,jj ) ) & … … 1148 1147 END SELECT 1149 1148 CALL lbc_lnk( 'dynspg_ts', zwz, 'F', 1._wp ) 1149 END SELECT 1150 ! 1151 SELECT CASE( nvor_scheme ) 1152 CASE( np_EEN ) 1150 1153 ! 1151 1154 ftne(1,:) = 0._wp ; ftnw(1,:) = 0._wp ; ftse(1,:) = 0._wp ; ftsw(1,:) = 0._wp … … 1157 1160 END_2D 1158 1161 ! 1159 CASE( np_EET ) != EEN scheme using e3t energy conserving scheme1162 CASE( np_EET ) != EEN scheme using e3t energy conserving scheme 1160 1163 ftne(1,:) = 0._wp ; ftnw(1,:) = 0._wp ; ftse(1,:) = 0._wp ; ftsw(1,:) = 0._wp 1161 1164 DO_2D( 0, 1, 0, 1 ) … … 1167 1170 END_2D 1168 1171 ! 1169 CASE( np_ENE, np_ENS , np_MIX ) != all other schemes (ENE, ENS, MIX) except ENT !1170 !1171 zwz(:,:) = 0._wp1172 zhf(:,:) = 0._wp1173 1174 !!gm assume 0 in both cases (which is almost surely WRONG ! ) as hvatf has been removed1175 !!gm A priori a better value should be something like :1176 !!gm zhf(i,j) = masked sum of ht(i,j) , ht(i+1,j) , ht(i,j+1) , (i+1,j+1)1177 !!gm divided by the sum of the corresponding mask1178 !!gm1179 !!1180 IF( .NOT.ln_sco ) THEN1181 1182 !!gm agree the JC comment : this should be done in a much clear way1183 1184 ! JC: It not clear yet what should be the depth at f-points over land in z-coordinate case1185 ! Set it to zero for the time being1186 ! IF( rn_hmin < 0._wp ) THEN ; jk = - INT( rn_hmin ) ! from a nb of level1187 ! ELSE ; jk = MINLOC( gdepw_0, mask = gdepw_0 > rn_hmin, dim = 1 ) ! from a depth1188 ! ENDIF1189 ! zhf(:,:) = gdepw_0(:,:,jk+1)1190 !1191 ELSE1192 !1193 !zhf(:,:) = hbatf(:,:)1194 DO_2D( 1, 0, 1, 0 )1195 zhf(ji,jj) = ( ht_0 (ji,jj ) + ht_0 (ji+1,jj ) &1196 & + ht_0 (ji,jj+1) + ht_0 (ji+1,jj+1) ) &1197 & / MAX( ssmask(ji,jj ) + ssmask(ji+1,jj ) &1198 & + ssmask(ji,jj+1) + ssmask(ji+1,jj+1) , 1._wp )1199 END_2D1200 ENDIF1201 !1202 DO jj = 1, jpjm1 ! keep only the value at the coast if sco1203 zhf(:,jj) = zhf(:,jj) * (1._wp- umask(:,jj,1) * umask(:,jj+1,1))1204 END DO1205 !1206 DO jk = 1, jpkm1 ! ocean point : sum of masked e3f1207 DO jj = 1, jpjm11208 zhf(:,jj) = zhf(:,jj) + e3f(:,jj,jk) * umask(:,jj,jk) * umask(:,jj+1,jk)1209 END DO1210 END DO1211 CALL lbc_lnk( 'dynspg_ts', zhf, 'F', 1._wp )1212 ! JC: TBC. hf should be greater than 01213 DO_2D( 1, 1, 1, 1 )1214 IF( zhf(ji,jj) /= 0._wp ) zwz(ji,jj) = 1._wp / zhf(ji,jj)1215 END_2D1216 zwz(:,:) = ff_f(:,:) * zwz(:,:)1217 1172 END SELECT 1218 1173 1219 1174 END SUBROUTINE dyn_cor_2d_init 1220 1221 1175 1222 1176 … … 1412 1366 END SUBROUTINE wad_spg 1413 1367 1414 1415 1368 1416 1369 SUBROUTINE dyn_drg_init( Kbb, Kmm, puu, pvv, puu_b ,pvv_b, pu_RHSi, pv_RHSi, pCdU_u, pCdU_v ) -
NEMO/branches/2020/dev_r13327_KERNEL-06_2_techene_e3/src/OCE/DYN/dynvor.F90
r13644 r13696 60 60 LOGICAL, PUBLIC :: ln_dynvor_eeT !: t-point energy conserving scheme (EET) 61 61 LOGICAL, PUBLIC :: ln_dynvor_een !: energy & enstrophy conserving scheme (EEN) 62 INTEGER, PUBLIC :: nn_een_e3f !: e3f=masked averaging of e3t divided by 4 (=0) or by the sum of mask (=1)63 62 LOGICAL, PUBLIC :: ln_dynvor_mix !: mixed scheme (MIX) 64 63 LOGICAL, PUBLIC :: ln_dynvor_msk !: vorticity multiplied by fmask (=T) or not (=F) (all vorticity schemes) 64 INTEGER, PUBLIC :: nn_e3f_typ !: e3f=masked averaging of e3t divided by 4 (=0) or by the sum of mask (=1) 65 65 66 66 INTEGER, PUBLIC :: nvor_scheme !: choice of the type of advection scheme … … 84 84 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: dj_e1v_2 ! = dj(e1v)/2 - - - - 85 85 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: di_e2v_2e1e2f ! = di(e2u)/(2*e1e2f) used in F-point metric term calculation 86 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: dj_e1u_2e1e2f ! = dj(e1v)/(2*e1e2f) - - - - 86 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: dj_e1u_2e1e2f ! = dj(e1v)/(2*e1e2f) - - - - 87 ! 88 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: e3f_0vor ! e3f used in EEN, ENE and ENS cases (key_qco only) 87 89 88 90 REAL(wp) :: r1_4 = 0.250_wp ! =1/4 … … 335 337 ! 336 338 INTEGER :: ji, jj, jk ! dummy loop indices 337 REAL(wp) :: zx1, zy1, zx2, zy2 ! local scalars339 REAL(wp) :: zx1, zy1, zx2, zy2, ze3f, zmsk ! local scalars 338 340 REAL(wp), DIMENSION(jpi,jpj) :: zwx, zwy, zwz ! 2D workspace 339 341 !!---------------------------------------------------------------------- … … 387 389 ! 388 390 #if defined key_qco 389 ! !== horizontal fluxes and potential vorticity ==! 391 DO_2D( 1, 0, 1, 0 ) !== potential vorticity ==! (key_qco) 392 zwz(ji,jj) = zwz(ji,jj) / e3f_vor(ji,jj,jk) 393 END_2D 394 #else 395 SELECT CASE( nn_e3f_typ ) !== potential vorticity ==! 396 CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4) 397 DO_2D( 1, 0, 1, 0 ) 398 ze3f = ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) & 399 & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) & 400 & + e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) & 401 & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) ) 402 IF( ze3f /= 0._wp ) THEN ; zwz(ji,jj) = zwz(ji,jj) * 4._wp / ze3f 403 ELSE ; zwz(ji,jj) = 0._wp 404 ENDIF 405 END_2D 406 CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask) 407 DO_2D( 1, 0, 1, 0 ) 408 ze3f = ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) & 409 & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) & 410 & + e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) & 411 & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) ) 412 zmsk = ( tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) & 413 & + tmask(ji,jj ,jk) + tmask(ji+1,jj ,jk) ) 414 IF( ze3f /= 0._wp ) THEN ; zwz(ji,jj) = zwz(ji,jj) * zmsk / ze3f 415 ELSE ; zwz(ji,jj) = 0._wp 416 ENDIF 417 END_2D 418 END SELECT 419 #endif 420 ! !== horizontal fluxes ==! 390 421 zwx(:,:) = e2u(:,:) * e3u(:,:,jk,Kmm) * pu(:,:,jk) 391 422 zwy(:,:) = e1v(:,:) * e3v(:,:,jk,Kmm) * pv(:,:,jk) 392 zwz(:,:) = zwz(:,:) / e3f(:,:,jk)393 #else394 ! !== horizontal fluxes and potential vorticity ==!395 IF( ln_sco ) THEN ! weighted by e3396 zwz(:,:) = zwz(:,:) / e3f(:,:,jk)397 zwx(:,:) = e2u(:,:) * e3u(:,:,jk,Kmm) * pu(:,:,jk)398 zwy(:,:) = e1v(:,:) * e3v(:,:,jk,Kmm) * pv(:,:,jk)399 ELSE ! ugly fix: unweighted by e3400 zwx(:,:) = e2u(:,:) * pu(:,:,jk)401 zwy(:,:) = e1v(:,:) * pv(:,:,jk)402 ENDIF403 #endif404 423 ! 405 424 ! !== compute and add the vorticity term trend =! … … 445 464 ! 446 465 INTEGER :: ji, jj, jk ! dummy loop indices 447 REAL(wp) :: zuav, zvau ! local scalars466 REAL(wp) :: zuav, zvau, ze3f, zmsk ! local scalars 448 467 REAL(wp), DIMENSION(jpi,jpj) :: zwx, zwy, zwz, zww ! 2D workspace 449 468 !!---------------------------------------------------------------------- … … 497 516 ! 498 517 #if defined key_qco 499 ! !== horizontal fluxes and potential vorticity ==! 518 DO_2D( 1, 0, 1, 0 ) !== potential vorticity ==! (key_qco) 519 zwz(ji,jj) = zwz(ji,jj) / e3f_vor(ji,jj,jk) 520 END_2D 521 #else 522 SELECT CASE( nn_e3f_typ ) !== potential vorticity ==! 523 CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4) 524 DO_2D( 1, 0, 1, 0 ) 525 ze3f = ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) & 526 & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) & 527 & + e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) & 528 & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) ) 529 IF( ze3f /= 0._wp ) THEN ; zwz(ji,jj) = zwz(ji,jj) * 4._wp / ze3f 530 ELSE ; zwz(ji,jj) = 0._wp 531 ENDIF 532 END_2D 533 CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask) 534 DO_2D( 1, 0, 1, 0 ) 535 ze3f = ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) & 536 & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) & 537 & + e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) & 538 & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) ) 539 zmsk = ( tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) & 540 & + tmask(ji,jj ,jk) + tmask(ji+1,jj ,jk) ) 541 IF( ze3f /= 0._wp ) THEN ; zwz(ji,jj) = zwz(ji,jj) * zmsk / ze3f 542 ELSE ; zwz(ji,jj) = 0._wp 543 ENDIF 544 END_2D 545 END SELECT 546 #endif 547 ! !== horizontal fluxes ==! 500 548 zwx(:,:) = e2u(:,:) * e3u(:,:,jk,Kmm) * pu(:,:,jk) 501 549 zwy(:,:) = e1v(:,:) * e3v(:,:,jk,Kmm) * pv(:,:,jk) 502 zwz(:,:) = zwz(:,:) / e3f(:,:,jk)503 #else504 ! !== horizontal fluxes and potential vorticity ==!505 IF( ln_sco ) THEN ! weighted by e3506 zwz(:,:) = zwz(:,:) / e3f(:,:,jk)507 zwx(:,:) = e2u(:,:) * e3u(:,:,jk,Kmm) * pu(:,:,jk)508 zwy(:,:) = e1v(:,:) * e3v(:,:,jk,Kmm) * pv(:,:,jk)509 ELSE ! ugly fix: unweighted by e3510 zwx(:,:) = e2u(:,:) * pu(:,:,jk)511 zwy(:,:) = e1v(:,:) * pv(:,:,jk)512 ENDIF513 #endif514 550 ! 515 551 ! !== compute and add the vorticity term trend =! … … 570 606 ! ! =============== 571 607 ! 572 SELECT CASE( nn_een_e3f ) ! == reciprocal of e3 at F-point 608 #if defined key_qco 609 DO_2D( 1, 0, 1, 0 ) ! == reciprocal of e3 at F-point (key_qco) 610 z1_e3f(ji,jj) = 1._wp / e3f_vor(ji,jj,jk) 611 END_2D 612 #else 613 SELECT CASE( nn_e3f_typ ) ! == reciprocal of e3 at F-point 573 614 CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4) 574 615 DO_2D( 1, 0, 1, 0 ) … … 594 635 END_2D 595 636 END SELECT 637 #endif 596 638 ! 597 639 SELECT CASE( kvor ) !== vorticity considered ==! … … 801 843 INTEGER :: ji, jj, jk ! dummy loop indices 802 844 INTEGER :: ioptio, ios ! local integer 845 REAL(wp) :: zmsk ! local scalars 803 846 !! 804 847 NAMELIST/namdyn_vor/ ln_dynvor_ens, ln_dynvor_ene, ln_dynvor_enT, ln_dynvor_eeT, & 805 & ln_dynvor_een, nn_e en_e3f, ln_dynvor_mix, ln_dynvor_msk848 & ln_dynvor_een, nn_e3f_typ , ln_dynvor_mix, ln_dynvor_msk 806 849 !!---------------------------------------------------------------------- 807 850 ! … … 825 868 WRITE(numout,*) ' energy conserving scheme (een using e3t) ln_dynvor_eeT = ', ln_dynvor_eeT 826 869 WRITE(numout,*) ' enstrophy and energy conserving scheme ln_dynvor_een = ', ln_dynvor_een 827 WRITE(numout,*) ' e3f = averaging /4 (=0) or /sum(tmask) (=1) nn_e en_e3f = ', nn_een_e3f870 WRITE(numout,*) ' e3f = averaging /4 (=0) or /sum(tmask) (=1) nn_e3f_typ = ', nn_e3f_typ 828 871 WRITE(numout,*) ' mixed enstrophy/energy conserving scheme ln_dynvor_mix = ', ln_dynvor_mix 829 872 WRITE(numout,*) ' masked (=T) or unmasked(=F) vorticity ln_dynvor_msk = ', ln_dynvor_msk … … 891 934 ! 892 935 END SELECT 893 936 !!st ADD !! pense a recalculer le hf_0 937 #if defined key_qco 938 SELECT CASE( nvor_scheme ) ! qco case: pre-computed a specific e3f_0 for some vorticity schemes 939 CASE( np_ENS , np_ENE , np_EEN , np_MIX ) 940 ! 941 ALLOCATE( e3f_0vor(jpi,jpj,jpk) ) 942 ! 943 SELECT CASE( nn_e3f_typ ) 944 CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4) 945 DO_3D( 0, 0, 0, 0, 1, jpk ) 946 e3f_0vor(ji,jj,jk) = ( e3t_0(ji ,jj+1,jk)*tmask(ji ,jj+1,jk) & 947 & + e3t_0(ji+1,jj+1,jk)*tmask(ji+1,jj+1,jk) & 948 & + e3t_0(ji ,jj ,jk)*tmask(ji ,jj ,jk) & 949 & + e3t_0(ji+1,jj ,jk)*tmask(ji+1,jj ,jk) ) * 0.25_wp 950 END_3D 951 CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask) 952 DO_3D( 0, 0, 0, 0, 1, jpk ) 953 zmsk = (tmask(ji,jj+1,jk) +tmask(ji+1,jj+1,jk) & 954 & + tmask(ji,jj ,jk) +tmask(ji+1,jj ,jk) ) 955 ! 956 IF( zmsk /= 0._wp ) THEN 957 e3f_0vor(ji,jj,jk) = ( e3t_0(ji ,jj+1,jk)*tmask(ji ,jj+1,jk) & 958 & + e3t_0(ji+1,jj+1,jk)*tmask(ji+1,jj+1,jk) & 959 & + e3t_0(ji ,jj ,jk)*tmask(ji ,jj ,jk) & 960 & + e3t_0(ji+1,jj ,jk)*tmask(ji+1,jj ,jk) ) / zmsk 961 ENDIF 962 END_3D 963 END SELECT 964 ! 965 CALL lbc_lnk( 'dynvor', e3f_0vor, 'F', 1._wp ) 966 ! ! insure e3f_0vor /= 0 967 WHERE( e3f_0vor(:,:,:) == 0._wp ) e3f_0vor(:,:,:) = e3f_0(:,:,:) 968 ! 969 END SELECT 970 ! 971 #endif 972 !!st END 894 973 IF(lwp) THEN ! Print the choice 895 974 WRITE(numout,*)
Note: See TracChangeset
for help on using the changeset viewer.
