From a65fbed2abb6c3590143972b074faf750cd16de7 Mon Sep 17 00:00:00 2001 From: Robert Hallberg Date: Wed, 24 Nov 2021 15:10:25 -0500 Subject: [PATCH] +Rescale fluxes passed to KPP_NonLocalTransport Applied appropriate rescaling for dimensional consistency testing to the net heat and salt fluxes calculated in calculateBuoyancyFlux and used in the two KPP_NonLocalTransport_...() routines. Several comments describing variables were also either corrected or added. All answers are bitwise identical. --- src/core/MOM_forcing_type.F90 | 52 ++++++++--------- .../vertical/MOM_CVMix_KPP.F90 | 57 +++++++++++-------- .../vertical/MOM_diabatic_driver.F90 | 28 ++++----- .../vertical/MOM_opacity.F90 | 7 +-- 4 files changed, 75 insertions(+), 69 deletions(-) diff --git a/src/core/MOM_forcing_type.F90 b/src/core/MOM_forcing_type.F90 index dced9537d9..9cec8587db 100644 --- a/src/core/MOM_forcing_type.F90 +++ b/src/core/MOM_forcing_type.F90 @@ -925,34 +925,32 @@ subroutine calculateBuoyancyFlux1d(G, GV, US, fluxes, optics, nsw, h, Temp, Salt real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: Salt !< salinity [ppt] type(thermo_var_ptrs), intent(inout) :: tv !< thermodynamics type integer, intent(in) :: j !< j-row to work on - real, dimension(SZI_(G),SZK_(GV)+1), intent(inout) :: buoyancyFlux !< buoyancy fluxes [L2 T-3 ~> m2 s-3] - real, dimension(SZI_(G)), intent(inout) :: netHeatMinusSW !< Surface heat flux excluding shortwave - !! [degC H s-1 ~> degC m s-1 or degC kg m-2 s-1] - real, dimension(SZI_(G)), intent(inout) :: netSalt !< surface salt flux - !! [ppt H s-1 ~> ppt m s-1 or ppt kg m-2 s-1] + real, dimension(SZI_(G),SZK_(GV)+1), intent(out) :: buoyancyFlux !< buoyancy fluxes [L2 T-3 ~> m2 s-3] + real, dimension(SZI_(G)), intent(out) :: netHeatMinusSW !< Surface heat flux excluding shortwave + !! [degC H T-1 ~> degC m s-1 or degC kg m-2 s-1] + real, dimension(SZI_(G)), intent(out) :: netSalt !< surface salt flux + !! [ppt H T-1 ~> ppt m s-1 or ppt kg m-2 s-1] ! local variables - integer :: k - real, parameter :: dt = 1. ! to return a rate from extractFluxes1d - real, dimension(SZI_(G)) :: netH ! net FW flux [H s-1 ~> m s-1 or kg m-2 s-1] + real, dimension(SZI_(G)) :: netH ! net FW flux [H T-1 ~> m s-1 or kg m-2 s-1] real, dimension(SZI_(G)) :: netEvap ! net FW flux leaving ocean via evaporation - ! [H s-1 ~> m s-1 or kg m-2 s-1] - real, dimension(SZI_(G)) :: netHeat ! net temp flux [degC H s-1 ~> degC m s-1 or degC kg m-2 s-1] + ! [H T-1 ~> m s-1 or kg m-2 s-1] + real, dimension(SZI_(G)) :: netHeat ! net temp flux [degC H T-1 ~> degC m s-1 or degC kg m-2 s-1] real, dimension(max(nsw,1), SZI_(G)) :: penSWbnd ! penetrating SW radiation by band - ! [degC H s-1 ~> degC m s-1 or degC kg m-2 s-1] + ! [degC H T-1 ~> degC m s-1 or degC kg m-2 s-1] real, dimension(SZI_(G)) :: pressure ! pressure at the surface [R L2 T-2 ~> Pa] real, dimension(SZI_(G)) :: dRhodT ! density partial derivative wrt temp [R degC-1 ~> kg m-3 degC-1] real, dimension(SZI_(G)) :: dRhodS ! density partial derivative wrt saln [R ppt-1 ~> kg m-3 ppt-1] real, dimension(SZI_(G),SZK_(GV)+1) :: netPen ! The net penetrating shortwave radiation at each level - ! [degC H s-1 ~> degC m s-1 or degC kg m-2 s-1] + ! [degC H T-1 ~> degC m s-1 or degC kg m-2 s-1] logical :: useRiverHeatContent logical :: useCalvingHeatContent real :: depthBeforeScalingFluxes ! A depth scale [H ~> m or kg m-2] - real :: GoRho ! The gravitational acceleration divided by mean density times some - ! unit conversion factors [L2 H-1 s R-1 T-3 ~> m4 kg-1 s-2 or m7 kg-2 s-2] + real :: GoRho ! The gravitational acceleration divided by mean density times a + ! unit conversion factor [L2 H-1 R-1 T-2 ~> m4 kg-1 s-2 or m7 kg-2 s-2] real :: H_limit_fluxes ! Another depth scale [H ~> m or kg m-2] - integer :: i + integer :: i, k ! smg: what do we do when have heat fluxes from calving and river? useRiverHeatContent = .False. @@ -961,25 +959,25 @@ subroutine calculateBuoyancyFlux1d(G, GV, US, fluxes, optics, nsw, h, Temp, Salt depthBeforeScalingFluxes = max( GV%Angstrom_H, 1.e-30*GV%m_to_H ) pressure(:) = 0. if (associated(tv%p_surf)) then ; do i=G%isc,G%iec ; pressure(i) = tv%p_surf(i,j) ; enddo ; endif - GoRho = (GV%g_Earth * GV%H_to_Z*US%T_to_s) / GV%Rho0 + GoRho = (GV%g_Earth * GV%H_to_Z) / GV%Rho0 H_limit_fluxes = depthBeforeScalingFluxes ! The surface forcing is contained in the fluxes type. ! We aggregate the thermodynamic forcing for a time step into the following: - ! netH = water added/removed via surface fluxes [H s-1 ~> m s-1 or kg m-2 s-1] - ! netHeat = heat via surface fluxes [degC H s-1 ~> degC m s-1 or degC kg m-2 s-1] - ! netSalt = salt via surface fluxes [ppt H s-1 ~> ppt m s-1 or gSalt m-2 s-1] + ! netH = water added/removed via surface fluxes [H T-1 ~> m s-1 or kg m-2 s-1] + ! netHeat = heat via surface fluxes [degC H T-1 ~> degC m s-1 or degC kg m-2 s-1] + ! netSalt = salt via surface fluxes [ppt H T-1 ~> ppt m s-1 or gSalt m-2 s-1] ! Note that unlike other calls to extractFLuxes1d() that return the time-integrated flux - ! this call returns the rate because dt=1 - call extractFluxes1d(G, GV, US, fluxes, optics, nsw, j, dt*US%s_to_T, & + ! this call returns the rate because dt=1 (in arbitrary time units) + call extractFluxes1d(G, GV, US, fluxes, optics, nsw, j, 1.0, & depthBeforeScalingFluxes, useRiverHeatContent, useCalvingHeatContent, & h(:,j,:), Temp(:,j,:), netH, netEvap, netHeatMinusSW, & netSalt, penSWbnd, tv, .false.) ! Sum over bands and attenuate as a function of depth ! netPen is the netSW as a function of depth - call sumSWoverBands(G, GV, US, h(:,j,:), optics_nbands(optics), optics, j, dt*US%s_to_T, & + call sumSWoverBands(G, GV, US, h(:,j,:), optics_nbands(optics), optics, j, 1.0, & H_limit_fluxes, .true., penSWbnd, netPen) ! Density derivatives @@ -987,12 +985,14 @@ subroutine calculateBuoyancyFlux1d(G, GV, US, fluxes, optics, nsw, h, Temp, Salt tv%eqn_of_state, EOS_domain(G%HI)) ! Adjust netSalt to reflect dilution effect of FW flux - netSalt(G%isc:G%iec) = netSalt(G%isc:G%iec) - Salt(G%isc:G%iec,j,1) * netH(G%isc:G%iec) ! ppt H/s + ! [ppt H T-1 ~> ppt m s-1 or ppt kg m-2 s-1] + netSalt(G%isc:G%iec) = netSalt(G%isc:G%iec) - Salt(G%isc:G%iec,j,1) * netH(G%isc:G%iec) ! Add in the SW heating for purposes of calculating the net ! surface buoyancy flux affecting the top layer. + ! [degC H T-1 ~> degC m s-1 or degC kg m-2 s-1] !netHeat(:) = netHeatMinusSW(:) + sum( penSWbnd, dim=1 ) - netHeat(G%isc:G%iec) = netHeatMinusSW(G%isc:G%iec) + netPen(G%isc:G%iec,1) ! K H/s + netHeat(G%isc:G%iec) = netHeatMinusSW(G%isc:G%iec) + netPen(G%isc:G%iec,1) ! Convert to a buoyancy flux, excluding penetrating SW heating buoyancyFlux(G%isc:G%iec,1) = - GoRho * ( dRhodS(G%isc:G%iec) * netSalt(G%isc:G%iec) + & @@ -1020,9 +1020,9 @@ subroutine calculateBuoyancyFlux2d(G, GV, US, fluxes, optics, h, Temp, Salt, tv, type(thermo_var_ptrs), intent(inout) :: tv !< thermodynamics type real, dimension(SZI_(G),SZJ_(G),SZK_(GV)+1), intent(inout) :: buoyancyFlux !< buoyancy fluxes [L2 T-3 ~> m2 s-3] real, dimension(SZI_(G),SZJ_(G)), intent(inout) :: netHeatMinusSW !< surface heat flux excluding shortwave - !! [degC H s-1 ~> degC m s-1 or degC kg m-2 s-1] + !! [degC H T-1 ~> degC m s-1 or degC kg m-2 s-1] real, dimension(SZI_(G),SZJ_(G)), intent(inout) :: netSalt !< Net surface salt flux - !! [ppt H s-1 ~> ppt m s-1 or ppt kg m-2 s-1] + !! [ppt H T-1 ~> ppt m s-1 or ppt kg m-2 s-1] ! local variables integer :: j diff --git a/src/parameterizations/vertical/MOM_CVMix_KPP.F90 b/src/parameterizations/vertical/MOM_CVMix_KPP.F90 index 5342c621dd..2ff0a21196 100644 --- a/src/parameterizations/vertical/MOM_CVMix_KPP.F90 +++ b/src/parameterizations/vertical/MOM_CVMix_KPP.F90 @@ -538,9 +538,9 @@ logical function KPP_init(paramFile, G, GV, US, diag, Time, CS, passive) CS%id_buoyFlux = register_diag_field('ocean_model', 'KPP_buoyFlux', diag%axesTi, Time, & 'Surface (and penetrating) buoyancy flux, as used by [CVMix] KPP', 'm2/s3', conversion=US%L_to_m**2*US%s_to_T**3) CS%id_QminusSW = register_diag_field('ocean_model', 'KPP_QminusSW', diag%axesT1, Time, & - 'Net temperature flux ignoring short-wave, as used by [CVMix] KPP', 'K m/s') + 'Net temperature flux ignoring short-wave, as used by [CVMix] KPP', 'K m/s', conversion=GV%H_to_m*US%s_to_T) CS%id_netS = register_diag_field('ocean_model', 'KPP_netSalt', diag%axesT1, Time, & - 'Effective net surface salt flux, as used by [CVMix] KPP', 'ppt m/s') + 'Effective net surface salt flux, as used by [CVMix] KPP', 'ppt m/s', conversion=GV%H_to_m*US%s_to_T) CS%id_Kt_KPP = register_diag_field('ocean_model', 'KPP_Kheat', diag%axesTi, Time, & 'Heat diffusivity due to KPP, as calculated by [CVMix] KPP', 'm2/s') CS%id_Kd_in = register_diag_field('ocean_model', 'KPP_Kd_in', diag%axesTi, Time, & @@ -554,13 +554,17 @@ logical function KPP_init(paramFile, G, GV, US, diag, Time, CS, passive) CS%id_NLTs = register_diag_field('ocean_model', 'KPP_NLtransport_salt', diag%axesTi, Time, & 'Non-local tranpsort (Cs*G(sigma)) for scalars, as calculated by [CVMix] KPP', 'nondim') CS%id_NLT_dTdt = register_diag_field('ocean_model', 'KPP_NLT_dTdt', diag%axesTL, Time, & - 'Temperature tendency due to non-local transport of heat, as calculated by [CVMix] KPP', 'K/s') + 'Temperature tendency due to non-local transport of heat, as calculated by [CVMix] KPP', & + 'K/s', conversion=US%s_to_T) CS%id_NLT_dSdt = register_diag_field('ocean_model', 'KPP_NLT_dSdt', diag%axesTL, Time, & - 'Salinity tendency due to non-local transport of salt, as calculated by [CVMix] KPP', 'ppt/s') + 'Salinity tendency due to non-local transport of salt, as calculated by [CVMix] KPP', & + 'ppt/s', conversion=US%s_to_T) CS%id_NLT_temp_budget = register_diag_field('ocean_model', 'KPP_NLT_temp_budget', diag%axesTL, Time, & - 'Heat content change due to non-local transport, as calculated by [CVMix] KPP', 'W/m^2') + 'Heat content change due to non-local transport, as calculated by [CVMix] KPP', & + 'W/m^2', conversion=US%QRZ_T_to_W_m2) CS%id_NLT_saln_budget = register_diag_field('ocean_model', 'KPP_NLT_saln_budget', diag%axesTL, Time, & - 'Salt content change due to non-local transport, as calculated by [CVMix] KPP', 'kg/(sec*m^2)') + 'Salt content change due to non-local transport, as calculated by [CVMix] KPP', & + 'kg/(sec*m^2)', conversion=US%RZ_T_to_kg_m2s) CS%id_Tsurf = register_diag_field('ocean_model', 'KPP_Tsurf', diag%axesT1, Time, & 'Temperature of surface layer (10% of OBL depth) as passed to [CVMix] KPP', 'C') CS%id_Ssurf = register_diag_field('ocean_model', 'KPP_Ssurf', diag%axesT1, Time, & @@ -1179,7 +1183,7 @@ subroutine KPP_compute_BLD(CS, G, GV, US, h, Temp, Salt, u, v, tv, uStar, buoyFl ! Calculate Bulk Richardson number from eq (21) of LMD94 BulkRi_1d = CVmix_kpp_compute_bulk_Richardson( & zt_cntr = cellHeight(1:GV%ke), & ! Depth of cell center [m] - delta_buoy_cntr=GoRho*deltaRho, & ! Bulk buoyancy difference, Br-B(z) [s-1] + delta_buoy_cntr=GoRho*deltaRho, & ! Bulk buoyancy difference, Br-B(z) [m s-2] delta_Vsqr_cntr=deltaU2, & ! Square of resolved velocity difference [m2 s-2] ws_cntr=Ws_1d, & ! Turbulent velocity scale profile [m s-1] N_iface=CS%N(i,j,:), & ! Buoyancy frequency [s-1] @@ -1285,12 +1289,12 @@ subroutine KPP_smooth_BLD(CS,G,GV,h) real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h !< Layer/level thicknesses [H ~> m or kg m-2] ! local - real, dimension(SZI_(G),SZJ_(G)) :: OBLdepth_prev ! OBLdepth before s.th smoothing iteration + real, dimension(SZI_(G),SZJ_(G)) :: OBLdepth_prev ! OBLdepth before s.th smoothing iteration [m] real, dimension( GV%ke ) :: cellHeight ! Cell center heights referenced to surface [m] ! (negative in the ocean) real, dimension( GV%ke+1 ) :: iFaceHeight ! Interface heights referenced to surface [m] ! (negative in the ocean) - real :: wc, ww, we, wn, ws ! averaging weights for smoothing + real :: wc, ww, we, wn, ws ! averaging weights for smoothing [nondim] real :: dh ! The local thickness used for calculating interface positions [m] real :: hcorr ! A cumulative correction arising from inflation of vanished layers [m] integer :: i, j, k, s @@ -1390,13 +1394,14 @@ subroutine KPP_NonLocalTransport_temp(CS, G, GV, h, nonLocalTrans, surfFlux, & real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h !< Layer/level thickness [H ~> m or kg m-2] real, dimension(SZI_(G),SZJ_(G),SZK_(GV)+1), intent(in) :: nonLocalTrans !< Non-local transport [nondim] real, dimension(SZI_(G),SZJ_(G)), intent(in) :: surfFlux !< Surface flux of temperature - !! [degC H s-1 ~> degC m s-1 or degC kg m-2 s-1] - real, intent(in) :: dt !< Time-step [s] - real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(inout) :: scalar !< temperature - real, intent(in) :: C_p !< Seawater specific heat capacity [J kg-1 degC-1] + !! [degC H T-1 ~> degC m s-1 or degC kg m-2 s-1] + real, intent(in) :: dt !< Time-step [T ~> s] + real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(inout) :: scalar !< temperature [degC] + real, intent(in) :: C_p !< Seawater specific heat capacity + !! [Q degC-1 ~> J kg-1 degC-1] integer :: i, j, k - real, dimension( SZI_(G), SZJ_(G),SZK_(GV) ) :: dtracer + real, dimension( SZI_(G), SZJ_(G),SZK_(GV) ) :: dtracer ! Rate of tracer change [degC T-1 ~> degC s-1] dtracer(:,:,:) = 0.0 @@ -1431,8 +1436,9 @@ subroutine KPP_NonLocalTransport_temp(CS, G, GV, h, nonLocalTrans, surfFlux, & do k = 1, GV%ke do j = G%jsc, G%jec do i = G%isc, G%iec + ! Here dtracer has units of [Q R Z T-1 ~> W m-2]. dtracer(i,j,k) = (nonLocalTrans(i,j,k) - nonLocalTrans(i,j,k+1)) * & - surfFlux(i,j) * C_p * GV%H_to_kg_m2 + surfFlux(i,j) * C_p * GV%H_to_RZ enddo enddo enddo @@ -1446,18 +1452,18 @@ end subroutine KPP_NonLocalTransport_temp !> This routine is a useful prototype for other material tracers. subroutine KPP_NonLocalTransport_saln(CS, G, GV, h, nonLocalTrans, surfFlux, dt, scalar) - type(KPP_CS), intent(in) :: CS !< Control structure - type(ocean_grid_type), intent(in) :: G !< Ocean grid - type(verticalGrid_type), intent(in) :: GV !< Ocean vertical grid - real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h !< Layer/level thickness [H ~> m or kg m-2] + type(KPP_CS), intent(in) :: CS !< Control structure + type(ocean_grid_type), intent(in) :: G !< Ocean grid + type(verticalGrid_type), intent(in) :: GV !< Ocean vertical grid + real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h !< Layer/level thickness [H ~> m or kg m-2] real, dimension(SZI_(G),SZJ_(G),SZK_(GV)+1), intent(in) :: nonLocalTrans !< Non-local transport [nondim] - real, dimension(SZI_(G),SZJ_(G)), intent(in) :: surfFlux !< Surface flux of salt - !! [ppt H s-1 ~> ppt m s-1 or ppt kg m-2 s-1] - real, intent(in) :: dt !< Time-step [s] - real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(inout) :: scalar !< Scalar (scalar units [conc]) + real, dimension(SZI_(G),SZJ_(G)), intent(in) :: surfFlux !< Surface flux of salt + !! [ppt H T-1 ~> ppt m s-1 or ppt kg m-2 s-1] + real, intent(in) :: dt !< Time-step [T ~> s] + real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(inout) :: scalar !< Salinity [ppt] integer :: i, j, k - real, dimension( SZI_(G), SZJ_(G),SZK_(GV) ) :: dtracer + real, dimension( SZI_(G), SZJ_(G),SZK_(GV) ) :: dtracer ! Rate of tracer change [ppt T-1 ~> ppt s-1] dtracer(:,:,:) = 0.0 @@ -1492,8 +1498,9 @@ subroutine KPP_NonLocalTransport_saln(CS, G, GV, h, nonLocalTrans, surfFlux, dt, do k = 1, GV%ke do j = G%jsc, G%jec do i = G%isc, G%iec + ! Here dtracer has units of [ppt R Z T-1 ~> ppt kg m-2 s-1] dtracer(i,j,k) = (nonLocalTrans(i,j,k) - nonLocalTrans(i,j,k+1)) * & - surfFlux(i,j) * GV%H_to_kg_m2 + surfFlux(i,j) * GV%H_to_RZ enddo enddo enddo diff --git a/src/parameterizations/vertical/MOM_diabatic_driver.F90 b/src/parameterizations/vertical/MOM_diabatic_driver.F90 index bac311bd6d..c0ea3aff53 100644 --- a/src/parameterizations/vertical/MOM_diabatic_driver.F90 +++ b/src/parameterizations/vertical/MOM_diabatic_driver.F90 @@ -248,9 +248,9 @@ module MOM_diabatic_driver real, allocatable, dimension(:,:,:) :: KPP_NLTscalar !< KPP non-local transport for scalars [nondim] real, allocatable, dimension(:,:,:) :: KPP_buoy_flux !< KPP forcing buoyancy flux [L2 T-3 ~> m2 s-3] real, allocatable, dimension(:,:) :: KPP_temp_flux !< KPP effective temperature flux - !! [degC H s-1 ~> degC m s-1 or degC kg m-2 s-1] + !! [degC H T-1 ~> degC m s-1 or degC kg m-2 s-1] real, allocatable, dimension(:,:) :: KPP_salt_flux !< KPP effective salt flux - !! [ppt H s-1 ~> ppt m s-1 or ppt kg m-2 s-1] + !! [ppt H T-1 ~> ppt m s-1 or ppt kg m-2 s-1] type(time_type), pointer :: Time !< Pointer to model time (needed for sponges) end type diabatic_CS @@ -680,17 +680,17 @@ subroutine diabatic_ALE_legacy(u, v, h, tv, Hml, fluxes, visc, ADp, CDp, dt, Tim call MOM_thermovar_chksum("after KPP", tv, G) call hchksum(Kd_heat, "after KPP Kd_heat", G%HI, haloshift=0, scale=US%Z2_T_to_m2_s) call hchksum(Kd_salt, "after KPP Kd_salt", G%HI, haloshift=0, scale=US%Z2_T_to_m2_s) - call hchksum(CS%KPP_temp_flux, "before KPP_applyNLT netHeat", G%HI, haloshift=0, scale=GV%H_to_m) - call hchksum(CS%KPP_salt_flux, "before KPP_applyNLT netSalt", G%HI, haloshift=0, scale=GV%H_to_m) + call hchksum(CS%KPP_temp_flux, "before KPP_applyNLT netHeat", G%HI, haloshift=0, scale=GV%H_to_m*US%s_to_T) + call hchksum(CS%KPP_salt_flux, "before KPP_applyNLT netSalt", G%HI, haloshift=0, scale=GV%H_to_m*US%s_to_T) call hchksum(CS%KPP_NLTheat, "before KPP_applyNLT NLTheat", G%HI, haloshift=0) call hchksum(CS%KPP_NLTscalar, "before KPP_applyNLT NLTscalar", G%HI, haloshift=0) endif ! Apply non-local transport of heat and salt ! Changes: tv%T, tv%S call KPP_NonLocalTransport_temp(CS%KPP_CSp, G, GV, h, CS%KPP_NLTheat, CS%KPP_temp_flux, & - US%T_to_s*dt, tv%T, US%Q_to_J_kg*tv%C_p) + dt, tv%T, tv%C_p) call KPP_NonLocalTransport_saln(CS%KPP_CSp, G, GV, h, CS%KPP_NLTscalar, CS%KPP_salt_flux, & - US%T_to_s*dt, tv%S) + dt, tv%S) call cpu_clock_end(id_clock_kpp) if (showCallTree) call callTree_waypoint("done with KPP_applyNonLocalTransport (diabatic)") if (CS%debugConservation) call MOM_state_stats('KPP_applyNonLocalTransport', u, v, h, tv%T, tv%S, G, GV, US) @@ -1255,17 +1255,17 @@ subroutine diabatic_ALE(u, v, h, tv, Hml, fluxes, visc, ADp, CDp, dt, Time_end, call MOM_thermovar_chksum("after KPP", tv, G) call hchksum(Kd_heat, "after KPP Kd_heat", G%HI, haloshift=0, scale=US%Z2_T_to_m2_s) call hchksum(Kd_salt, "after KPP Kd_salt", G%HI, haloshift=0, scale=US%Z2_T_to_m2_s) - call hchksum(CS%KPP_temp_flux, "before KPP_applyNLT netHeat", G%HI, haloshift=0, scale=GV%H_to_m) - call hchksum(CS%KPP_salt_flux, "before KPP_applyNLT netSalt", G%HI, haloshift=0, scale=GV%H_to_m) + call hchksum(CS%KPP_temp_flux, "before KPP_applyNLT netHeat", G%HI, haloshift=0, scale=GV%H_to_m*US%s_to_T) + call hchksum(CS%KPP_salt_flux, "before KPP_applyNLT netSalt", G%HI, haloshift=0, scale=GV%H_to_m*US%s_to_T) call hchksum(CS%KPP_NLTheat, "before KPP_applyNLT NLTheat", G%HI, haloshift=0) call hchksum(CS%KPP_NLTscalar, "before KPP_applyNLT NLTscalar", G%HI, haloshift=0) endif ! Apply non-local transport of heat and salt ! Changes: tv%T, tv%S call KPP_NonLocalTransport_temp(CS%KPP_CSp, G, GV, h, CS%KPP_NLTheat, CS%KPP_temp_flux, & - US%T_to_s*dt, tv%T, US%Q_to_J_kg*tv%C_p) + dt, tv%T, tv%C_p) call KPP_NonLocalTransport_saln(CS%KPP_CSp, G, GV, h, CS%KPP_NLTscalar, CS%KPP_salt_flux, & - US%T_to_s*dt, tv%S) + dt, tv%S) call cpu_clock_end(id_clock_kpp) if (showCallTree) call callTree_waypoint("done with KPP_applyNonLocalTransport (diabatic)") if (CS%debugConservation) call MOM_state_stats('KPP_applyNonLocalTransport', u, v, h, tv%T, tv%S, G, GV, US) @@ -1877,17 +1877,17 @@ subroutine layered_diabatic(u, v, h, tv, Hml, fluxes, visc, ADp, CDp, dt, Time_e if (CS%useKPP) then call cpu_clock_begin(id_clock_kpp) if (CS%debug) then - call hchksum(CS%KPP_temp_flux, "before KPP_applyNLT netHeat", G%HI, haloshift=0, scale=GV%H_to_m) - call hchksum(CS%KPP_salt_flux, "before KPP_applyNLT netSalt", G%HI, haloshift=0, scale=GV%H_to_m) + call hchksum(CS%KPP_temp_flux, "before KPP_applyNLT netHeat", G%HI, haloshift=0, scale=GV%H_to_m*US%s_to_T) + call hchksum(CS%KPP_salt_flux, "before KPP_applyNLT netSalt", G%HI, haloshift=0, scale=GV%H_to_m*US%s_to_T) call hchksum(CS%KPP_NLTheat, "before KPP_applyNLT NLTheat", G%HI, haloshift=0) call hchksum(CS%KPP_NLTscalar, "before KPP_applyNLT NLTscalar", G%HI, haloshift=0) endif ! Apply non-local transport of heat and salt ! Changes: tv%T, tv%S call KPP_NonLocalTransport_temp(CS%KPP_CSp, G, GV, h, CS%KPP_NLTheat, CS%KPP_temp_flux, & - US%T_to_s*dt, tv%T, US%Q_to_J_kg*tv%C_p) + dt, tv%T, tv%C_p) call KPP_NonLocalTransport_saln(CS%KPP_CSp, G, GV, h, CS%KPP_NLTscalar, CS%KPP_salt_flux, & - US%T_to_s*dt, tv%S) + dt, tv%S) call cpu_clock_end(id_clock_kpp) if (showCallTree) call callTree_waypoint("done with KPP_applyNonLocalTransport (diabatic)") if (CS%debugConservation) call MOM_state_stats('KPP_applyNonLocalTransport', u, v, h, tv%T, tv%S, G, GV, US) diff --git a/src/parameterizations/vertical/MOM_opacity.F90 b/src/parameterizations/vertical/MOM_opacity.F90 index 5dec767f5b..02d49d024d 100644 --- a/src/parameterizations/vertical/MOM_opacity.F90 +++ b/src/parameterizations/vertical/MOM_opacity.F90 @@ -773,7 +773,7 @@ end subroutine absorbRemainingSW !> This subroutine calculates the total shortwave heat flux integrated over !! bands as a function of depth. This routine is only called for computing -!! buoyancy fluxes for use in KPP. This routine does not updat e the state. +!! buoyancy fluxes for use in KPP. This routine does not update the state. subroutine sumSWoverBands(G, GV, US, h, nsw, optics, j, dt, & H_limit_fluxes, absorbAllSW, iPen_SW_bnd, netPen) type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure. @@ -793,9 +793,8 @@ subroutine sumSWoverBands(G, GV, US, h, nsw, optics, j, dt, & logical, intent(in) :: absorbAllSW !< If true, ensure that all shortwave !! radiation is absorbed in the ocean water column. real, dimension(max(nsw,1),SZI_(G)), intent(in) :: iPen_SW_bnd !< The incident penetrating shortwave - !! heating in each band that hits the bottom and - !! will be redistributed through the water column - !! [degC H ~> degC m or degC kg m-2]; size nsw x SZI_(G). + !! in each band at the sea surface; size nsw x SZI_(G) + !! [degC H ~> degC m or degC kg m-2]. real, dimension(SZI_(G),SZK_(GV)+1), & intent(inout) :: netPen !< Net penetrating shortwave heat flux at each !! interface, summed across all bands