Skip to content

Commit

Permalink
Modifications needed for when enthalpy is via CPL
Browse files Browse the repository at this point in the history 
* Introduces a new constant (EnthalpyConst = 1.0, by default) which
is set to 0.0 when fluxes%heat_content_evap is associated. This
constant is used in the expression that accounts for the temperature
 of the mass exchange (dTemp) to avoid double-couting for the enthalpy
terms when they are provided via coupler.

* Use heat_content_evap to determine if the diagostics heat_content_massin,
heat_content_massout, and TempxPmE should be calculated.
  • Loading branch information
@gustavo-marques @alperaltuntas
gustavo-marques authored and alperaltuntas committed Apr 21, 2022
1 parent 7924fba commit fc05abe
Showing 1 changed file with 30 additions and 20 deletions.
50 changes: 30 additions & 20 deletions src/parameterizations/vertical/MOM_diabatic_aux.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1031,7 +1031,9 @@ subroutine applyBoundaryFluxesInOut(CS, G, GV, US, dt, fluxes, optics, nsw, h, t
real :: dThickness, dTemp, dSalt
real :: fractionOfForcing, hOld, Ithickness
real :: RivermixConst ! A constant used in implementing river mixing [R Z2 T-1 ~> Pa s].

real :: EnthalpyConst ! A constant used to control the enthalpy calculation
! By default EnthalpyConst = 1.0. If fluxes%heat_content_evap
! is associated enthalpy is provided via coupler and EnthalpyConst = 0.0.
real, dimension(SZI_(G)) :: &
d_pres, & ! pressure change across a layer [R L2 T-2 ~> Pa]
p_lay, & ! average pressure in a layer [R L2 T-2 ~> Pa]
Expand Down Expand Up @@ -1095,6 +1097,9 @@ subroutine applyBoundaryFluxesInOut(CS, G, GV, US, dt, fluxes, optics, nsw, h, t
! Only apply forcing if fluxes%sw is associated.
if (.not.associated(fluxes%sw) .and. .not.calculate_energetics) return

EnthalpyConst = 1.0
if (associated(fluxes%heat_content_evap)) EnthalpyConst = 0.0

if (calculate_buoyancy) then
SurfPressure(:) = 0.0
GoRho = US%L_to_Z**2*GV%g_Earth / GV%Rho0
Expand All @@ -1116,7 +1121,8 @@ subroutine applyBoundaryFluxesInOut(CS, G, GV, US, dt, fluxes, optics, nsw, h, t
!$OMP nonPenSW,hGrounding,CS,Idt,aggregate_FW_forcing, &
!$OMP minimum_forcing_depth,evap_CFL_limit,dt,EOSdom, &
!$OMP calculate_buoyancy,netPen_rate,SkinBuoyFlux,GoRho, &
!$OMP calculate_energetics,dSV_dT,dSV_dS,cTKE,g_Hconv2) &
!$OMP calculate_energetics,dSV_dT,dSV_dS,cTKE,g_Hconv2, &
!$OMP EnthalpyConst) &
!$OMP private(opacityBand,h2d,T2d,netMassInOut,netMassOut, &
!$OMP netHeat,netSalt,Pen_SW_bnd,fractionOfForcing, &
!$OMP IforcingDepthScale, &
Expand Down Expand Up @@ -1258,17 +1264,19 @@ subroutine applyBoundaryFluxesInOut(CS, G, GV, US, dt, fluxes, optics, nsw, h, t
! This line accounts for the temperature of the mass exchange
Temp_in = T2d(i,k)
Salin_in = 0.0
dTemp = dTemp + dThickness*Temp_in
dTemp = dTemp + dThickness*Temp_in*EnthalpyConst

! Diagnostics of heat content associated with mass fluxes
if (associated(fluxes%heat_content_massin)) &
fluxes%heat_content_massin(i,j) = fluxes%heat_content_massin(i,j) + &
T2d(i,k) * max(0.,dThickness) * GV%H_to_RZ * fluxes%C_p * Idt
if (associated(fluxes%heat_content_massout)) &
fluxes%heat_content_massout(i,j) = fluxes%heat_content_massout(i,j) + &
T2d(i,k) * min(0.,dThickness) * GV%H_to_RZ * fluxes%C_p * Idt
if (associated(tv%TempxPmE)) tv%TempxPmE(i,j) = tv%TempxPmE(i,j) + &
T2d(i,k) * dThickness * GV%H_to_RZ
if (.not. associated(fluxes%heat_content_evap)) then
if (associated(fluxes%heat_content_massin)) &
fluxes%heat_content_massin(i,j) = fluxes%heat_content_massin(i,j) + &
T2d(i,k) * max(0.,dThickness) * GV%H_to_RZ * fluxes%C_p * Idt
if (associated(fluxes%heat_content_massout)) &
fluxes%heat_content_massout(i,j) = fluxes%heat_content_massout(i,j) + &
T2d(i,k) * min(0.,dThickness) * GV%H_to_RZ * fluxes%C_p * Idt
if (associated(tv%TempxPmE)) tv%TempxPmE(i,j) = tv%TempxPmE(i,j) + &
T2d(i,k) * dThickness * GV%H_to_RZ
endif

! Determine the energetics of river mixing before updating the state.
if (calculate_energetics .and. associated(fluxes%lrunoff) .and. CS%do_rivermix) then
Expand Down Expand Up @@ -1341,17 +1349,19 @@ subroutine applyBoundaryFluxesInOut(CS, G, GV, US, dt, fluxes, optics, nsw, h, t
netSalt(i) = netSalt(i) - dSalt

! This line accounts for the temperature of the mass exchange
dTemp = dTemp + dThickness*T2d(i,k)
dTemp = dTemp + dThickness*T2d(i,k)*EnthalpyConst

! Diagnostics of heat content associated with mass fluxes
if (associated(fluxes%heat_content_massin)) &
fluxes%heat_content_massin(i,j) = fluxes%heat_content_massin(i,j) + &
T2d(i,k) * max(0.,dThickness) * GV%H_to_RZ * fluxes%C_p * Idt
if (associated(fluxes%heat_content_massout)) &
fluxes%heat_content_massout(i,j) = fluxes%heat_content_massout(i,j) + &
T2d(i,k) * min(0.,dThickness) * GV%H_to_RZ * fluxes%C_p * Idt
if (associated(tv%TempxPmE)) tv%TempxPmE(i,j) = tv%TempxPmE(i,j) + &
T2d(i,k) * dThickness * GV%H_to_RZ
if (.not. associated(fluxes%heat_content_evap)) then
if (associated(fluxes%heat_content_massin)) &
fluxes%heat_content_massin(i,j) = fluxes%heat_content_massin(i,j) + &
T2d(i,k) * max(0.,dThickness) * GV%H_to_RZ * fluxes%C_p * Idt
if (associated(fluxes%heat_content_massout)) &
fluxes%heat_content_massout(i,j) = fluxes%heat_content_massout(i,j) + &
T2d(i,k) * min(0.,dThickness) * GV%H_to_RZ * fluxes%C_p * Idt
if (associated(tv%TempxPmE)) tv%TempxPmE(i,j) = tv%TempxPmE(i,j) + &
T2d(i,k) * dThickness * GV%H_to_RZ
endif

! Update state by the appropriate increment.
hOld = h2d(i,k) ! Keep original thickness in hand
Expand Down

0 comments on commit fc05abe

Please sign in to comment.