Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Adding Integrated Canopy Effects to CCPP/PBL Scheme #253

Open
wants to merge 55 commits into
base: ufs/dev
Choose a base branch
from
+346 −4
Open

Adding Integrated Canopy Effects to CCPP/PBL Scheme #253

Show file tree
Hide file tree
Changes from 54 commits
Commits
Show all changes
55 commits
Select commit Hold shift + click to select a range
67b693e
Testing hedmf change.
drnimbusrain Sep 28, 2022
e2c7f47
Rolled back hedmf test change.
drnimbusrain Sep 29, 2022
1fb9060
Merge branch 'ufs-community:ufs/dev' into feature/aqm_canopy
drnimbusrain Dec 5, 2022
d825523
First draft commit of canopy turbulence effects in hedmf.f
drnimbusrain Dec 5, 2022
574ad4e
Fixed non-ascii character
drnimbusrain Dec 5, 2022
1665715
Changed from canopy in hedmf to satmedmvdifq for GFSv16.
drnimbusrain Dec 6, 2022
1c587a3
Added more arguments for canopy varibles.
drnimbusrain Dec 6, 2022
fe808c6
Fixed non-ascii character.
drnimbusrain Dec 6, 2022
634ebb1
Added canopy height grid cell conidtion.
drnimbusrain Dec 12, 2022
08dd5a6
Updated to modify TKE instead of K for canopy.
drnimbusrain Feb 7, 2023
22d3fd7
Rolled back to modifying K directly for canopy effects.
drnimbusrain Feb 7, 2023
936a28e
Updated Meta data
drnimbusrain Feb 7, 2023
989ed53
Added lai and vegtype to the meta
drnimbusrain Feb 9, 2023
99ba757
Fixed bugs.
drnimbusrain Feb 9, 2023
de1c0b4
Added canopy_utils_mod and noahmp_tables to satmedmfvdifq meta depend…
drnimbusrain Feb 10, 2023
90762fb
Changed some comments.
drnimbusrain Feb 10, 2023
07bafc0
Fixed IF statement bug in canopy conditions.
drnimbusrain Feb 10, 2023
6b686e5
Fixed bug in canopy variable inputs.
drnimbusrain Feb 10, 2023
ea79a91
Merge remote-tracking branch 'origin/develop' into feature/aqm_canopy
drnimbusrain Feb 11, 2023
be76560
Merge branch 'ufs-community:ufs/dev' into feature/aqm_canopy
drnimbusrain Feb 11, 2023
8122173
Updated bugs.
drnimbusrain Feb 12, 2023
f99c2bb
Fixed bugs (again...).
drnimbusrain Feb 12, 2023
e68de59
Merge branch 'ufs-community:ufs/dev' into feature/aqm_canopy
drnimbusrain Mar 3, 2023
92e33b2
Added "do_canopy" flag for eddy effects in TKE-EDMF
drnimbusrain Mar 8, 2023
d9d6465
Fixed bug in syntax
drnimbusrain Mar 8, 2023
da70f9b
Added canopy effects for multiple model layers...if necessary
drnimbusrain Mar 23, 2023
3f5a1d5
Merge branch 'ufs-community:ufs/dev' into feature/aqm_canopy
drnimbusrain Mar 23, 2023
03cbe6f
Updated canopy effect to use interface levels, zi.
drnimbusrain Mar 24, 2023
b8c61a3
Merge branch 'feature/aqm_canopy' of https://github.com/noaa-oar-arl/…
drnimbusrain Mar 24, 2023
b0465a5
Merge branch 'ufs-community:ufs/dev' into feature/aqm_canopy
drnimbusrain Mar 25, 2023
a962bdc
Merge branch 'ufs-community:ufs/dev' into feature/aqm_canopy
drnimbusrain Mar 28, 2023
9490726
Merge branch 'ufs-community:ufs/dev' into feature/aqm_canopy
drnimbusrain May 31, 2023
97d0417
Merge branch 'ufs-community:ufs/dev' into feature/aqm_canopy
drnimbusrain Jun 30, 2023
df68b5c
Merge branch 'ufs-community:ufs/dev' into feature/aqm_canopy
drnimbusrain Jul 17, 2023
29033e2
Merge branch 'ufs-community:ufs/dev' into feature/aqm_canopy
drnimbusrain Nov 16, 2023
911e4f9
Added TODO comments for input canopy variables and draft meta.
drnimbusrain Nov 16, 2023
7aa45c4
Merge branch 'feature/aqm_canopy' of https://github.com/noaa-oar-arl/…
drnimbusrain Nov 16, 2023
3a9338f
Fixed canopy height table and changed ffrac to frt.
drnimbusrain Dec 19, 2023
1e28ccd
Removed noah_mp_table dependency for aqm_canopy
drnimbusrain Dec 19, 2023
b0706c2
Pass via interface 5 AQM canopy inputs and build-in diagnostic arrays…
iri01 Feb 28, 2024
7726000
Corrected canopy arrays meta names and dimensions for consistency wit…
iri01 Feb 28, 2024
4d171ce
Remove unnecessary definitions of canopy arrays
iri01 Feb 28, 2024
caee56d
Merge pull request #1 from noaa-oar-arl/feature/aqm_canopy_new
iri01 Feb 28, 2024
7a0e1e2
Fix to canopy arrays definitions: add intent
iri01 Mar 1, 2024
336447d
Merge pull request #2 from noaa-oar-arl/feature/aqm_canopy_new_fix
iri01 Mar 1, 2024
fc1fb3f
Fixed cfch and cpopu units in the canopy meta variables.
drnimbusrain Mar 3, 2024
cbcb544
Fixed standard_name for naux2d and naux3d in meta file.
drnimbusrain Mar 3, 2024
8b82cf9
Bug fix to remove FCH and divide by ZFL.
drnimbusrain Apr 4, 2024
d3d17e9
Merge branch 'feature/aqm_canopy' into feature/aqm_canopy2
drnimbusrain Feb 6, 2025
eb85a48
Removed Files.
drnimbusrain Feb 6, 2025
4065eb7
Commenting out aux arrays.
drnimbusrain Feb 10, 2025
5016165
Fixed units of cpopu for canopy work.
drnimbusrain Feb 10, 2025
d98a4d7
Removed old lai/vegtype from canopy satmedmfvdifq
drnimbusrain Feb 10, 2025
aedd1c4
Use pi physical constant from ccpp, correct units of cpopu.
iri01 Mar 1, 2025
5abb270
Bug fix on the use of the canopy module
iri01 Mar 2, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
235 changes: 232 additions & 3 deletions physics/PBL/SATMEDMF/satmedmfvdifq.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,7 @@ module satmedmfvdifq
use mfpbltq_mod
use tridi_mod
use mfscuq_mod

contains

!> \defgroup module_satmedmfvdifq GFS TKE-EDMF PBL Module
Expand Down Expand Up @@ -74,7 +75,7 @@ end subroutine satmedmfvdifq_init
!! (mfscuq.f).
!! \section detail_satmedmfvidfq GFS satmedmfvdifq Detailed Algorithm
subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
& ntiw,ntke,grav,rd,cp,rv,hvap,hfus,fv,eps,epsm1, &
& ntiw,ntke,grav,pi,rd,cp,rv,hvap,hfus,fv,eps,epsm1, &
& dv,du,tdt,rtg,u1,v1,t1,q1,usfco,vsfco,icplocn2atm, &
& swh,hlw,xmu,garea,zvfun,sigmaf, &
& psk,rbsoil,zorl,u10m,v10m,fm,fh, &
Expand All @@ -83,9 +84,16 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
& dspheat,dusfc,dvsfc,dtsfc,dqsfc,hpbl,dkt,dku, &
& kinver,xkzm_m,xkzm_h,xkzm_s,dspfac,bl_upfr,bl_dnfr, &
& rlmx,elmx,sfc_rlm,tc_pbl, &
!IVAI: canopy inputs from AQM
& do_canopy, claie, cfch, cfrt, cclu, cpopu, &
!IVAI
!TODO -Canopy Inputs from UFS
& ntqv,dtend,dtidx,index_of_temperature,index_of_x_wind, &
& index_of_y_wind,index_of_process_pbl,gen_tend,ldiag3d, &
& errmsg,errflg)
& errmsg,errflg)
!IVAI: aux arrays
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

These aux arrays are purely for diagnosing the implemention, correct? I believe they need to be commented out or removed before merging this into ufs/dev

Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thank you for your comments here. Yes, indeed the aux arrays were for diagnosing the implementation. We can comment out for now if that's acceptable path forward.

Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes please, thanks.

! & naux2d,naux3d,aux2d,aux3d)

!
use machine , only : kind_phys
use funcphys , only : fpvs
Expand All @@ -101,11 +109,17 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
integer, intent(out) :: kpbl(:)
logical, intent(in) :: gen_tend,ldiag3d
!
real(kind=kind_phys), intent(in) :: grav,rd,cp,rv,hvap,hfus,fv, &
real(kind=kind_phys), intent(in) :: grav,pi,rd,cp,rv,hvap,hfus,fv, &
& eps,epsm1
real(kind=kind_phys), intent(in) :: delt, xkzm_m, xkzm_h, xkzm_s
real(kind=kind_phys), intent(in) :: dspfac, bl_upfr, bl_dnfr
real(kind=kind_phys), intent(in) :: rlmx, elmx
!PCC CANOPY------------------------------------
logical, intent(in) :: do_canopy
!IVAI: canopy inputs
real(kind=kind_phys), intent(in) :: claie(:), cfch(:), cfrt(:),
& cclu(:), cpopu(:)
!----------------------------------------------
real(kind=kind_phys), intent(inout) :: dv(:,:), du(:,:), &
& tdt(:,:), rtg(:,:,:)
real(kind=kind_phys), intent(in) :: &
Expand Down Expand Up @@ -254,11 +268,44 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
!
real(kind=kind_phys) ck0, ck1, ch0, ch1, ce0, rchck
!
!PCC_CANOPY------------------------------------
real(kind=kind_phys) qlcr, zstblmax, hcrinv
!
real(kind=kind_phys) h1

real(kind=kind_phys) bfac, mffac

!PCC_CANOPY------------------------------------
integer COUNTCAN,KCAN
integer kount !IVAI
real(kind=kind_phys) FCH, MOL, HOL, TLCAN,
& SIGMACAN, RRCAN, BBCAN,
& AACAN, ZCAN, ZFL, BOTCAN,
& EDDYVEST1, EDDYVEST_INT

! in canopy eddy diffusivity [ m**2/s ]
real(kind=kind_phys), allocatable :: EDDYVESTX ( : )
! in canopy layer [m]
real(kind=kind_phys), allocatable :: ZCANX ( : )
! Declare local maximum canopy layers
integer, parameter :: MAXCAN = 1000
integer, parameter :: mvt = 30 ! use 30 instead of 27
!Based on MODIS IGBP 20 Category Dataset
real :: fch_table(mvt) !< top of canopy (m)
data ( fch_table(i),i=1,mvt) /
& 20.0, 20.0, 18.0, 16.0, 16.0, 1.10,
& 1.10, 13.0, 10.0, 1.00, 5.00, 2.00,
& 15.0, 1.50, 0.00, 0.00, 0.00, 4.00,
& 2.00, 0.50, 0.00, 0.00, 0.00, 0.00,
& 0.00, 0.00, 0.00, 0.00, 0.00, 0.00 /
!----------------------------------------------

!IVAI
! integer, intent(in) :: naux2d,naux3d
! real(kind_phys), intent(inout) :: aux2d(:,:)
! real(kind_phys), intent(inout) :: aux3d(:,:,:)
!IVAI

!!
parameter(bfac=100.)
parameter(wfac=7.0,cfac=4.5)
Expand All @@ -282,6 +329,14 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
parameter(ck1=0.15,ch1=0.15)
parameter(cs0=0.4,csmf=0.5)
parameter(rchck=1.5,ndt=20)
!PCC_CANOPY------------------------------------
if (do_canopy) then
if(.not.allocated(EDDYVESTX))
& allocate( EDDYVESTX ( MAXCAN ) )
if(.not.allocated(ZCANX))
& allocate( ZCANX ( MAXCAN ) )
endif
!----------------------------------------------

if (tc_pbl == 0) then
ck0 = 0.4
Expand Down Expand Up @@ -1293,6 +1348,179 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
!
enddo
enddo
!PCC_CANOPY------------------------------------
kount=0 !IVAI
if (do_canopy) then

!IVAI
! print*, 'SATMEDMFVDIFQ_RUN: CLAIE = ', claie(:)
! print*, 'SATMEDMFVDIFQ_RUN: CFCH = ' , cfch (:)
! print*, 'SATMEDMFVDIFQ_RUN: CFRT = ' , cfrt (:)
! print*, 'SATMEDMFVDIFQ_RUN: CCLU = ' , cclu (:)
! print*, 'SATMEDMFVDIFQ_RUN: CPOPU= ' , cpopu(:)
! 2D aux arrays: canopy data in diffusion
! aux2d(:,1) = cfch (:)
! aux2d(:,2) = claie(:)
! aux2d(:,3) = cfrt(:)

! 3D aux arrays: before canopy correction
! aux3d(:,:,1) = dkq(:,:)
! aux3d(:,:,2) = dkt(:,:)
! aux3d(:,:,3) = dku(:,:)
!IVAI
do k = 1, km1-1
do i = 1, im

!IVAI: AQM canopy Inputs
! FCH = fch_table(vegtype(i)) !top of canopy from look-up table
FCH = cfch(i) !top of canopy from AQM canopy inputs
IF (k .EQ. 1) THEN !use model layer interfaces
KCAN = 1
ELSE
IF ( cfch(i) .GT. zi(i,k)
& .AND. cfch(i) .LE. zi(i,k+1) ) THEN
KCAN = 1
ELSE
KCAN = 0
END IF
END IF

IF (KCAN .EQ. 1) THEN !canopy inside model layer
! Check for other Contiguous Canopy Grid Cell Conditions

! Not a contigous canopy cell
IF ( claie(i) .LT. 0.1
& .OR. cfch (i) .LT. 0.5
!IVAI: modified contiguous canopy condition
! & .OR. MAX(0.0, 1.0 - cfrt(i)) .GT. 0.5
& .OR. MAX(0.0, 1.0 - cfrt(i)) .GT. 0.75
!IVAI
& .OR. cpopu(i) .GT. 10000.0
& .OR. (EXP(-0.5*claie(i)*cclu(i)) .GT. 0.45
& .AND. cfch(i) .LT. 18.) ) THEN


!TODO: Canopy Inputs
! IF ( XCANOPYLAI .LT. 0.1 !from canopy inputs
! IF ( lai(i) .LT. 0.1 !from LSM
! & .OR. FCH .LT. 0.5 ) THEN
! & .OR. MAX(0.0, 1.0 - XCANOPYFRT) .GT. 0.5
! & .OR. POPU .GT. 10000.0
! & .OR. EXP(-0.5*XCANOPYLAI*XCANOPYCLU).GT. 0.45
! & .AND. FCH .LT. 18.0 ) THEN

dkt(i,k)= dkt(i,k)
dkq(i,k)= dkq(i,k)
dku(i,k)= dku(i,k)

ELSE ! There is a contiguous forest canopy, apply correction over canopy layers

! Output contiguous canopy mask
! if (kount .EQ. 0 ) aux2d(i,5) = aux2d(i,5) + 1

!Raupauch M. R. A Practical Lagrangian method for relating scalar
!concentrations to
! source distributions in vegetation canopies. Q. J. R. Meteor. Soc.
! (1989), 115, pp 609-632
MOL = zol(i)/zl(i,k) !Monin-Obukhov Length in layer
HOL = FCH/MOL !local canopy stability parameter (hc/MOL)
ZCAN = zi(i,k+1) ! Initialize each model layer top that contains canopy (m)
! Integrate across total model interface
ZFL = ZCAN ! Set ZFL = ZCAN
COUNTCAN = 0 ! Initialize canopy layers

IF (k .EQ. 1) THEN !Find bottom in each model layer
BOTCAN = 0.5
ELSE
BOTCAN = zi(i,k)
END IF

DO WHILE (ZCAN.GE.BOTCAN)
! TLCAN = Lagrangian timescale
TLCAN = (FCH/ustar(i)) * (
& (0.256 * (ZCAN-(0.75*FCH))/FCH ) +
& (0.492*EXP((-0.256*ZCAN/FCH)/0.492)) )
IF ( HOL .LT. -0.1 ) THEN !STRONG UNSTABLE
IF ( ZCAN/FCH .GT. 1.25 ) THEN !SIGMACAN = Eulerian vertical velocity variance
SIGMACAN = 1.25*ustar(i)
END IF
IF ( ZCAN/FCH .GE. 0.175
& .AND. ZCAN/FCH .LE. 1.25 ) THEN
SIGMACAN = ustar(i) * ( 0.75 +
& (0.5 * COS((PI/1.06818) *
& (1.25 - (ZCAN/FCH)))) )
END IF
IF ( ZCAN/FCH .LT. 0.175 ) THEN
SIGMACAN = 0.25*ustar(i)
END IF
END IF
IF ( HOL .GE. -0.1 .AND. HOL .LT. 0.1 ) THEN !WEAKLY UNSTABLE to NEUTRAL
IF ( ZCAN/FCH .GT. 1.25 ) THEN
SIGMACAN = 1.0*ustar(i)
END IF
IF ( ZCAN/FCH .GE. 0.175
& .AND. ZCAN/FCH .LE. 1.25 ) THEN
SIGMACAN = ustar(i) * ( 0.625 +
& (0.375* COS((PI/1.06818) *
& (1.25 - (ZCAN/FCH)))) )
END IF
IF ( ZCAN/FCH .LT. 0.175 ) THEN
SIGMACAN = 0.25*ustar(i)
END IF
END IF
IF ( HOL .GE. 0.1 .AND. HOL .LT. 0.9 ) THEN !STABLE
IF ( ZCAN/FCH .GT. 1.25 ) THEN
SIGMACAN = 0.25*(4.375 - (3.75*HOL))*ustar(i)
END IF
IF ( ZCAN/FCH .GE. 0.175
& .AND. ZCAN/FCH .LE. 1.25 ) THEN
RRCAN=4.375-(3.75*HOL)
AACAN=(0.125*RRCAN) + 0.125
BBCAN=(0.125*RRCAN) - 0.125
SIGMACAN = ustar(i) * ( AACAN +
& (BBCAN * COS((PI/1.06818) *
& (1.25 - (ZCAN/FCH)))) )
END IF
IF ( ZCAN/FCH .LT. 0.175 ) THEN
SIGMACAN = 0.25*ustar(i)
END IF
END IF
IF ( HOL .GE. 0.9 ) THEN !VERY STABLE
SIGMACAN = 0.25*ustar(i)
END IF
IF ( ZCAN .EQ. ZFL ) THEN ! Each model layer that includes canopy
EDDYVEST1 = (SIGMACAN*SIGMACAN)*TLCAN
ELSE IF ( ZCAN .LE. FCH ) THEN !in-canopy layers and set arrays
COUNTCAN = COUNTCAN + 1
ZCANX (COUNTCAN) = ZCAN
EDDYVESTX (COUNTCAN) = (SIGMACAN*SIGMACAN)*TLCAN
END IF
ZCAN = ZCAN-0.5 !step down in-canopy resolution of 0.5m
END DO !end loop on canopy layers
EDDYVEST_INT = IntegrateTrapezoid((ZCANX(COUNTCAN:1:-1)
& ),EDDYVESTX(COUNTCAN:1:-1)) / ZFL
dkt(i,k)= (dkt(i,k)/EDDYVEST1) * EDDYVEST_INT !Scale dkt to resolved eddy diffusivity
dkq(i,k)= (dkq(i,k)/EDDYVEST1) * EDDYVEST_INT !Scale dkq to resolved eddy diffusivity
dku(i,k)= (dku(i,k)/EDDYVEST1) * EDDYVEST_INT !Scale dku to resolved eddy diffusivity

!IVAI: Output contiguos canopy correction bottom layer and 3D
! if ( kount .EQ. 0)
! & aux2d(i,4) = 1./EDDYVEST1 * EDDYVEST_INT
! aux3d(i,k,4) = 1./EDDYVEST1 * EDDYVEST_INT
!IVAI

END IF ! contigous canopy conditions

END IF ! (KCAN .EQ. 1) model layer(s) containing canopy

enddo !i

kount = kount + 1 !IVAI

enddo !k

endif !do_canopy

!> ## Compute TKE.
!! - Compute a minimum TKE deduced from background diffusivity for momentum.
!
Expand Down Expand Up @@ -1564,6 +1792,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
tkeh(i,k) = 0.5 * (tke(i,k)+tke(i,k+1))
enddo
enddo

do k=1,kps
do i=1,im
e_diff(i,k) = tke(i,k) - tke(i,k+1)
Expand Down
Loading