From b13d25ad5fae9968756b392d71cd961cc21dc517 Mon Sep 17 00:00:00 2001 From: Milan Date: Sat, 3 Jun 2023 09:41:06 -0400 Subject: [PATCH] precip units in [m] to [mm] for output --- src/SpeedyWeather.jl | 2 +- src/output/output.jl | 6 ++--- src/physics/column_variables.jl | 34 ++++++++++++------------- src/physics/define_column.jl | 30 +++++++++++----------- src/physics/large_scale_condensation.jl | 14 ++++++---- 5 files changed, 44 insertions(+), 42 deletions(-) diff --git a/src/SpeedyWeather.jl b/src/SpeedyWeather.jl index 3ca257fcb..194890462 100644 --- a/src/SpeedyWeather.jl +++ b/src/SpeedyWeather.jl @@ -168,6 +168,7 @@ include("physics/thermodynamics.jl") include("physics/boundary_layer.jl") include("physics/temperature_relaxation.jl") include("physics/vertical_diffusion.jl") +include("physics/large_scale_condensation.jl") include("physics/pretty_printing.jl") # MODELS @@ -175,7 +176,6 @@ include("dynamics/models.jl") # # PHYSICS # include("physics/convection.jl") -include("physics/large_scale_condensation.jl") # include("physics/longwave_radiation.jl") # include("physics/shortwave_radiation.jl") diff --git a/src/output/output.jl b/src/output/output.jl index f8e486f1e..7a62bad00 100644 --- a/src/output/output.jl +++ b/src/output/output.jl @@ -442,9 +442,9 @@ function write_netcdf_variables!( output::OutputWriter, precip_large_scale .= 0 precip_convection .= 0 - # convert from [m/s] to [mm/Δt], where Δt is the output time step - precip_cond *= 1000*output.output_dt_sec - precip_conv *= 1000*output.output_dt_sec + # convert from [m] to [mm] within output time step (e.g. 6hours) + precip_cond *= 1000 + precip_conv *= 1000 if Model <: PrimitiveEquation @. pres = exp(pres)/100 # convert from log(pₛ) to surface pressure pₛ [hPa] diff --git a/src/physics/column_variables.jl b/src/physics/column_variables.jl index 0f66464c2..33edc17f5 100644 --- a/src/physics/column_variables.jl +++ b/src/physics/column_variables.jl @@ -1,6 +1,5 @@ """ - get_column!(C,D,ij,G) - +$(TYPEDSIGNATURES) Update `C::ColumnVariables` by copying the prognostic variables from `D::DiagnosticVariables` at gridpoint index `ij`. Provide `G::Geometry` for coordinate information.""" function get_column!( C::ColumnVariables, @@ -31,6 +30,7 @@ function get_column!( C::ColumnVariables, C.humid[k] = layer.grid_variables.humid_grid[ij] # as well as geopotential (not actually prognostic though) + # TODO geopot on the grid is currently not computed in dynamics C.geopot[k] = layer.grid_variables.geopot_grid[ij] end end @@ -47,8 +47,7 @@ function get_column!( C::ColumnVariables, end """ - write_column_tendencies!(D,C,ij,G) - +$(TYPEDSIGNATURES) Write the parametrization tendencies from `C::ColumnVariables` into the horizontal fields of tendencies stored in `D::DiagnosticVariables` at gridpoint index `ij`.""" function write_column_tendencies!( D::DiagnosticVariables, @@ -72,27 +71,26 @@ function write_column_tendencies!( D::DiagnosticVariables, end """ - reset_column!(column::ColumnVariables) - +$(TYPEDSIGNATURES) Set the accumulators (tendencies but also vertical sums and similar) back to zero for `column` to be reused at other grid points.""" function reset_column!(column::ColumnVariables{NF}) where NF # set tendencies to 0 for += accumulation - fill!(column.u_tend,0) - fill!(column.v_tend,0) - fill!(column.temp_tend,0) - fill!(column.humid_tend,0) + column.u_tend .= 0 + column.v_tend .= 0 + column.temp_tend .= 0 + column.humid_tend .= 0 # set fluxes to 0 for += accumulation - fill!(column.flux_u_upward,0) - fill!(column.flux_u_downward,0) - fill!(column.flux_v_upward,0) - fill!(column.flux_v_downward,0) - fill!(column.flux_humid_upward,0) - fill!(column.flux_humid_downward,0) - fill!(column.flux_temp_upward,0) - fill!(column.flux_temp_downward,0) + column.flux_u_upward .= 0 + column.flux_u_downward .= 0 + column.flux_v_upward .= 0 + column.flux_v_downward .= 0 + column.flux_humid_upward .= 0 + column.flux_humid_downward .= 0 + column.flux_temp_upward .= 0 + column.flux_temp_downward .= 0 # # Convection # column.cloud_top = column.nlev+1 diff --git a/src/physics/define_column.jl b/src/physics/define_column.jl index 548f101c2..36ef159bf 100644 --- a/src/physics/define_column.jl +++ b/src/physics/define_column.jl @@ -17,14 +17,14 @@ Base.@kwdef mutable struct ColumnVariables{NF<:AbstractFloat} <: AbstractColumnV latd::NF = 0 # latitude, needed for shortwave radiation # PROGNOSTIC VARIABLES - const u::Vector{NF} = zeros(NF,nlev) # zonal velocity - const v::Vector{NF} = zeros(NF,nlev) # meridional velocity - const temp::Vector{NF} = zeros(NF,nlev) # temperature - const humid::Vector{NF} = zeros(NF,nlev) # specific humidity + const u::Vector{NF} = zeros(NF,nlev) # zonal velocity [m/s] + const v::Vector{NF} = zeros(NF,nlev) # meridional velocity [m/s] + const temp::Vector{NF} = zeros(NF,nlev) # temperature [K] + const humid::Vector{NF} = zeros(NF,nlev) # specific humidity [kg/kg] # (log) pressure per layer, surface is prognostic, last element here, but precompute other layers too - const ln_pres::Vector{NF} = zeros(NF,nlev+1) # logarithm of pressure [log(hPa)] - const pres::Vector{NF} = zeros(NF,nlev+1) # pressure [hPa] + const ln_pres::Vector{NF} = zeros(NF,nlev+1) # logarithm of pressure [log(Pa)] + const pres::Vector{NF} = zeros(NF,nlev+1) # pressure [Pa] # TENDENCIES to accumulate the parametrizations into const u_tend::Vector{NF} = zeros(NF,nlev) # zonal velocity [m] @@ -49,8 +49,8 @@ Base.@kwdef mutable struct ColumnVariables{NF<:AbstractFloat} <: AbstractColumnV const flux_humid_downward::Vector{NF} = zeros(NF,nlev+1) # THERMODYNAMICS - const sat_humid::Vector{NF} = zeros(NF,nlev) # Saturation specific humidity - const sat_vap_pres::Vector{NF} = zeros(NF,nlev) # Saturation vapour pressure + const sat_humid::Vector{NF} = zeros(NF,nlev) # Saturation specific humidity [kg/kg] + const sat_vap_pres::Vector{NF} = zeros(NF,nlev) # Saturation vapour pressure [Pa] const dry_static_energy::Vector{NF} = zeros(NF,nlev) # Dry static energy const moist_static_energy::Vector{NF} = zeros(NF,nlev) # Moist static energy @@ -62,18 +62,18 @@ Base.@kwdef mutable struct ColumnVariables{NF<:AbstractFloat} <: AbstractColumnV sat_moist_static_energy_half::Vector{NF} = zeros(NF,nlev) # Saturation moist static energy interpolated to half-levels # Convection - conditional_instability::Bool = false # Whether a conditional instability exists in this column (condition 1) - activate_convection::Bool = false # Whether convection should be activated in this column (condition 2) - cloud_top::Int = nlev+1 # Top-of-convection layer - excess_humidity::NF = 0 # Excess humidity due to convection - cloud_base_mass_flux::NF = 0 # Mass flux at the top of the PBL - precip_convection::NF = 0 # Precipitation due to convection + conditional_instability::Bool = false # in this column? (condition 1) + activate_convection::Bool = false # in this column? (condition 2) + cloud_top::Int = nlev+1 # Top-of-convection layer + excess_humidity::NF = 0 # due to convection + cloud_base_mass_flux::NF = 0 # Mass flux at the top of the PBL + precip_convection::NF = 0 # Precipitation due to convection [m] net_flux_humid::Vector{NF} = zeros(NF,nlev) # Net fluxes of moisture in this column net_flux_dry_static_energy::Vector{NF} = zeros(NF,nlev) # Net fluxes of dry static energy in this column entrainment_profile::Vector{NF} = zeros(NF,nlev) # Entrainment coefficients # Large-scale condensation - precip_large_scale::NF = 0 # Precipitation due to large-scale condensation + precip_large_scale::NF = 0 # precipitation [m] # Longwave radiation ## New vars in radlw_down! diff --git a/src/physics/large_scale_condensation.jl b/src/physics/large_scale_condensation.jl index d092ccf27..8533c4978 100644 --- a/src/physics/large_scale_condensation.jl +++ b/src/physics/large_scale_condensation.jl @@ -71,7 +71,7 @@ function large_scale_condensation!( model::PrimitiveWet, ) large_scale_condensation!(column,model.large_scale_condensation, - model.geometry,model.constants,model.atmosphere) + model.geometry,model.constants,model.atmosphere,model.time_stepping) end @@ -81,6 +81,7 @@ function large_scale_condensation!( geometry::Geometry, constants::DynamicsConstants, atmosphere::AbstractAtmosphere, + time_stepping::TimeStepper, ) where NF (;relative_threshold,humid_tend_max) = scheme @@ -94,6 +95,9 @@ function large_scale_condensation!( pₛ = pres[end] # surface pressure (;gravity, water_density) = constants + (;Δt_sec) = time_stepping + Δt_gρ = Δt_sec/gravity/water_density + (;σ_levels_thick) = geometry latent_heat = convert(NF, atmosphere.latent_heat_condensation/atmosphere.cₚ) @@ -106,16 +110,16 @@ function large_scale_condensation!( if humid[k] > humid_threshold # accumulate in tendencies (nothing is added if humidity not above threshold) humid_tend_k = -(humid[k] - humid_threshold) / time_scale # Formula 22 - temp_tend[k] += -latent_heat * min(humid_tend_k, humid_tend_max[k]*pres[k]) # Formula 23 + temp_tend[k] += -latent_heat * min(humid_tend_k, humid_tend_max[k]*pres[k]) # Formula 23 # If there is large-scale condensation at a level higher (i.e. smaller k) than # the cloud-top previously diagnosed due to convection, then increase the cloud-top - column.cloud_top = min(column.cloud_top, k) # Page 7 (last sentence) + column.cloud_top = min(column.cloud_top, k) # Page 7 (last sentence) # 2. Precipitation due to large-scale condensation [kg/m²/s] /ρ for [m/s] # += for vertical integral - Δpₖ_g = pₛ*σ_levels_thick[k]/gravity/water_density # Formula 4 - column.precip_large_scale += -Δpₖ_g * humid_tend_k # Formula 25 + Δpₖ_g = pₛ*σ_levels_thick[k]*Δt_gρ # Formula 4 *Δt for [m] of rain during Δt + column.precip_large_scale += -Δpₖ_g * humid_tend_k # Formula 25, unit [m] # only write into humid_tend now to allow humid_tend != 0 before this scheme is called humid_tend[k] += humid_tend_k