From 7f01bce57ac35e4478b5306e2fe1fb7b075b4860 Mon Sep 17 00:00:00 2001
From: VIE Benoit <vie@sxphynh>
Date: Thu, 30 Nov 2023 12:20:15 +0100
Subject: [PATCH] Turbulence now changes rr, rs, rg, rh - 1 test case ok
---
src/MNH/ini_budget.f90 | 70 +++++++++++++++++++++++++++++++----------
src/PHYEX/turb/turb.f90 | 64 ++++++++++++++++++++++++++++---------
2 files changed, 104 insertions(+), 30 deletions(-)
diff --git a/src/MNH/ini_budget.f90 b/src/MNH/ini_budget.f90
index a5ed1e2d4..93bc1ab1e 100644
--- a/src/MNH/ini_budget.f90
+++ b/src/MNH/ini_budget.f90
@@ -1964,9 +1964,20 @@ if ( tbudgets(NBUDGET_RR)%lenabled ) then
tzsource%lavailable = ohorelax_rr
call Budget_source_add( tbudgets(NBUDGET_RR), tzsource )
+ tzsource%cmnhname = 'VTURB'
+ tzsource%clongname = 'vertical turbulent diffusion'
+ tzsource%lavailable = hturb == 'TKEL'
+ call Budget_source_add( tbudgets(NBUDGET_RR), tzsource )
+
+ tzsource%cmnhname = 'HTURB'
+ tzsource%clongname = 'horizontal turbulent diffusion'
+ tzsource%lavailable = hturb == 'TKEL' .and. HTURBDIM == '3DIM'
+ call Budget_source_add( tbudgets(NBUDGET_RR), tzsource )
+
tzsource%cmnhname = 'NETUR'
tzsource%clongname = 'negativity correction induced by turbulence'
- tzsource%lavailable = hturb == 'TKEL' .and. ( hcloud == 'KHKO' .or. hcloud == 'C2R2' .or. hcloud == 'LIMA' )
+ tzsource%lavailable = hturb == 'TKEL' .and. ( hcloud == 'ICE3' .or. hcloud == 'ICE4' &
+ .or. hcloud == 'KHKO' .or. hcloud == 'C2R2' .or. hcloud == 'LIMA' )
call Budget_source_add( tbudgets(NBUDGET_RR), tzsource )
tzsource%cmnhname = 'VISC'
@@ -2451,11 +2462,20 @@ if ( tbudgets(NBUDGET_RS)%lenabled ) then
tzsource%lavailable = ohorelax_rs
call Budget_source_add( tbudgets(NBUDGET_RS), tzsource )
-! tzsource%cmnhname = 'NETUR'
-! tzsource%clongname = 'negativity correction induced by turbulence'
-! tzsource%lavailable = hturb == 'TKEL' .and. ( hcloud == 'ICE3' .or. hcloud == 'ICE4' &
-! .or. hcloud == 'KHKO' .or. hcloud == 'C2R2' .or. hcloud == 'LIMA' )
-! call Budget_source_add( tbudgets(NBUDGET_RS), tzsource nneturrs )
+ tzsource%cmnhname = 'VTURB'
+ tzsource%clongname = 'vertical turbulent diffusion'
+ tzsource%lavailable = hturb == 'TKEL'
+ call Budget_source_add( tbudgets(NBUDGET_RS), tzsource )
+
+ tzsource%cmnhname = 'HTURB'
+ tzsource%clongname = 'horizontal turbulent diffusion'
+ tzsource%lavailable = hturb == 'TKEL' .and. HTURBDIM == '3DIM'
+ call Budget_source_add( tbudgets(NBUDGET_RS), tzsource )
+
+ tzsource%cmnhname = 'NETUR'
+ tzsource%clongname = 'negativity correction induced by turbulence'
+ tzsource%lavailable = hturb == 'TKEL' .and. ( hcloud == 'ICE3' .or. hcloud == 'ICE4' .or. hcloud == 'LIMA' )
+ call Budget_source_add( tbudgets(NBUDGET_RS), tzsource )
tzsource%cmnhname = 'VISC'
tzsource%clongname = 'viscosity'
@@ -2642,11 +2662,20 @@ if ( tbudgets(NBUDGET_RG)%lenabled ) then
tzsource%lavailable = ohorelax_rg
call Budget_source_add( tbudgets(NBUDGET_RG), tzsource )
-! tzsource%cmnhname = 'NETUR'
-! tzsource%clongname = 'negativity correction induced by turbulence'
-! tzsource%lavailable = hturb == 'TKEL' .and. ( hcloud == 'ICE3' .or. hcloud == 'ICE4' &
-! .or. hcloud == 'KHKO' .or. hcloud == 'C2R2' .or. hcloud == 'LIMA' )
-! call Budget_source_add( tbudgets(NBUDGET_RG), tzsource nneturrg )
+ tzsource%cmnhname = 'VTURB'
+ tzsource%clongname = 'vertical turbulent diffusion'
+ tzsource%lavailable = hturb == 'TKEL'
+ call Budget_source_add( tbudgets(NBUDGET_RG), tzsource )
+
+ tzsource%cmnhname = 'HTURB'
+ tzsource%clongname = 'horizontal turbulent diffusion'
+ tzsource%lavailable = hturb == 'TKEL' .and. HTURBDIM == '3DIM'
+ call Budget_source_add( tbudgets(NBUDGET_RG), tzsource )
+
+ tzsource%cmnhname = 'NETUR'
+ tzsource%clongname = 'negativity correction induced by turbulence'
+ tzsource%lavailable = hturb == 'TKEL' .and. ( hcloud == 'ICE3' .or. hcloud == 'ICE4' .or. hcloud == 'LIMA' )
+ call Budget_source_add( tbudgets(NBUDGET_RG), tzsource )
tzsource%cmnhname = 'VISC'
tzsource%clongname = 'viscosity'
@@ -2851,11 +2880,20 @@ if ( tbudgets(NBUDGET_RH)%lenabled ) then
tzsource%lavailable = ohorelax_rh
call Budget_source_add( tbudgets(NBUDGET_RH), tzsource )
-! tzsource%cmnhname = 'NETUR'
-! tzsource%clongname = 'negativity correction induced by turbulence'
-! tzsource%lavailable = hturb == 'TKEL' .and. ( hcloud == 'ICE3' .or. hcloud == 'ICE4' &
-! .or. hcloud == 'KHKO' .or. hcloud == 'C2R2' .or. hcloud == 'LIMA' )
-! call Budget_source_add( tbudgets(NBUDGET_RH), tzsource nneturrh )
+ tzsource%cmnhname = 'VTURB'
+ tzsource%clongname = 'vertical turbulent diffusion'
+ tzsource%lavailable = hturb == 'TKEL'
+ call Budget_source_add( tbudgets(NBUDGET_RH), tzsource )
+
+ tzsource%cmnhname = 'HTURB'
+ tzsource%clongname = 'horizontal turbulent diffusion'
+ tzsource%lavailable = hturb == 'TKEL' .and. HTURBDIM == '3DIM'
+ call Budget_source_add( tbudgets(NBUDGET_RH), tzsource )
+
+ tzsource%cmnhname = 'NETUR'
+ tzsource%clongname = 'negativity correction induced by turbulence'
+ tzsource%lavailable = hturb == 'TKEL' .and. ( hcloud == 'ICE3' .or. hcloud == 'ICE4' .or. hcloud == 'LIMA' )
+ call Budget_source_add( tbudgets(NBUDGET_RH), tzsource )
tzsource%cmnhname = 'VISC'
tzsource%clongname = 'viscosity'
diff --git a/src/PHYEX/turb/turb.f90 b/src/PHYEX/turb/turb.f90
index 42b1d1cb8..79eb73fe6 100644
--- a/src/PHYEX/turb/turb.f90
+++ b/src/PHYEX/turb/turb.f90
@@ -502,6 +502,11 @@ REAL :: ZALPHA ! work coefficient :
REAL :: ZTIME1, ZTIME2
TYPE(TFIELDMETADATA) :: TZFIELD
!
+REAL, DIMENSION(D%NIJT,D%NKT,KSV+KRR) :: ZWORKT
+REAL, DIMENSION(D%NIJT,D%NKT,KSV+KRR) :: ZWORKS
+REAL, DIMENSION(D%NIJT, KSV+KRR) :: ZWORKSFSV
+REAL, DIMENSION(D%NIJT,D%NKT,KSV+KRR) :: ZWORKWSV
+!
!* 1.PRELIMINARIES
! -------------
!
@@ -541,6 +546,13 @@ END IF
!Save LIMA scalar variables sources
ZRSVS(:,:,1:KSV)=PRSVS(:,:,1:KSV)
!
+ZWORKT(:,:,1:KSV)=PSVT(:,:,:)
+ZWORKS(:,:,1:KSV)=PRSVS(:,:,:)
+ZWORKT(:,:,KSV+1:KSV+KRR)=PRT(:,:,:)
+ZWORKS(:,:,KSV+1:KSV+KRR)=PRRS(:,:,:)
+ZWORKSFSV(:,:)=0.
+ZWORKWSV(:,:,:)=0.
+ZWORKSFSV(:,1:KSV)=PSFSV(:,:)
!
!----------------------------------------------------------------------------
!
@@ -991,33 +1003,37 @@ IF( BUCONF%LBUDGET_RV ) THEN
END IF
IF( BUCONF%LBUDGET_RC ) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RC), 'VTURB', PRRS (:,:, 2) )
+IF( BUCONF%LBUDGET_RR ) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RR), 'VTURB', ZWORKS (:,:, KSV + 3) )
IF( BUCONF%LBUDGET_RI ) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RI), 'VTURB', PRRS (:,:, 4) )
+IF( BUCONF%LBUDGET_RS ) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RS), 'VTURB', ZWORKS (:,:, KSV + 5) )
+IF( BUCONF%LBUDGET_RG ) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RG), 'VTURB', ZWORKS (:,:, KSV + 6) )
+IF( BUCONF%LBUDGET_RH ) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RH), 'VTURB', ZWORKS (:,:, KSV + 7) )
IF( BUCONF%LBUDGET_SV ) THEN
DO JSV = 1, KSV
- CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_SV1 - 1 + JSV), 'VTURB', PRSVS(:,:, JSV) )
+ CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_SV1 - 1 + JSV), 'VTURB', ZWORKS(:,:, JSV) )
END DO
END IF
CALL TURB_VER(D,CST,CSTURB,TURBN,NEBN,TLES, &
KRR,KRRL,KRRI,KGRADIENTS, &
OOCEAN, ODEEPOC, OCOMPUTE_SRC, &
- KSV,KSV_LGBEG,KSV_LGEND, &
+ KSV+KRR,KSV_LGBEG,KSV_LGEND, &
ZEXPL, O2D, ONOMIXLG, OFLAT, &
OCOUPLES,OBLOWSNOW,OFLYER, PRSNOW, &
PTSTEP,TPFILE, &
PDXX,PDYY,PDZZ,PDZX,PDZY,PDIRCOSZW,PZZ, &
PCOSSLOPE,PSINSLOPE, &
PRHODJ,PTHVREF,PSFU,PSFV, &
- PSFTH,PSFRV,PSFSV,PSFTH,PSFRV,PSFSV, &
+ PSFTH,PSFRV,ZWORKSFSV,PSFTH,PSFRV,ZWORKSFSV, &
ZCDUEFF,ZTAU11M,ZTAU12M,ZTAU33M, &
- PUT,PVT,PWT,ZUSLOPE,ZVSLOPE,PTHLT,PRT,PSVT, &
+ PUT,PVT,PWT,ZUSLOPE,ZVSLOPE,PTHLT,PRT,ZWORKT, &
PTKET,ZLM,PLENGTHM,PLENGTHH,ZLEPS,MFMOIST, &
ZLOCPEXNM,ZATHETA,ZAMOIST,PSRCT,ZFRAC_ICE, &
ZFWTH,ZFWR,ZFTH2,ZFR2,ZFTHR,PBL_DEPTH, &
PSBL_DEPTH,ZLMO,PHGRAD,PZS, &
- PRUS,PRVS,PRWS,PRTHLS,PRRS,PRSVS, &
- PDP,PTP,PSIGS,PWTH,PWRC,PWSV, &
+ PRUS,PRVS,PRWS,PRTHLS,PRRS,ZWORKS, &
+ PDP,PTP,PSIGS,PWTH,PWRC,ZWORKWSV, &
PSSTFL, PSSTFL_C, PSSRFL_C,PSSUFL_C,PSSVFL_C, &
PSSUFL,PSSVFL )
@@ -1054,11 +1070,15 @@ IF( BUCONF%LBUDGET_RV ) THEN
END IF
IF( BUCONF%LBUDGET_RC ) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RC), 'VTURB', PRRS(:,:, 2) )
+IF( BUCONF%LBUDGET_RR ) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RR), 'VTURB', ZWORKS(:,:, KSV + 3) )
IF( BUCONF%LBUDGET_RI ) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RI), 'VTURB', PRRS(:,:, 4) )
+IF( BUCONF%LBUDGET_RS ) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RS), 'VTURB', ZWORKS(:,:, KSV + 5) )
+IF( BUCONF%LBUDGET_RG ) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RG), 'VTURB', ZWORKS(:,:, KSV + 6) )
+IF( BUCONF%LBUDGET_RH ) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RH), 'VTURB', ZWORKS(:,:, KSV + 7) )
IF( BUCONF%LBUDGET_SV ) THEN
DO JSV = 1, KSV
- CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_SV1 - 1 + JSV), 'VTURB', PRSVS(:,:, JSV) )
+ CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_SV1 - 1 + JSV), 'VTURB', ZWORKS(:,:, JSV) )
END DO
END IF
!
@@ -1089,15 +1109,20 @@ IF( TURBN%CTURBDIM == '3DIM' ) THEN
END IF
IF( BUCONF%LBUDGET_RC ) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RC), 'HTURB', PRRS(:,:, 2) )
+ IF( BUCONF%LBUDGET_RR ) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RR), 'HTURB', ZWORKS(:,:, KSV+3) )
IF( BUCONF%LBUDGET_RI ) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RI), 'HTURB', PRRS(:,:, 4) )
+ IF( BUCONF%LBUDGET_RS ) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RS), 'HTURB', ZWORKS(:,:, KSV+5) )
+ IF( BUCONF%LBUDGET_RG ) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RG), 'HTURB', ZWORKS(:,:, KSV+6) )
+ IF( BUCONF%LBUDGET_RH ) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RH), 'HTURB', ZWORKS(:,:, KSV+7) )
IF( BUCONF%LBUDGET_SV ) THEN
DO JSV = 1, KSV
- CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_SV1 - 1 + JSV), 'HTURB', PRSVS(:,:, JSV) )
+ CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_SV1 - 1 + JSV), 'HTURB', ZWORKS(:,:, JSV) )
END DO
END IF
CALL TURB_HOR_SPLT(D,CST,CSTURB, TURBN, NEBN, TLES, &
KSPLIT, KRR, KRRL, KRRI, KSV,KSV_LGBEG,KSV_LGEND, &
+ KSPLIT, KRR, KRRL, KRRI, KSV+KRR,KSV_LGBEG,KSV_LGEND,&
PTSTEP,HLBCX,HLBCY, OFLAT,O2D, ONOMIXLG, &
OOCEAN,OCOMPUTE_SRC,OBLOWSNOW,PRSNOW, &
TPFILE, KHALO, &
@@ -1105,14 +1130,14 @@ IF( TURBN%CTURBDIM == '3DIM' ) THEN
PDIRCOSXW,PDIRCOSYW,PDIRCOSZW, &
PCOSSLOPE,PSINSLOPE, &
PRHODJ,PTHVREF, &
- PSFTH,PSFRV,PSFSV, &
+ PSFTH,PSFRV,ZWORKSFSV, &
ZCDUEFF,ZTAU11M,ZTAU12M,ZTAU22M,ZTAU33M, &
- PUT,PVT,PWT,ZUSLOPE,ZVSLOPE,PTHLT,PRT,PSVT, &
+ PUT,PVT,PWT,ZUSLOPE,ZVSLOPE,PTHLT,PRT,ZWORKT, &
PTKET,ZLM,ZLEPS, &
ZLOCPEXNM,ZATHETA,ZAMOIST,PSRCT,ZFRAC_ICE, &
PDP,PTP,PSIGS, &
ZTRH, &
- PRUS,PRVS,PRWS,PRTHLS,PRRS,PRSVS )
+ PRUS,PRVS,PRWS,PRTHLS,PRRS,ZWORKS )
!
! IF (HCLOUD == 'LIMA') THEN
! IF (KSV_LIMA_NR.GT.0) PRSVS(:,:,KSV_LIMA_NR) = ZRSVS(:,:,KSV_LIMA_NR)
@@ -1147,11 +1172,15 @@ IF( TURBN%CTURBDIM == '3DIM' ) THEN
END IF
IF( BUCONF%LBUDGET_RC ) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RC), 'HTURB', PRRS(:,:, 2) )
+ IF( BUCONF%LBUDGET_RR ) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RR), 'HTURB', ZWORKS(:,:, KSV+3) )
IF( BUCONF%LBUDGET_RI ) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RI), 'HTURB', PRRS(:,:, 4) )
+ IF( BUCONF%LBUDGET_RS ) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RS), 'HTURB', ZWORKS(:,:, KSV+5) )
+ IF( BUCONF%LBUDGET_RG ) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RG), 'HTURB', ZWORKS(:,:, KSV+6) )
+ IF( BUCONF%LBUDGET_RH ) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RH), 'HTURB', ZWORKS(:,:, KSV+7) )
IF( BUCONF%LBUDGET_SV ) THEN
DO JSV = 1, KSV
- CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_SV1 - 1 + JSV), 'HTURB', PRSVS(:,:, JSV) )
+ CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_SV1 - 1 + JSV), 'HTURB', ZWORKS(:,:, JSV) )
END DO
END IF
END IF
@@ -1314,8 +1343,15 @@ IF ( KRRL >= 1 ) THEN
* PRRS(:,:,2)
!$mnh_end_expand_array(JIJ=IIJB:IIJE,JK=1:IKT)
END IF
-END IF
-
+END IF!
+!
+PRSVS(:,:,:)=ZWORKS(:,:,1:KSV)
+IF (KRR.GE.3) PRRS(:,:,3)=ZWORKS(:,:,KSV+3)
+IF (KRR.GE.5) PRRS(:,:,5)=ZWORKS(:,:,KSV+5)
+IF (KRR.GE.6) PRRS(:,:,6)=ZWORKS(:,:,KSV+6)
+IF (KRR.GE.7) PRRS(:,:,7)=ZWORKS(:,:,KSV+7)
+IF (OFLYER) PWSV(:,:,:)=ZWORKWSV(:,:,1:KSV)
+!
! Remove non-physical negative values (unnecessary in a perfect world) + corresponding budgets
CALL SOURCES_NEG_CORRECT_PHY(D,KSV,HCLOUD,HELEC,'NETUR',KRR,PTSTEP,PPABST,PTHLT,PRT,PRTHLS,PRRS,PRSVS)
!----------------------------------------------------------------------------
--
GitLab