diff --git a/docs/bibliography.bib b/docs/bibliography.bib index bd226707d..fb88ca5eb 100644 --- a/docs/bibliography.bib +++ b/docs/bibliography.bib @@ -46,6 +46,14 @@ @Article{Alpert2022 doi = "10.1039/D1EA00077B", } +@article{BarklieGokhale1959, + title = {The freezing of supercooled wter drops. Alberta hal, 1958, and related studies}, + author = {Barklie, R. H. D. and Gokhale, N. R.}, + journal = {McGill University Stormy Weather Group Sci. Rep. MW-30}, + pages = {43-64}, + year = {1959}, +} + @Article{Baumgartner2022, author = {Baumgartner, M. and Rolf, C. and Groo{\ss}, J.-U. and Schneider, J. and Schorr, T. and M\"ohler, O. and Spichtinger, P. and Kr\"amer, M.}, title = {New investigations on homogeneous ice nucleation: the effects of water activity and water saturation formulations}, @@ -69,6 +77,16 @@ @article{Beheng1994 doi = {10.1016/0169-8095(94)90020-5} } +@article{Bigg1953, + title = {The supercooling of water.}, + author = {Bigg, E.K.}, + journal = {Proc. Phys. Soc.}, + volume = {66B}, + pages = {688-694}, + year = {1953}, + doi = {10.1088/0370-1301/66/8/309} +} + @article{BrownFrancis1995, title = {Improved Measurements of the Ice Water Content in Cirrus Using a Total-Water Probe}, author = {Philip R. A. Brown and Peter N. Francis}, @@ -599,6 +617,17 @@ @article{Spichtinger2023 DOI = {10.5194/acp-23-2035-2023} } +@article{Thompson2004, + author = {Gregory Thompson and Roy M. Rasmussen and Kevin Manning}, + title = {Explicit Forecasts of Winter Precipitation Using an Improved Bulk Microphysics Scheme. Part I: Description and Sensitivity Analysis}, + journal = {Monthly Weather Review}, + year = {2004}, + volume = {132}, + number = {2}, + doi = {10.1175/1520-0493(2004)132<0519:EFOWPU>2.0.CO;2}, + pages = {519 - 542}, +} + @article{TripoliCotton1980, title = {A Numerical Investigation of Several Factors Contributing to the Observed Variable Intensity of Deep Convection over South Florida}, author = {Tripoli, G.J. and Cotton, W.R.}, diff --git a/docs/src/API.md b/docs/src/API.md index 2b1a970a3..6b7e209f9 100644 --- a/docs/src/API.md +++ b/docs/src/API.md @@ -90,6 +90,8 @@ HetIceNucleation.dust_activated_number_fraction HetIceNucleation.MohlerDepositionRate HetIceNucleation.deposition_J HetIceNucleation.ABIFM_J +HetIceNucleation.P3_deposition_N_i +HetIceNucleation.P3_het_N_i HetIceNucleation.INP_concentration_frequency ``` diff --git a/docs/src/IceNucleation.md b/docs/src/IceNucleation.md index 12890c908..65eee3764 100644 --- a/docs/src/IceNucleation.md +++ b/docs/src/IceNucleation.md @@ -63,7 +63,7 @@ Limited experimental values for the free parameters ``a`` and ``S_0`` can be fou freezing occurs in a different ice nucleation mode (either a second deposition or other condensation type mode). -## Water Activity Based Deposition Nucleation +### Water Activity Based Deposition Nucleation The water activity based deposition nucleation model is analagous to ABIFM for immersion freezing (see [ABIFM for Sulphuric Acid Containing Droplets](https://clima.github.io/CloudMicrophysics.jl/dev/IceNucleation/#ABIFM-for-Sulphuric-Acid-Containing-Droplets) section below). It calculates a nucleation rate coefficient, ``J``, which @@ -94,7 +94,22 @@ include("plots/activity_based_deposition.jl") ![](water_activity_depo_nuc.svg) -## ABIFM for Sulphuric Acid Containing Droplets +### Cooper (1986) Ice Crystal Number +The P3 scheme as described in [MorrisonMilbrandt2015](@cite) follows + the implementation of [Thompson2004](@cite) where a parameterization from + Cooper 1986 was used. Number of ice crystals formed is derived from + ambient temperature. +```math +\begin{equation} + N_i = 0.005 exp[0.304(T_0 - T)] +\end{equation} +``` +Where ``T_0 = 273.15K`` and ``T`` is the ambient temperature. +Because the parameterization is prone to overpredict at low temperatures, + ``N_i = N_i(T = 233K)`` for all values of ``N_i`` at which ``T < 233K``. + +## Immersion Freezing +### ABIFM for Sulphuric Acid Containing Droplets Water Activity-Based Immersion Freezing Model (ABFIM) is a method of parameterizing immersion freezing inspired by the time-dependent classical nucleation theory (CNT). More on CNT can be found in [Karthika2016](@cite). @@ -169,7 +184,23 @@ It is also important to note that this plot is reflective of cirrus clouds and shows only a very small temperature range. The two curves are slightly off because of small differences in parameterizations for vapor pressures. -## Homogeneous Freezing for Sulphuric Acid Containing Droplets +### Bigg (1953) Volume and Time Dependent Heterogeneous Freezing +Heterogeneous freezing in the P3 scheme as described in [MorrisonMilbrandt2015](@cite) + follows the parameterization from [Bigg1953](@cite) with parameters from + [BarklieGokhale1959](@cite). The number of ice nucleated in a timestep via + heterogeneous freezing is determined by +```math +\begin{equation} + N_{ice} = N_{liq} \left[ 1 - exp(-B V_{l} \Delta t exp(aT_s)) \right] +\end{equation} +``` +where `a` and `B` are parameters taken from [BarklieGokhale1959](@cite) for + rainwater as 0.65 and 2e-4 respectively. `T_s` is the difference between + 273.15K and ambient temperature. `V_l` is the volume of droplets to be + frozen and `\Delta t` is timestep in which freezing occurs. + +## Homogeneous Freezing +### Homogeneous Freezing for Sulphuric Acid Containing Droplets Homogeneous freezing occurs when supercooled liquid droplets freeze on their own. Closly based off [Koop2000](@cite), this parameterization determines a homoegneous nucleation rate coefficient, ``J_{hom}``, using water activity. The change in water activity, @@ -185,7 +216,7 @@ The nucleation rate coefficient is determined with the cubic function from [Koop ``` This parameterization is valid only when ``0.26 < \Delta a_w < 0.36`` and ``185K < T < 235K``. -### Homogeneous Ice Nucleation Rate Coefficient +### Homogeneous Freezing Example Figures Here is a comparison of our parameterization of ``J_{hom}`` compared to Koop 2000 as plotted in figure 1 of [Spichtinger2023](@cite). Our parameterization differs in the calculation of ``\Delta a_w``. We define water activity to be a ratio of saturated vapor pressures whereas diff --git a/docs/src/IceNucleationParcel0D.md b/docs/src/IceNucleationParcel0D.md index bf70908a4..5aeb800d5 100644 --- a/docs/src/IceNucleationParcel0D.md +++ b/docs/src/IceNucleationParcel0D.md @@ -315,3 +315,21 @@ The following plots show the parcel model running with homogeneous freezing and include("../../parcel/Jensen_et_al_2022.jl") ``` ![](Jensen_et_al_2022.svg) + +## P3 Ice Nucleation Parameterizations +The parcel also includes ice nucleation parameterizations used in + the P3 scheme as described in [MorrisonMilbrandt2015](@cite). + Deposition nucleation is based on the ice crystal number parameterization + from Cooper (1986). The heterogeneous freezing parameterization, which + follows Bigg(1953) with parameters from Barklie aand Gokhale (1959), is + treated as immersion freezing in the parcel. Homogeneous freezing happens + instantaneously at 233.15K. +Shown below are three separate parcel simulations for deposition nucleation, + immersion freezing, and homogeneous freezing. Note that initial temperature + varies for each run. The deposition nucleation run does not conserve + INP number, while the other two freezing modes do. Updraft velocity is + set to 0.5 m/s. +```@example +include("../../parcel/P3_ice_nuc.jl") +``` +![](P3_ice_nuc.svg) diff --git a/parcel/P3_ice_nuc.jl b/parcel/P3_ice_nuc.jl new file mode 100644 index 000000000..032a465e5 --- /dev/null +++ b/parcel/P3_ice_nuc.jl @@ -0,0 +1,125 @@ +import OrdinaryDiffEq as ODE +import CairoMakie as MK +import Thermodynamics as TD +import CloudMicrophysics as CM +import CLIMAParameters as CP + +# definition of the ODE problem for parcel model +include(joinpath(pkgdir(CM), "parcel", "parcel.jl")) +FT = Float32 +# get free parameters +tps = TD.Parameters.ThermodynamicsParameters(FT) +aps = CMP.AirProperties(FT) +wps = CMP.WaterProperties(FT) +ip = CMP.IceNucleationParameters(FT) + +# Constants +ρₗ = wps.ρw +R_v = TD.Parameters.R_v(tps) +R_d = TD.Parameters.R_d(tps) + +# Initial conditions +Nₐ = FT(2000) +Nₗ = FT(2000) +Nᵢ = FT(0) +rₗ = FT(1.25e-6) +p₀ = FT(20000) +T₀_dep = FT(235) +T₀_het = FT(235) +T₀_hom = FT(233.2) +qᵥ = FT(8.3e-4) +qₗ = FT(Nₗ * 4 / 3 * π * rₗ^3 * ρₗ / 1.2) +qᵢ = FT(0) +x_sulph = FT(0) + +# Moisture dependent initial conditions +q = TD.PhasePartition(qᵥ + qₗ + qᵢ, qₗ, qᵢ) +ts_dep = TD.PhaseNonEquil_pTq(tps, p₀, T₀_dep, q) +ts_het = TD.PhaseNonEquil_pTq(tps, p₀, T₀_het, q) +ts_hom = TD.PhaseNonEquil_pTq(tps, p₀, T₀_hom, q) +ρₐ_dep = TD.air_density(tps, ts_dep) +ρₐ_het = TD.air_density(tps, ts_het) +ρₐ_hom = TD.air_density(tps, ts_hom) +Rₐ = TD.gas_constant_air(tps, q) +eₛ_dep = TD.saturation_vapor_pressure(tps, T₀_dep, TD.Liquid()) +eₛ_het = TD.saturation_vapor_pressure(tps, T₀_het, TD.Liquid()) +eₛ_hom = TD.saturation_vapor_pressure(tps, T₀_hom, TD.Liquid()) +e = qᵥ * p₀ * R_v / Rₐ + +Sₗ_dep = FT(e / eₛ_dep) +Sₗ_het = FT(e / eₛ_het) +Sₗ_hom = FT(e / eₛ_hom) +IC_dep = [Sₗ_dep, p₀, T₀_dep, qᵥ, qₗ, qᵢ, Nₐ, Nₗ, Nᵢ, x_sulph] +IC_het = [Sₗ_het, p₀, T₀_het, qᵥ, qₗ, qᵢ, Nₐ, Nₗ, Nᵢ, x_sulph] +IC_hom = [Sₗ_hom, p₀, T₀_hom, qᵥ, qₗ, qᵢ, Nₐ, Nₗ, Nᵢ, x_sulph] + +ξ(T) = + TD.saturation_vapor_pressure(tps, T, TD.Liquid()) / + TD.saturation_vapor_pressure(tps, T, TD.Ice()) +S_i(T, S_liq) = ξ(T) * S_liq + +# Simulation parameters passed into ODE solver +r_nuc = FT(1.25e-6) # assumed size of nucleated particles +w = FT(0.5) # updraft speed +α_m = FT(0.5) # accomodation coefficient +const_dt = FT(0.1) # model timestep +t_max = FT(50) +aerosol = [] +ice_nucleation_modes_list = [["P3_Deposition"], ["P3_het"], ["P3_hom"]] +growth_modes = ["Deposition"] +droplet_size_distribution_list = [["Monodisperse"]] + +# Plotting +fig = MK.Figure(resolution = (900, 600)) +ax1 = MK.Axis(fig[1, 1], ylabel = "Ice Saturation [-]") +ax2 = MK.Axis(fig[1, 2], ylabel = "ICNC [cm^-3]") +ax7 = MK.Axis(fig[1, 3], ylabel = "Temperature [K]") +ax3 = MK.Axis(fig[2, 1], ylabel = "Ice Saturation [-]") +ax4 = MK.Axis(fig[2, 2], ylabel = "ICNC [cm^-3]") +ax8 = MK.Axis(fig[2, 3], ylabel = "Temperature [K]") +ax5 = MK.Axis(fig[3, 1], ylabel = "Ice Saturation [-]", xlabel = "Time [s]") +ax6 = MK.Axis(fig[3, 2], ylabel = "ICNC [cm^-3]", xlabel = "Time [s]") +ax9 = MK.Axis(fig[3, 3], ylabel = "Temperature [K]", xlabel = "Time [s]") + +for ice_nucleation_modes in ice_nucleation_modes_list + nuc_mode = ice_nucleation_modes[1] + droplet_size_distribution = droplet_size_distribution_list[1] + p = (; + wps, + aps, + tps, + ip, + const_dt, + r_nuc, + w, + α_m, + aerosol, + ice_nucleation_modes, + growth_modes, + droplet_size_distribution, + ) + # solve ODE + #! format: off + if "P3_Deposition" in ice_nucleation_modes + sol = run_parcel(IC_dep, FT(0), t_max, p) + MK.lines!(ax1, sol.t, S_i.(sol[3, :], (sol[1, :])), label = nuc_mode) # saturation + MK.lines!(ax2, sol.t, sol[9, :] * 1e-6) # ICNC + MK.lines!(ax7, sol.t, sol[3, :]) # Temperature + MK.axislegend(ax1, framevisible = false, labelsize = 12, orientation = :horizontal, position = :rt) + elseif "P3_het" in ice_nucleation_modes + sol = run_parcel(IC_het, FT(0), t_max, p) + MK.lines!(ax3, sol.t, S_i.(sol[3, :], (sol[1, :])), label = nuc_mode, color = :red) # saturation + MK.lines!(ax4, sol.t, sol[9, :] * 1e-6, color = :red) # ICNC + MK.lines!(ax8, sol.t, sol[3, :], color = :red) # Temperature + MK.axislegend(ax3, framevisible = false, labelsize = 12, orientation = :horizontal, position = :lt) + elseif "P3_hom" in ice_nucleation_modes + sol = run_parcel(IC_hom, FT(0), t_max, p) + MK.lines!(ax5, sol.t, S_i.(sol[3, :], (sol[1, :])), label = nuc_mode, color = :orange) # saturation + MK.lines!(ax6, sol.t, sol[9, :] * 1e-6, color = :orange) # ICNC + MK.lines!(ax9, sol.t, sol[3, :], color = :orange) # Temperature + MK.axislegend(ax5, framevisible = false, labelsize = 12, orientation = :horizontal, position = :lt) + end + #! format: on +end + +MK.save("P3_ice_nuc.svg", fig) diff --git a/parcel/parcel.jl b/parcel/parcel.jl index b85a8bb17..56dc3b863 100644 --- a/parcel/parcel.jl +++ b/parcel/parcel.jl @@ -107,6 +107,11 @@ function parcel_model(dY, Y, p, t) dqi_dt_new_depo = dN_act_dt_depo * 4 / 3 * FT(π) * r_nuc^3 * ρ_ice / ρ_air end end + if "P3_Deposition" in ice_nucleation_modes + Nᵢ_depo = CMI_het.P3_deposition_N_i(ip.p3, T) + dN_act_dt_depo = max(FT(0), (Nᵢ_depo - N_ice) / const_dt) + dqi_dt_new_depo = dN_act_dt_depo * 4 / 3 * FT(π) * r_nuc^3 * ρ_ice / ρ_air + end dN_act_dt_immersion = FT(0) dqi_dt_new_immers = FT(0) @@ -126,6 +131,19 @@ function parcel_model(dY, Y, p, t) dN_act_dt_immersion = max(FT(0), J_immersion * N_liq * A_l) dqi_dt_new_immers = dN_act_dt_immersion * 4 / 3 * FT(π) * r_l^3 * ρ_ice / ρ_air end + if "P3_het" in ice_nucleation_modes + if "Monodisperse" in droplet_size_distribution + r_l = cbrt(q_liq / N_liq / (4 / 3 * FT(π)) / ρ_liq * ρ_air) + elseif "Gamma" in droplet_size_distribution + λ = cbrt(32 * FT(π) * N_liq / q_liq * ρ_liq / ρ_air) + #A = N_liq* λ^2 + r_l = 2 / λ + end + V_l = FT(4 / 3 * π * r_l^3) + Nᵢ_het = CMI_het.P3_het_N_i(ip.p3, T, N_liq, V_l, const_dt) + dN_act_dt_immersion = max(FT(0), (Nᵢ_het - N_ice) / const_dt) + dqi_dt_new_immers = dN_act_dt_immersion * 4 / 3 * FT(π) * r_l^3 * ρ_ice / ρ_air + end dN_act_dt_hom = FT(0) dqi_dt_new_hom = FT(0) @@ -159,6 +177,23 @@ function parcel_model(dY, Y, p, t) dqi_dt_new_hom = dN_act_dt_hom * 4 / 3 * FT(π) * r_aero^3 * ρ_ice / ρ_air end end + if "P3_hom" in ice_nucleation_modes + if "Monodisperse" in droplet_size_distribution + r_l = cbrt(q_liq / N_liq / (4 / 3 * FT(π)) / ρ_liq * ρ_air) + elseif "Gamma" in droplet_size_distribution + λ = cbrt(32 * FT(π) * N_liq / q_liq * ρ_liq / ρ_air) + #A = N_liq* λ^2 + r_l = 2 / λ + elseif "Lognormal" in droplet_size_distribution + r_l = R_mode_liq * exp(σ) + end + if T < FT(233.15) && N_liq > FT(0) + dN_act_dt_hom = N_liq / const_dt + else + dN_act_dt_hom = FT(0) + end + dqi_dt_new_hom = dN_act_dt_hom * 4 / 3 * FT(π) * r_l^3 * ρ_ice / ρ_air + end dN_ice_dt = dN_act_dt_depo + dN_act_dt_immersion + dN_act_dt_hom if "Jensen_Example" in droplet_size_distribution @@ -306,12 +341,25 @@ function run_parcel(IC, t_0, t_end, p) print("\n Water activity based deposition nucleation ") print("(Note that this mode only runs for monodisperse size distributions.) ") end + if "P3_Deposition" in ice_nucleation_modes + print("\n P3 deposition nucleation ") + timestepper = ODE.Euler() + end if "ImmersionFreezing" in ice_nucleation_modes - print("\n Immersion freezing ") + print("\n Immersion freezing (ABIFM) ") + print("(Note that this mode only runs for monodisperse and gamma size distributions.) ") + end + if "P3_het" in ice_nucleation_modes + print("\n P3 heterogeneous freezing ") print("(Note that this mode only runs for monodisperse and gamma size distributions.) ") + timestepper = ODE.Euler() end if "HomogeneousFreezing" in ice_nucleation_modes - print("\n Homogeneous freezing ") + print("\n Water activity based homogeneous freezing ") + end + if "P3_hom" in ice_nucleation_modes + print("\n P3 homogeneous freezing ") + timestepper = ODE.Euler() end print("\n") @@ -338,8 +386,14 @@ function run_parcel(IC, t_0, t_end, p) if "Lognormal" in droplet_size_distribution print("\n Lognormal size distribution") end - if "MohlerAF_Deposition" in ice_nucleation_modes || "MohlerRate_Deposition" in ice_nucleation_modes || "ActivityBasedDeposition" in ice_nucleation_modes - print("\nWARNING: Multiple deposition nucleation parameterizations chosen to run in one simulation. Parcel will only run MohlerAF_Deposition for now.") + if "MohlerAF_Deposition" in ice_nucleation_modes || + "MohlerRate_Deposition" in ice_nucleation_modes || + "ActivityBasedDeposition" in ice_nucleation_modes || + "P3_Deposition" in ice_nucleation_modes + print("\nWARNING: Multiple deposition nucleation parameterizations chosen to run in one simulation. Parcel can only properly handle one for now.") + end + if "P3_het" in ice_nucleation_modes + print("\nWARNING: Assuming P3_het parameterization is a parameterization for immeresion freezing.") end println(" ") diff --git a/src/IceNucleation.jl b/src/IceNucleation.jl index 78724fcdc..f702515b1 100644 --- a/src/IceNucleation.jl +++ b/src/IceNucleation.jl @@ -10,6 +10,8 @@ export dust_activated_number_fraction export MohlerDepositionRate export deposition_J export ABIFM_J +export P3_deposition_N_i +export P3_het_N_i export INP_concentration_frequency """ @@ -109,6 +111,59 @@ function ABIFM_J( return max(FT(0), FT(10)^logJ * FT(1e4)) # converts cm^-2 s^-1 to m^-2 s^-1 end +""" + P3_deposition_N_i(ip, T) + + - `ip` - a struct with ice nucleation parameters, + - `T` - air temperature [K]. + +Returns the number of ice nucleated via deposition nucleation with units of m^-3. +From Thompson et al 2004 eqn 2 as used in Morrison & Milbrandt 2015. +""" +function P3_deposition_N_i(ip::CMP.MorrisonMilbrandt2014, T::FT) where {FT} + + T₀ = ip.T₀ # 0°C + T_thres = ip.T_dep_thres # cutoff temperature + + Nᵢ = ifelse( + T < T_thres, + max(FT(0), FT(ip.c₁ * exp(ip.c₂ * (T₀ - T_thres)))), + max(FT(0), FT(ip.c₁ * exp(ip.c₂ * (T₀ - T)))), + ) + return Nᵢ * 1e3 # converts L^-1 to m^-3 +end + +""" + P3_het_N_i(ip, T, N_l, B, V_l, a, Δt) + + - `ip` - a struct with ice nucleation parameters, + - `T` - air temperature [K], + - `N_l` - number of droplets [m^-3], + - `B` - water-type dependent parameter [cm^-3 s^-1], + - `V_l` - volume of droplets to be frozen [m^3], + - `a` - empirical parameter [C^-1], + - `Δt` - timestep. + +Returns the number of ice nucleated within Δt via heterogeneous freezing +with units of m^-3. From Pruppacher & Klett 1997 eqn (9-51) as used in +Morrison & Milbrandt 2015. +""" +function P3_het_N_i( + ip::CMP.MorrisonMilbrandt2014, + T::FT, + N_l::FT, + V_l::FT, + Δt::FT, +) where {FT} + + a = ip.het_a # (celcius)^-1 + B = ip.het_B # cm^-3 s^-1 for rain water + V_l_converted = V_l * 1e6 # converted from m^3 to cm^3 + Tₛ = ip.T₀ - T + + return N_l * (1 - exp(-B * V_l_converted * Δt * exp(a * Tₛ))) +end + """ INP_concentration_frequency(INPC,T) diff --git a/src/parameters/IceNucleation.jl b/src/parameters/IceNucleation.jl index 801918c2f..02a97572a 100644 --- a/src/parameters/IceNucleation.jl +++ b/src/parameters/IceNucleation.jl @@ -66,16 +66,55 @@ function Koop2000(td::CP.AbstractTOMLDict) end """ - IceNucleationParameters{FT, DEP, HOM} + MorrisonMilbrandt2014{FT} + +Parameters for ice nucleation from Morrison & Milbrandt 2014 +DOI: 10.1175/JAS-D-14-0065.1 + +# Fields +$(DocStringExtensions.FIELDS) +""" +Base.@kwdef struct MorrisonMilbrandt2014{FT} <: ParametersType{FT} + "Cutoff temperature for deposition nucleation [K]" + T_dep_thres::FT + "coefficient [-]" + c₁::FT + "coefficient [-]" + c₂::FT + "T₀" + T₀::FT + "heterogeneous freezing parameter a [°C^-1]" + het_a::FT + "heterogeneous freezing parameter B [cm^-3 s^-1]" + het_B::FT +end + +function MorrisonMilbrandt2014(td::CP.AbstractTOMLDict) + name_map = (; + :temperature_homogenous_nucleation => :T_dep_thres, + :Thompson2004_c1_Cooper => :c₁, + :Thompson2004_c2_Cooper => :c₂, + :temperature_water_freeze => :T₀, + :BarklieGokhale1959_a_parameter => :het_a, + :BarklieGokhale1959_B_parameter => :het_B, + ) + parameters = CP.get_parameter_values(td, name_map, "CloudMicrophysics") + FT = CP.float_type(td) + return MorrisonMilbrandt2014{FT}(; parameters...) +end + +""" + IceNucleationParameters{FT, DEP, HOM, P3_type} Parameters for ice nucleation # Fields $(DocStringExtensions.FIELDS) """ -struct IceNucleationParameters{FT, DEP, HOM} <: ParametersType{FT} +struct IceNucleationParameters{FT, DEP, HOM, P3_type} <: ParametersType{FT} deposition::DEP homogeneous::HOM + p3::P3_type end IceNucleationParameters(::Type{FT}) where {FT <: AbstractFloat} = @@ -84,10 +123,16 @@ IceNucleationParameters(::Type{FT}) where {FT <: AbstractFloat} = function IceNucleationParameters(toml_dict::CP.AbstractTOMLDict) deposition = Mohler2006(toml_dict) homogeneous = Koop2000(toml_dict) + p3 = MorrisonMilbrandt2014(toml_dict) DEP = typeof(deposition) HOM = typeof(homogeneous) + P3_type = typeof(p3) FT = CP.float_type(toml_dict) - return IceNucleationParameters{FT, DEP, HOM}(deposition, homogeneous) + return IceNucleationParameters{FT, DEP, HOM, P3_type}( + deposition, + homogeneous, + p3, + ) end diff --git a/test/gpu_tests.jl b/test/gpu_tests.jl index 0330d0518..acb3a3313 100644 --- a/test/gpu_tests.jl +++ b/test/gpu_tests.jl @@ -568,6 +568,35 @@ end end end +@kernel function IceNucleation_P3_deposition_N_i_kernel!( + output::AbstractArray{FT}, + ip, + T, +) where {FT} + + i = @index(Group, Linear) + + @inbounds begin + output[1] = CMI_het.P3_deposition_N_i(ip, T[1]) + end +end + +@kernel function IceNucleation_P3_het_N_i_kernel!( + output::AbstractArray{FT}, + ip, + T, + N_l, + V_l, + Δt, +) where {FT} + + i = @index(Group, Linear) + + @inbounds begin + output[1] = CMI_het.P3_het_N_i(ip, T[1], N_l[1], V_l[1], Δt[1]) + end +end + @kernel function IceNucleation_INPC_frequency_kernel!( output::AbstractArray{FT}, ip_frostenberg, @@ -989,6 +1018,32 @@ function test_gpu(FT) dims = (1, 1) (; output, ndrange) = setup_output(dims, FT) + ip_P3 = ip.p3 + T = ArrayType([FT(240)]) + + kernel! = IceNucleation_P3_deposition_N_i_kernel!(backend, work_groups) + kernel!(output, ip_P3, T; ndrange) + + # test if P3_deposition_N_i is callable and returns reasonable values + @test Array(output)[1] ≈ FT(119018.93920746) + + dims = (1, 1) + (; output, ndrange) = setup_output(dims, FT) + + T = ArrayType([FT(240)]) + N_l = ArrayType([FT(2000)]) + V_l = ArrayType([FT(3e-18)]) + Δt = ArrayType([FT(0.1)]) + + kernel! = IceNucleation_P3_het_N_i_kernel!(backend, work_groups) + kernel!(output, ip_P3, T, N_l, V_l, Δt; ndrange) + + # test if P3_het_N_i is callable and returns reasonable values + @test Array(output)[1] ≈ FT(0.0002736160475969029) + + dims = (1, 1) + (; output, ndrange) = setup_output(dims, FT) + T = FT(233) INPC = FT(220000) diff --git a/test/heterogeneous_ice_nucleation_tests.jl b/test/heterogeneous_ice_nucleation_tests.jl index aaf5d3f5d..cc411112d 100644 --- a/test/heterogeneous_ice_nucleation_tests.jl +++ b/test/heterogeneous_ice_nucleation_tests.jl @@ -146,6 +146,22 @@ function test_heterogeneous_ice_nucleation(FT) end end + TT.@testset "P3 Deposition Nᵢ" begin + + T_warm = FT(235) + T_cold = FT(234) + + T_too_cold = FT(232) + + # higher ice concentration at colder temperatures + TT.@test CMI_het.P3_deposition_N_i(ip.p3, T_cold) > + CMI_het.P3_deposition_N_i(ip.p3, T_warm) + + # if colder than threshold T, use threshold T + TT.@test CMI_het.P3_deposition_N_i(ip.p3, T_too_cold) == + CMI_het.P3_deposition_N_i(ip.p3, ip.p3.T_dep_thres) + end + TT.@testset "ABIFM J" begin T_warm_1 = FT(229.2) @@ -179,6 +195,20 @@ function test_heterogeneous_ice_nucleation(FT) end end + TT.@testset "P3 Heterogeneous Nᵢ" begin + + T_warm = FT(235) + T_cold = FT(234) + N_liq = FT(2e5) + r_l = FT(2e-5) + V_l = FT(4 / 3 * FT(π) * r_l^3) + Δt = FT(0.1) + + # higher ice concentration at colder temperatures + TT.@test CMI_het.P3_het_N_i(ip.p3, T_cold, N_liq, V_l, Δt) > + CMI_het.P3_het_N_i(ip.p3, T_warm, N_liq, V_l, Δt) + end + TT.@testset "Frostenberg" begin temperatures = FT.([233, 257]) diff --git a/test/performance_tests.jl b/test/performance_tests.jl index 9f88a9015..c7fc5d5bb 100644 --- a/test/performance_tests.jl +++ b/test/performance_tests.jl @@ -68,7 +68,7 @@ function benchmark_test(FT) mixed_nuc = CMP.MixedNucleationParameters(FT) h2so4_nuc = CMP.H2S04NucleationParameters(FT) organ_nuc = CMP.OrganicNucleationParameters(FT) - # air and thermodunamics parameters + # air and thermodynamics parameters aps = CMP.AirProperties(FT) tps = TD.Parameters.ThermodynamicsParameters(FT) @@ -113,6 +113,9 @@ function benchmark_test(FT) x_sulph = FT(0.1) Delta_a_w = FT(0.27) + r_liq = FT(1e-6) + V_liq = FT(4 / 3 * π * r_liq^3) + Δt = FT(25) INPC = FT(1e5) @@ -148,11 +151,13 @@ function benchmark_test(FT) 80, ) bench_press(CMI_het.deposition_J, (kaolinite, Delta_a_w), 230) + bench_press(CMI_het.P3_deposition_N_i, (ip.p3, T_air_cold), 230) bench_press(CMI_het.ABIFM_J, (desert_dust, Delta_a_w), 230) + bench_press(CMI_het.P3_het_N_i, (ip.p3, T_air_cold, N_liq, V_liq, Δt), 230) bench_press( CMI_het.INP_concentration_frequency, (ip_frostenberg, INPC, T_air_cold), - 110, + 150, ) bench_press(CMI_hom.homogeneous_J, (ip.homogeneous, Delta_a_w), 230)