Failure on saturation adjustment for certain namelist parameters & cases

Description: Certain combinations of namelist parameters lead to DomainError or maxiter reached in saturation_adjustment_given_pθq, stemming from values of H<0 being fed into the PhaseEquil_pθq function when performing saturation adjustment. Occasionally θ_liq_ice=NaN, q_tot=NaN, T=NaN as well.

https://github.com/CliMA/TurbulenceConvection.jl/blob/971a265708beabbd0a1d359ba1ab3a0ae07994d2/src/update_aux.jl#L82

This is encountered in the first few iterations of the simulation.

TODO: Investigate root cause.

Stacktrace:

  DomainError with -25.191093523420196:
  Exponentiation yielding a complex result requires a complex argument.
  Replace x^y with (x+0im)^y, Complex(x)^y, or similar.
  Stacktrace:
    [1] throw_exp_domainerror(x::Float64)
      @ Base.Math ./math.jl:37
    [2] ^
      @ ./math.jl:901 [inlined]
    [3] saturation_vapor_pressure(param_set::EarthParameterSet{NamedTuple{(:MSLP, :cp_d, :cp_v, :R_d, :R_v, :molmass_ratio, :T_freeze, :τ_precip, :τ_acnv, :c_ε, :c_div, :α_b, :α_a, :α_d, :H_up_min, :ω_pr, :c_δ, :β, :χ, :c_t, :c_λ, :w_min, :μ_0, :μ, :c_m, :c_d, :c_b, :Pr_n, :smin_ub, :smin_rm, :l_max, :stoch_ε_lognormal_var, :stoch_δ_lognormal_var, :sde_ϵ_θ, :sde_ϵ_σ, :sde_δ_θ, :sde_δ_σ), NTuple{37, Float64}}}, T::Float64, LH_0::Float64, Δcp::Float64)
      @ Thermodynamics ~/.julia/packages/Thermodynamics/OfVhK/src/relations.jl:878
    [4] saturation_vapor_pressure
      @ ~/.julia/packages/Thermodynamics/OfVhK/src/relations.jl:864 [inlined]
    [5] q_vap_saturation_from_pressure
      @ ~/.julia/packages/Thermodynamics/OfVhK/src/relations.jl:1038 [inlined]
    [6] PhasePartition_equil_given_p(param_set::EarthParameterSet{NamedTuple{(:MSLP, :cp_d, :cp_v, :R_d, :R_v, :molmass_ratio, :T_freeze, :τ_precip, :τ_acnv, :c_ε, :c_div, :α_b, :α_a, :α_d, :H_up_min, :ω_pr, :c_δ, :β, :χ, :c_t, :c_λ, :w_min, :μ_0, :μ, :c_m, :c_d, :c_b, :Pr_n, :smin_ub, :smin_rm, :l_max, :stoch_ε_lognormal_var, :stoch_δ_lognormal_var, :sde_ϵ_θ, :sde_ϵ_σ, :sde_δ_θ, :sde_δ_σ), NTuple{37, Float64}}}, T::Float64, p::Float64, q_tot::Float64, phase_type::Type{Thermodynamics.PhaseEquil})
      @ Thermodynamics ~/.julia/packages/Thermodynamics/OfVhK/src/relations.jl:1260

. . .

      @ TurbulenceConvection /central~/clima/TurbulenceConvection.jl/src/update_aux.jl:84
   [14] update(edmf::TurbulenceConvection.EDMF_PrognosticTKE{EarthParameterSet{NamedTuple{(:MSLP, :cp_d, :cp_v, :R_d, :R_v, :molmass_ratio, :T_freeze, :τ_precip, :τ_acnv, :c_ε, :c_div, :α_b, :α_a, :α_d, :H_up_min, :ω_pr, :c_δ, :β, :χ, :c_t, :c_λ, :w_min, :μ_0, :μ, :c_m, :c_d, :c_b, :Pr_n, :smin_ub, :smin_rm, :l_max, :stoch_ε_lognormal_var, :stoch_δ_lognormal_var, :sde_ϵ_θ, :sde_ϵ_σ, :sde_δ_θ, :sde_δ_σ), NTuple{37, Float64}}}, Vector{Float64}, Matrix{Float64}, NamedTuple{(:A_θq_gm, :A_uv_gm, :A_TKE, :A_Hvar, :A_QTvar, :A_HQTcov, :b_θ_liq_ice_gm, :b_q_tot_gm, :b_u_gm, :b_v_gm, :b_TKE, :b_Hvar, :b_QTvar, :b_HQTcov), Tuple{LinearAlgebra.Tridiagonal{Float64, Vector{Float64}}, LinearAlgebra.Tridiagonal{Float64, Vector{Float64}}, LinearAlgebra.Tridiagonal{Float64, Vector{Float64}}, LinearAlgebra.Tridiagonal{Float64, Vector{Float64}}, LinearAlgebra.Tridiagonal{Float64, Vector{Float64}}, LinearAlgebra.Tridiagonal{Float64, Vector{Float64}}, Vector{Float64}, Vector{Float64}, Vector{Float64}, Vector{Float64}, Vector{Float64}, Vector{Float64}, Vector{Float64}, Vector{Float64}}}}, GMV::TurbulenceConvection.GridMeanVariables{EarthParameterSet{NamedTuple{(:MSLP, :cp_d, :cp_v, :R_d, :R_v, :molmass_ratio, :T_freeze, :τ_precip, :τ_acnv, :c_ε, :c_div, :α_b, :α_a, :α_d, :H_up_min, :ω_pr, :c_δ, :β, :χ, :c_t, :c_λ, :w_min, :μ_0, :μ, :c_m, :c_d, :c_b, :Pr_n, :smin_ub, :smin_rm, :l_max, :stoch_ε_lognormal_var, :stoch_δ_lognormal_var, :sde_ϵ_θ, :sde_ϵ_σ, :sde_δ_θ, :sde_δ_σ), NTuple{37, Float64}}}}, Case::TurbulenceConvection.CasesBase{Main.Cases.Bomex}, TS::TurbulenceConvection.TimeStepping)
      @ TurbulenceConvection /central~/clima/TurbulenceConvection.jl/src/Turbulence_PrognosticTKE.jl:431
   [15] run(self::Simulation1d)
      @ Main /central~/clima/TurbulenceConvection.jl/integration_tests/utils/main.jl:148
   [16] main1d(namelist::Dict{Any, Any}; time_run::Bool)
      @ Main /central~/clima/TurbulenceConvection.jl/integration_tests/utils/main.jl:226

. . .

CliMA / TurbulenceConvection.jl

Failure on saturation adjustment for certain namelist parameters & cases #356