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.

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

clm lake fix: pressure is not density + #142 and #148 #139

Merged
merged 15 commits into from
Dec 21, 2023
Merged
Show file tree
Hide file tree
Changes from 6 commits
Commits
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
23 changes: 14 additions & 9 deletions physics/clm_lake.f90
Original file line number Diff line number Diff line change
Expand Up @@ -238,7 +238,10 @@ subroutine calculate_z_dz_lake(i,input_lakedepth,clm_lakedepth,z_lake,dz_lake)
real(kind_lake) :: dz_lake(nlevlake) ! layer thickness for lake (m)
real(kind_lake) :: depthratio

if (input_lakedepth(i) == spval) then
if (input_lakedepth(i) == spval .or. input_lakedepth(i) < 0.1) then
SamuelTrahanNOAA marked this conversation as resolved.
Show resolved Hide resolved
! This is a safeguard against:
! 1. missing in the lakedepth database (== spval)
! 2. errors in model cycling or unexpected changes in the orography database (< 0.1)
clm_lakedepth(i) = zlak(nlevlake) + 0.5_kind_lake*dzlak(nlevlake)
z_lake(1:nlevlake) = zlak(1:nlevlake)
dz_lake(1:nlevlake) = dzlak(1:nlevlake)
Expand Down Expand Up @@ -267,7 +270,7 @@ SUBROUTINE clm_lake_run( &

! Atmospheric model state inputs:
tg3, pgr, zlvl, gt0, prsi, phii, qvcurr, gu0, gv0, xlat_d, xlon_d, &
ch, cm, dlwsfci, dswsfci, oro_lakedepth, wind, rho0, tsfc, &
ch, cm, dlwsfci, dswsfci, oro_lakedepth, wind, tsfc, &
flag_iter, ISLTYP, rainncprv, raincprv, &

! Feedback to atmosphere:
Expand All @@ -283,7 +286,7 @@ SUBROUTINE clm_lake_run( &

salty, savedtke12d, snowdp2d, h2osno2d, snl2d, t_grnd2d, t_lake3d, &
lake_icefrac3d, t_soisno3d, h2osoi_ice3d, h2osoi_liq3d, h2osoi_vol3d, &
z3d, dz3d, zi3d, &
z3d, dz3d, zi3d, t1, qv1, prsl1, &
input_lakedepth, clm_lakedepth, cannot_freeze, &

! Error reporting:
Expand Down Expand Up @@ -321,8 +324,8 @@ SUBROUTINE clm_lake_run( &
!
REAL(KIND_PHYS), DIMENSION(:), INTENT(IN):: &
tg3, pgr, zlvl, qvcurr, xlat_d, xlon_d, ch, cm, &
dlwsfci, dswsfci, oro_lakedepth, wind, rho0, &
rainncprv, raincprv
dlwsfci, dswsfci, oro_lakedepth, wind, &
rainncprv, raincprv, t1, qv1, prsl1
REAL(KIND_PHYS), DIMENSION(:,:), INTENT(in) :: gu0, gv0, prsi, gt0, phii
LOGICAL, DIMENSION(:), INTENT(IN) :: flag_iter
INTEGER, DIMENSION(:), INTENT(IN) :: ISLTYP
Expand Down Expand Up @@ -450,6 +453,7 @@ SUBROUTINE clm_lake_run( &
logical, parameter :: feedback_to_atmosphere = .true. ! FIXME: REMOVE

real(kind_lake) :: to_radians, lat_d, lon_d, qss, tkm, bd
real(kind_lake) :: rho0 ! lowest model level air density

integer :: month,num1,num2,day_of_month,isl
real(kind_lake) :: wght1,wght2,Tclim,depthratio
Expand Down Expand Up @@ -693,12 +697,13 @@ SUBROUTINE clm_lake_run( &

!-- The CLM output is combined for fractional ice and water
if( t_grnd(c) >= tfrz ) then
qfx = eflx_lh_tot(c)*invhvap
qfx = eflx_lh_tot(c)*invhvap
else
qfx = eflx_lh_tot(c)*invhsub ! heat flux (W/m^2)=>mass flux(kg/(sm^2))
qfx = eflx_lh_tot(c)*invhsub ! heat flux (W/m^2)=>mass flux(kg/(sm^2))
endif
evap_wat(i) = qfx/rho0(i) ! kinematic_surface_upward_latent_heat_flux_over_water
hflx_wat(i)=eflx_sh_tot(c)/(rho0(i)*cpair) ! kinematic_surface_upward_sensible_heat_flux_over_water
rho0 = prsl1(i) / (rair*t1(i)*(1.0 + con_fvirt*qv1(i)))
evap_wat(i) = qfx/rho0 ! kinematic_surface_upward_latent_heat_flux_over_water
hflx_wat(i) = eflx_sh_tot(c)/(rho0*cpair) ! kinematic_surface_upward_sensible_heat_flux_over_water
gflx_wat(I) = eflx_gnet(c) ![W/m/m] upward_heat_flux_in_soil_over_water
ep1d_water(i) = eflx_lh_tot(c) ![W/m/m] surface_upward_potential_latent_heat_flux_over_water
tsurf_water(I) = t_grnd(c) ![K] surface skin temperature after iteration over water
Expand Down
32 changes: 24 additions & 8 deletions physics/clm_lake.meta
Original file line number Diff line number Diff line change
Expand Up @@ -305,14 +305,6 @@
type = real
kind = kind_phys
intent = in
[rho0]
standard_name = air_pressure_at_surface_adjacent_layer
long_name = mean pressure at lowest model layer
units = Pa
dimensions = (horizontal_loop_extent)
type = real
kind = kind_phys
intent = in
[tsfc]
standard_name = surface_skin_temperature
long_name = surface skin temperature
Expand Down Expand Up @@ -732,6 +724,30 @@
type = real
kind = kind_phys
intent = in
[t1]
standard_name = air_temperature_at_surface_adjacent_layer
long_name = mean temperature at lowest model layer
units = K
dimensions = (horizontal_loop_extent)
type = real
kind = kind_phys
intent = in
[qv1]
standard_name = specific_humidity_at_surface_adjacent_layer
long_name = water vapor specific humidity at lowest model layer
units = kg kg-1
dimensions = (horizontal_loop_extent)
type = real
kind = kind_phys
intent = in
[prsl1]
standard_name = air_pressure_at_surface_adjacent_layer
long_name = mean pressure at lowest model layer
units = Pa
dimensions = (horizontal_loop_extent)
type = real
kind = kind_phys
intent = in
[errmsg]
standard_name = ccpp_error_message
long_name = error message for error handling in CCPP
Expand Down
22 changes: 6 additions & 16 deletions physics/module_sf_ruclsm.F90
Original file line number Diff line number Diff line change
Expand Up @@ -1703,7 +1703,8 @@ SUBROUTINE SFCTMP (debug_print, delt,ktau,conflx,i,j, & !--- input varia
IF (NEWSN > zero .and. snowfracnewsn > 0.99_kind_phys .and. rhosnfall < 450._kind_phys) THEN
! new snow
KEEP_SNOW_ALBEDO = one
!snow_mosaic=0. ! ???
! turn off separate treatment of snow covered and snow-free portions of the grid cell
snow_mosaic=0. ! ???
ENDIF

IF (debug_print ) THEN
Expand Down Expand Up @@ -2076,7 +2077,6 @@ SUBROUTINE SFCTMP (debug_print, delt,ktau,conflx,i,j, & !--- input varia
hfx = hfxs*(one-snowfrac) + hfx*snowfrac
s = ss*(one-snowfrac) + s*snowfrac
evapl = evapls*(one-snowfrac)
sublim = sublim*snowfrac
prcpl = prcpls*(one-snowfrac) + prcpl*snowfrac
fltot = fltots*(one-snowfrac) + fltot*snowfrac
ALB = MAX(keep_snow_albedo*alb, &
Expand All @@ -2088,10 +2088,6 @@ SUBROUTINE SFCTMP (debug_print, delt,ktau,conflx,i,j, & !--- input varia

runoff1 = runoff1s*(one-snowfrac) + runoff1*snowfrac
runoff2 = runoff2s*(one-snowfrac) + runoff2*snowfrac
smelt = smelt * snowfrac
snoh = snoh * snowfrac
snflx = snflx * snowfrac
snom = snom * snowfrac
mavail = mavails*(one-snowfrac) + one*snowfrac
infiltr = infiltrs*(one-snowfrac) + infiltr*snowfrac

Expand All @@ -2115,7 +2111,7 @@ SUBROUTINE SFCTMP (debug_print, delt,ktau,conflx,i,j, & !--- input varia
qvg = qvgs*(one-snowfrac) + qvg*snowfrac
qsg = qsgs*(one-snowfrac) + qsg*snowfrac
qcg = qcgs*(one-snowfrac) + qcg*snowfrac
sublim = eeta*snowfrac
sublim = eeta
eeta = eetas*(one-snowfrac) + eeta*snowfrac
qfx = qfxs*(one-snowfrac) + qfx*snowfrac
hfx = hfxs*(one-snowfrac) + hfx*snowfrac
Expand All @@ -2129,10 +2125,6 @@ SUBROUTINE SFCTMP (debug_print, delt,ktau,conflx,i,j, & !--- input varia
(emissn - emiss_snowfree) * snowfrac), emissn))
runoff1 = runoff1s*(one-snowfrac) + runoff1*snowfrac
runoff2 = runoff2s*(one-snowfrac) + runoff2*snowfrac
smelt = smelt * snowfrac
snoh = snoh * snowfrac
snflx = snflx * snowfrac
snom = snom * snowfrac
IF (debug_print ) THEN
print *,'SOILT combined on ice', soilt
ENDIF
Expand Down Expand Up @@ -2215,15 +2207,13 @@ SUBROUTINE SFCTMP (debug_print, delt,ktau,conflx,i,j, & !--- input varia
IF (debug_print ) then
!if (abs(xlat-testptlat).lt.0.2 .and. abs(xlon-testptlon).lt.0.2)then
print *,'Snowfallac xlat, xlon',xlat,xlon
print *,'newsn,rhonewsn,newsnowratio=',newsn,rhonewsn,newsnowratio
print *,'newsn [m],rhonewsn,newsnowratio=',newsn,rhonewsn,newsnowratio
print *,'Time-step newsn depth [m], swe [m]',newsn,newsn*rhonewsn
print *,'Time-step smelt: swe [m]' ,smelt*delt
print *,'Time-step sublim: swe,[kg m-2]',sublim*delt
endif

snowfallac = snowfallac + max(zero,(newsn*rhonewsn - & ! source of snow (swe) [m]
(smelt+sublim*1.e-3_kind_phys)*delt*newsnowratio) & ! sink: melting and sublimation, (swe) [m]
/rhonewsn)*rhowater ! snow accumulation in snow depth [mm]
snowfallac = snowfallac + newsn * 1.e3_kind_phys ! accumulated snow depth [mm], using variable snow density

IF (debug_print ) THEN
!if (abs(xlat-testptlat).lt.0.2 .and. abs(xlon-testptlon).lt.0.2)then
Expand Down Expand Up @@ -5596,7 +5586,7 @@ SUBROUTINE SNOWTEMP( debug_print,xlat,xlon, &
nmelt = 1
soiltfrac=snowfrac*tfrz+(one-snowfrac)*SOILT
QSG=min(QSG, QSN(soiltfrac,TBQ)/PP)
qvg=qsg
qvg=snowfrac*qsg+(one-snowfrac)*qvg
T3 = STBOLT*TN*TN*TN
UPFLUX = T3 * 0.5_kind_phys*(TN + SOILTfrac)
XINET = EMISS*(GLW-UPFLUX)
Expand Down
Loading