From fc05abe214722c798a4c59c9cd9d24f98debb2e4 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Tue, 29 Mar 2022 15:17:19 -0600 Subject: [PATCH] Modifications needed for when enthalpy is via CPL * 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. --- .../vertical/MOM_diabatic_aux.F90 | 50 +++++++++++-------- 1 file changed, 30 insertions(+), 20 deletions(-) diff --git a/src/parameterizations/vertical/MOM_diabatic_aux.F90 b/src/parameterizations/vertical/MOM_diabatic_aux.F90 index 13d25f06f5..ea7739450f 100644 --- a/src/parameterizations/vertical/MOM_diabatic_aux.F90 +++ b/src/parameterizations/vertical/MOM_diabatic_aux.F90 @@ -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] @@ -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 @@ -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, & @@ -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 @@ -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