Atmospheric temperature profiles

Here, we plot the atmospheric reference state profiles for a few different polynomial orders and number of elements.

using ClimateMachine
const clima_dir = dirname(dirname(pathof(ClimateMachine)));
using Plots
include(joinpath(clima_dir, "docs", "plothelpers.jl"));
include(joinpath(clima_dir, "test", "Atmos", "Model", "get_atmos_ref_states.jl"));

function export_ref_state_plot(nelem_vert, N_poly)
    solver_config = get_atmos_ref_states(nelem_vert, N_poly, 0.5)
    z = get_z(solver_config.dg.grid)
    all_data = dict_of_nodal_states(solver_config, ["z"])
    T = all_data["ref_state.T"]
    ρ = all_data["ref_state.ρ"]
    p = all_data["ref_state.p"]
    ρe = all_data["ref_state.ρe"]
    p1 = plot(T, z./10^3, xlabel="Temperature [K]");
    p2 = plot(ρ, z./10^3, xlabel="Density [kg/m^3]");
    p3 = plot(p./10^3, z./10^3, xlabel="Pressure [kPa]");
    p4 = plot(ρe./10^3, z./10^3, xlabel="Total energy [kJ]");
    plot(p1, p2, p3, p4, layout=(1,4), ylabel="z [km]")
    savefig("N_poly_$(N_poly).png")
end

export_ref_state_plot(80, 1)
export_ref_state_plot(40, 2)
export_ref_state_plot(20, 4)
┌ Info: Model composition
│     param_set = Main.##ex-#346.EarthParameterSet()
│     problem = ClimateMachine.Atmos.AtmosProblem{ClimateMachine.Atmos.AtmosBC{ClimateMachine.Atmos.Impenetrable{ClimateMachine.Atmos.FreeSlip},ClimateMachine.Atmos.Insulating,ClimateMachine.Atmos.Impermeable,ClimateMachine.Atmos.ImpermeableTracer,ClimateMachine.TurbulenceConvection.NoTurbConvBC},Main.##ex-#346.var"#7#8",typeof(ClimateMachine.Atmos.atmos_problem_init_state_auxiliary)}(ClimateMachine.Atmos.AtmosBC{ClimateMachine.Atmos.Impenetrable{ClimateMachine.Atmos.FreeSlip},ClimateMachine.Atmos.Insulating,ClimateMachine.Atmos.Impermeable,ClimateMachine.Atmos.ImpermeableTracer,ClimateMachine.TurbulenceConvection.NoTurbConvBC}(ClimateMachine.Atmos.Impenetrable{ClimateMachine.Atmos.FreeSlip}(ClimateMachine.Atmos.FreeSlip()), ClimateMachine.Atmos.Insulating(), ClimateMachine.Atmos.Impermeable(), ClimateMachine.Atmos.ImpermeableTracer(), ClimateMachine.TurbulenceConvection.NoTurbConvBC()), Main.##ex-#346.var"#7#8"(), ClimateMachine.Atmos.atmos_problem_init_state_auxiliary)
│     orientation = ClimateMachine.Orientations.FlatOrientation()
│     ref_state = ClimateMachine.Atmos.HydrostaticState{ClimateMachine.TemperatureProfiles.DecayingTemperatureProfile{Float64},Float64}(ClimateMachine.TemperatureProfiles.DecayingTemperatureProfile{Float64}(290.0, 220.0, 8484.271021693852), 0.5)
│     turbulence = ClimateMachine.TurbulenceClosures.SmagorinskyLilly{Float64}(0.21)
│     turbconv = ClimateMachine.TurbulenceConvection.NoTurbConv()
│     hyperdiffusion = ClimateMachine.TurbulenceClosures.NoHyperDiffusion()
│     viscoussponge = ClimateMachine.TurbulenceClosures.NoViscousSponge()
│     moisture = ClimateMachine.Atmos.EquilMoist{Float64}(3, 0.1)
│     precipitation = ClimateMachine.Atmos.NoPrecipitation()
│     radiation = ClimateMachine.Atmos.NoRadiation()
│     source = (ClimateMachine.Atmos.Gravity(), ClimateMachine.Atmos.Coriolis(), ClimateMachine.Atmos.GeostrophicForcing{Float64}(7.62e-5, 0.0, 0.0))
│     tracers = ClimateMachine.Atmos.NoTracers()
└     data_config = nothing
┌ Info: Establishing single stack configuration for ref_state
│     precision        = Float64
│     polynomial order = 1
│     domain_min       = 0.00 m x0.00 m x0.00 m
│     domain_max       = 1.00 m x1.00 m x25000.00 m
│     #vert elems      = 80
│     MPI ranks        = 1
│     min(Δ_horz)      = 1.00 m
└     min(Δ_vert)      = 312.50 m
[ Info: Initializing ref_state
┌ Info: Model composition
│     param_set = Main.##ex-#346.EarthParameterSet()
│     problem = ClimateMachine.Atmos.AtmosProblem{ClimateMachine.Atmos.AtmosBC{ClimateMachine.Atmos.Impenetrable{ClimateMachine.Atmos.FreeSlip},ClimateMachine.Atmos.Insulating,ClimateMachine.Atmos.Impermeable,ClimateMachine.Atmos.ImpermeableTracer,ClimateMachine.TurbulenceConvection.NoTurbConvBC},Main.##ex-#346.var"#7#8",typeof(ClimateMachine.Atmos.atmos_problem_init_state_auxiliary)}(ClimateMachine.Atmos.AtmosBC{ClimateMachine.Atmos.Impenetrable{ClimateMachine.Atmos.FreeSlip},ClimateMachine.Atmos.Insulating,ClimateMachine.Atmos.Impermeable,ClimateMachine.Atmos.ImpermeableTracer,ClimateMachine.TurbulenceConvection.NoTurbConvBC}(ClimateMachine.Atmos.Impenetrable{ClimateMachine.Atmos.FreeSlip}(ClimateMachine.Atmos.FreeSlip()), ClimateMachine.Atmos.Insulating(), ClimateMachine.Atmos.Impermeable(), ClimateMachine.Atmos.ImpermeableTracer(), ClimateMachine.TurbulenceConvection.NoTurbConvBC()), Main.##ex-#346.var"#7#8"(), ClimateMachine.Atmos.atmos_problem_init_state_auxiliary)
│     orientation = ClimateMachine.Orientations.FlatOrientation()
│     ref_state = ClimateMachine.Atmos.HydrostaticState{ClimateMachine.TemperatureProfiles.DecayingTemperatureProfile{Float64},Float64}(ClimateMachine.TemperatureProfiles.DecayingTemperatureProfile{Float64}(290.0, 220.0, 8484.271021693852), 0.5)
│     turbulence = ClimateMachine.TurbulenceClosures.SmagorinskyLilly{Float64}(0.21)
│     turbconv = ClimateMachine.TurbulenceConvection.NoTurbConv()
│     hyperdiffusion = ClimateMachine.TurbulenceClosures.NoHyperDiffusion()
│     viscoussponge = ClimateMachine.TurbulenceClosures.NoViscousSponge()
│     moisture = ClimateMachine.Atmos.EquilMoist{Float64}(3, 0.1)
│     precipitation = ClimateMachine.Atmos.NoPrecipitation()
│     radiation = ClimateMachine.Atmos.NoRadiation()
│     source = (ClimateMachine.Atmos.Gravity(), ClimateMachine.Atmos.Coriolis(), ClimateMachine.Atmos.GeostrophicForcing{Float64}(7.62e-5, 0.0, 0.0))
│     tracers = ClimateMachine.Atmos.NoTracers()
└     data_config = nothing
┌ Info: Establishing single stack configuration for ref_state
│     precision        = Float64
│     polynomial order = 2
│     domain_min       = 0.00 m x0.00 m x0.00 m
│     domain_max       = 1.00 m x1.00 m x25000.00 m
│     #vert elems      = 40
│     MPI ranks        = 1
│     min(Δ_horz)      = 0.50 m
└     min(Δ_vert)      = 312.50 m
[ Info: Initializing ref_state
┌ Info: Model composition
│     param_set = Main.##ex-#346.EarthParameterSet()
│     problem = ClimateMachine.Atmos.AtmosProblem{ClimateMachine.Atmos.AtmosBC{ClimateMachine.Atmos.Impenetrable{ClimateMachine.Atmos.FreeSlip},ClimateMachine.Atmos.Insulating,ClimateMachine.Atmos.Impermeable,ClimateMachine.Atmos.ImpermeableTracer,ClimateMachine.TurbulenceConvection.NoTurbConvBC},Main.##ex-#346.var"#7#8",typeof(ClimateMachine.Atmos.atmos_problem_init_state_auxiliary)}(ClimateMachine.Atmos.AtmosBC{ClimateMachine.Atmos.Impenetrable{ClimateMachine.Atmos.FreeSlip},ClimateMachine.Atmos.Insulating,ClimateMachine.Atmos.Impermeable,ClimateMachine.Atmos.ImpermeableTracer,ClimateMachine.TurbulenceConvection.NoTurbConvBC}(ClimateMachine.Atmos.Impenetrable{ClimateMachine.Atmos.FreeSlip}(ClimateMachine.Atmos.FreeSlip()), ClimateMachine.Atmos.Insulating(), ClimateMachine.Atmos.Impermeable(), ClimateMachine.Atmos.ImpermeableTracer(), ClimateMachine.TurbulenceConvection.NoTurbConvBC()), Main.##ex-#346.var"#7#8"(), ClimateMachine.Atmos.atmos_problem_init_state_auxiliary)
│     orientation = ClimateMachine.Orientations.FlatOrientation()
│     ref_state = ClimateMachine.Atmos.HydrostaticState{ClimateMachine.TemperatureProfiles.DecayingTemperatureProfile{Float64},Float64}(ClimateMachine.TemperatureProfiles.DecayingTemperatureProfile{Float64}(290.0, 220.0, 8484.271021693852), 0.5)
│     turbulence = ClimateMachine.TurbulenceClosures.SmagorinskyLilly{Float64}(0.21)
│     turbconv = ClimateMachine.TurbulenceConvection.NoTurbConv()
│     hyperdiffusion = ClimateMachine.TurbulenceClosures.NoHyperDiffusion()
│     viscoussponge = ClimateMachine.TurbulenceClosures.NoViscousSponge()
│     moisture = ClimateMachine.Atmos.EquilMoist{Float64}(3, 0.1)
│     precipitation = ClimateMachine.Atmos.NoPrecipitation()
│     radiation = ClimateMachine.Atmos.NoRadiation()
│     source = (ClimateMachine.Atmos.Gravity(), ClimateMachine.Atmos.Coriolis(), ClimateMachine.Atmos.GeostrophicForcing{Float64}(7.62e-5, 0.0, 0.0))
│     tracers = ClimateMachine.Atmos.NoTracers()
└     data_config = nothing
┌ Info: Establishing single stack configuration for ref_state
│     precision        = Float64
│     polynomial order = 4
│     domain_min       = 0.00 m x0.00 m x0.00 m
│     domain_max       = 1.00 m x1.00 m x25000.00 m
│     #vert elems      = 20
│     MPI ranks        = 1
│     min(Δ_horz)      = 0.17 m
└     min(Δ_vert)      = 215.84 m
[ Info: Initializing ref_state

Polynomial order 1, 80 elements

Polynomial order 2, 40 elements

Polynomial order 4, 20 elements