Skip to content

Commit 03f5828

Browse files
haakon-ehaakon-e
committed
feat: add p3 collisions
1 parent b311d9e commit 03f5828

File tree

7 files changed

+740
-6
lines changed

7 files changed

+740
-6
lines changed

‎Project.toml

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -12,6 +12,7 @@ LogExpFunctions = "2ab3a3ac-af41-5b50-aa03-7779005ae688"
1212
QuadGK = "1fd47b50-473d-5c70-9696-f719f8f3bcdc"
1313
RootSolvers = "7181ea78-2dcb-4de3-ab41-2b8ab5a31e74"
1414
SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b"
15+
StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
1516
Thermodynamics = "b60c26fb-14c3-4610-9d3e-2d17fe7ff00c"
1617

1718
[weakdeps]
@@ -32,5 +33,6 @@ MLJ = "0.20"
3233
QuadGK = "2.9"
3334
RootSolvers = "0.3, 0.4"
3435
SpecialFunctions = "1, 2"
36+
StaticArrays = "1.9"
3537
Thermodynamics = "0.12.4"
3638
julia = "1.6"

‎docs/make.jl

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -71,7 +71,7 @@ pages = Any[
7171
"References" => "References.md",
7272
]
7373

74-
mathengine = MathJax(
74+
mathengine = MathJax3(
7575
Dict(
7676
:TeX => Dict(
7777
:equationNumbers => Dict(:autoNumber => "AMS"),

‎docs/src/API.md

Lines changed: 6 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -182,10 +182,16 @@ P3Scheme.ice_melt
182182
```
183183

184184
#### Collisions with liquid droplets
185+
```@docs
186+
P3Scheme.bulk_liquid_ice_collision_sources
187+
```
185188
Supporting methods:
186189
```@docs
190+
P3Scheme.volumetric_collision_rate_integrand
187191
P3Scheme.compute_max_freeze_rate
188192
P3Scheme.compute_local_rime_density
193+
P3Scheme.get_liquid_integrals
194+
P3Scheme.∫liquid_ice_collisions
189195
```
190196

191197
### Supporting integral methods

‎docs/src/P3Scheme.md

Lines changed: 300 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -429,10 +429,306 @@ When modifying process rates, we now need to consider whether they are concerned
429429
in addition to whether they become sources and sinks of different prognostic variables with the inclusion of ``F_{liq}``.
430430
With the addition of liquid fraction, too, come new process rates.
431431

432-
!!! note
433-
TODO - Implement process rates, complete docs.
432+
## Microphysical Process Rates
434433

435-
## Heterogeneous Freezing
434+
The sources and sinks of the prognostic variables are given by the following equations:
435+
436+
!!! details "Symbol definitions"
437+
Table and text in this box is copied from [MorrisonMilbrandt2015](@cite), c.f. Appendix B.
438+
439+
| Symbol | Description |
440+
|:-------|:------------|
441+
| $\textcolor{cyan}{\textsf{WRM}}$ | Warm rain microphysics |
442+
| $\textcolor{magenta}{\textsf{NUC}}$ | Ice nucleation |
443+
| $\textcolor{lime}{\textsf{DEP}}$ | Ice deposition |
444+
| $\textcolor{lime}{\textsf{SUB}}$ | Ice sublimation |
445+
| $\textcolor{brown}{\textsf{COL}}$ | Collision/collection between liquid and ice |
446+
| $\textcolor{brown}{\textsf{SHD}}$ | Drop shedding |
447+
| $\textcolor{orange}{\textsf{HET}}$ | Heterogeneous freezing of cloud droplets and rain |
448+
| $\textcolor{orange}{\textsf{HOM}}$ | Homogeneous freezing of cloud droplets and rain |
449+
| $\textcolor{pink}{\textsf{MLT}}$ | Melting of ice |
450+
| $\textcolor{brown}{\textsf{WET}}$ | Particle densification due to wet growth in subfreezing temperatures |
451+
452+
The naming convention for the process rates below is as follows:
453+
- The first letter describes whether the process involves a change in mass (Q), number (N), or volume (B) mixing ratio.
454+
- we currently formulate the scheme in terms of mass content [kg/m³] instead of specific mass [kg/kg], so rates for $L_{rim}$ and $L_{ice}$ are given in [kg/m³/s] instead of [kg/kg/s]
455+
- For source and sink processes that do not involve multispecies interaction, the second letter indicates the species as cloud water (C), rain (R), or ice (I).
456+
- For source/sink processes that involve multispecies interaction (i.e., the same process acting as a source for one species but a sink for another), the second letter (C, R, or I) indicates the species that is reduced as a result of the process.
457+
- The remaining three letters indicate the type of microphysical process as defined in the table above.
458+
459+
460+
Liquid phase:
461+
```math
462+
\begin{align*}
463+
S_{q_c} &= \textcolor{cyan}{\text{QCWRM}}
464+
\textcolor{orange}{ - \text{QCHET} - \text{QCHOM}}
465+
\textcolor{brown}{ - \text{QCCOL}}
466+
\\
467+
S_{q_r} &= \textcolor{cyan}{\text{QRWRM}}
468+
\textcolor{orange}{ - \text{QRHET} - \text{QRHOM}}
469+
\textcolor{pink}{ + \text{QIMLT}}
470+
\textcolor{brown}{ - \text{QRCOL} + \text{QCSHD}}
471+
\\
472+
S_{N_c} &= \textcolor{cyan}{\text{NCWRM}}
473+
\textcolor{orange}{ - \text{NCHET} - \text{NCHOM}}
474+
\textcolor{brown}{ - \text{NCCOL}}
475+
\\
476+
S_{N_r} &= \textcolor{cyan}{\text{NRWRM}}
477+
\textcolor{orange}{ - \text{NRHET} - \text{NRHOM}}
478+
\textcolor{pink}{ + \text{NIMLT}}
479+
\textcolor{brown}{ - \text{NRCOL} + \text{NRSHD}}
480+
\\
481+
\end{align*}
482+
```
483+
Ice phase:
484+
```math
485+
\begin{align*}
486+
S_{L_{rim}} &=
487+
\textcolor{orange}{\text{QCHET} + \text{QRHET} + \text{QCHOM} + \text{QRHOM}}
488+
+ F_{rim}(\textcolor{lime}{ - \text{QISUB}}
489+
\textcolor{pink}{ - \text{QIMLT}})
490+
\\
491+
&\quad
492+
\textcolor{brown}{ + \text{QCCOL} + \text{QRCOL} + \text{QIWET}}
493+
\\
494+
S_{L_{ice}} &=
495+
\textcolor{orange}{\text{QCHET} + \text{QRHET} + \text{QCHOM} + \text{QRHOM}}
496+
\textcolor{lime}{ - \text{QISUB}}
497+
\textcolor{pink}{ - \text{QIMLT}}
498+
\\
499+
&\quad
500+
\textcolor{brown}{ + \text{QCCOL} + \text{QRCOL} + \text{QIWET}}
501+
\textcolor{lime}{ + \text{QIDEP}}
502+
\textcolor{magenta}{ + \text{QINUC}}
503+
\\
504+
S_{N_{ice}} &=
505+
\textcolor{orange}{ \text{NCHET} + \text{NRHET} + \text{NCHOM} + \text{NRHOM}}
506+
\textcolor{lime}{ - \text{NISUB}}
507+
\textcolor{magenta}{+ \text{NINUC} }
508+
\\
509+
S_{B_{rim}} &=
510+
\frac{1}{ρ^*} (
511+
\textcolor{orange}{\text{QCHET} + \text{QRHET} + \text{QCHOM} + \text{QRHOM}})
512+
+ \frac{1}{ρ_{rim}} F_{rim}(
513+
\textcolor{lime}{ - \text{QISUB}}
514+
\textcolor{pink}{ - \text{QIMLT}})
515+
\\
516+
&\quad
517+
\textcolor{brown}{+ \text{BCCOL} + \text{BRCOL} + \text{BIWET}}
518+
\\
519+
\end{align*}
520+
```
521+
522+
Terms that are related are colored the same:
523+
524+
- The $\textcolor{cyan}{\textsf{cyan}}$ colors → warm rain microphysics (i.e. WRM)
525+
- nucleation, condensation, autoconversion, accretion, evaporation
526+
- The $\textcolor{orange}{\textsf{orange}}$ colors → heterogeneous and homogeneous freezing (i.e. HET, HOM)
527+
- The $\textcolor{brown}{\textsf{brown}}$ colors → liquid-ice collisions (i.e. COL, SHD, WET)
528+
- The $\textcolor{pink}{\textsf{pink}}$ colors → melting (i.e. MLT)
529+
- The $\textcolor{lime}{\textsf{lime}}$ colors → deposition/sublimation (i.e. DEP, SUB)
530+
- The $\textcolor{magenta}{\textsf{magenta}}$ colors → nucleation (i.e. NUC)
531+
532+
!!! note "Differences with Morrison & Milbrandt (2015)"
533+
There are numerous differences between the P3 scheme in Morrison & Milbrandt (2015) and the P3 scheme in this package:
534+
- the maximum freezing rate (wet growth limit) is computed at the particle level, instead of as a bulk process
535+
- changes in rime volume due to (dry) collisions with cloud and rain are computed at the particle level, with a local rime density evaluated as a function of both the liquid and ice particle diameters
536+
- in the wet growth regime, the densification process is computed as a rapidly adjusting bulk process, instead of instantenous particle-level densification. We also scale by the fraction of the mass rate that undergoes wet growth.
537+
538+
539+
Below, we describe the different processes in more detail.
540+
541+
### Collisions with liquid droplets
542+
543+
Collisions between liquid droplets (cloud or rain) and ice particles are parameterized
544+
by integrating the volumetric collision rate $∂_t\mathcal{V}$ [$\text{m}^3/\text{s}$]
545+
over the particle size distributions. The volumetric collision rate is on the form
546+
```math
547+
∂_t\mathcal{V}_l(D_i, D_l) = E(D_i, D_l) K(D_i, D_l) |v_i(D_i) - v_l(D_l)|
548+
```
549+
where $l ∈ \{c, r\}$ denotes cloud or rain,
550+
- ``E(D_i, D_l)`` is the collision efficiency (which we assume to be 1),
551+
- ``K(D_i, D_l) = π (r_i + D_l/2)^2`` is the collision cross section,
552+
- ``r_i`` is the effective radius of the ice particle, which is the radius of a circle with the same cross-sectional area as the (in general, non-spherical) ice particle,
553+
- ``v_i(D_i)`` and ``v_l(D_l)`` are the terminal velocities of the ice particle and the liquid droplet, respectively.
554+
555+
!!! note "Number of collisions per unit time"
556+
If we multiply ``∂_t\mathcal{V}_l(D_i, D_l)`` by the liquid droplet size distribution ``N_l(D_l)``, we obtain the kernel
557+
```math
558+
∂_t\mathcal{V}_l(D_i, D_l) N_l(D_l),
559+
```
560+
which represents the number of liquid droplets of size ``D_l`` colliding with an ice particle of size ``D_i`` per unit time. If we know, say, the mass of each liquid droplet, and integrate the mass-kernel over all liquid drop sizes, we get the total mass of the collected liquid droplets by an ice particle of size ``D_i``.
561+
562+
The quantity $∂_t\mathcal{N}_{\text{col},l}(D_i)$ [$\text{s}^{-1}$],
563+
```math
564+
∂_t\mathcal{N}_{\text{col},l}(D_i) = ∫_0^∞ ∂_t\mathcal{V}_l(D_i, D_l) N_l(D_l) \mathrm{d}D_l
565+
```
566+
quantifies how many liquid droplets collide with an ice particle of size $D_i$ per unit time.
567+
568+
!!! note "Bulk loss of liquid droplets"
569+
The bulk loss of liquid droplets is then given by
570+
```math
571+
\textcolor{brown}{\text{NCCOL}} = ∫_0^∞ ∂_t\mathcal{N}_{\text{col},l}(D_i) N'_i(D_i) \mathrm{d}D_i,
572+
\quad
573+
\textcolor{brown}{\text{NRCOL}} = ∫_0^∞ ∂_t\mathcal{N}_{\text{col},r}(D_i) N'_i(D_i) \mathrm{d}D_i,
574+
```
575+
576+
The total mass of the collected liquid droplets $∂_t\mathcal{M}_{col,l}(D_i)$ [$\text{kg/s}$] can be expressed as
577+
```math
578+
∂_t\mathcal{M}_{\text{col},l}(D_i)
579+
= ∫_0^∞ ∂_t\mathcal{V}_l(D_i, D_l) N_l(D_l) m_l(D_l) \mathrm{d}D_l
580+
```
581+
582+
When liquid droplets collide with ice particles, we assume that they can either freeze or be shed as
583+
rain particles. Above freezing, all particles are shed. In subfreezing temperatures, normally all particles are frozen. However,
584+
close to the freezing temperature, there may not be enough time for all particles to freeze.
585+
We refer to this as the wet growth regime. In this case, the particles that are not frozen are shed as rain.
586+
The wet growth regime is determined by comparing the freezing rate to the collection rate.
587+
588+
The maximum freezing rate is based on Musil (1970) [Musil1970](@cite),
589+
which considers the heat transfer rate for a spherical wet hailstone.
590+
The maximum freezing rate [kg/s] is computed as
591+
```math
592+
∂_t\mathcal{M}_\text{max}(D_i)
593+
= 2π D_i F_v(D_i) \frac{- K_t ΔT + H_v D_v Δq_{v,\text{sat}}}{L_f + C_p ΔT}
594+
```
595+
where the first term in the numerator is the heat transfer rate to the air surrounding the hailstone,
596+
and the second term is evaporative cooling. The denominator is the heat that must be dissipated by the hailstone.
597+
598+
!!! details "Symbol definitions"
599+
| Symbol | Units | Description |
600+
|:-----------|:------|:------------|
601+
| $D_i$ | $\text{m}$ | diameter of the ice particle |
602+
| $F_v(D_i)$ | $-$ | ventilation factor |
603+
| $K_t$ | $\text{W/(m K)}$ | thermal conductivity |
604+
| $ΔT$ | $\text{K}$ | temperature difference between the surface of the ice particle and the surrounding air |
605+
| $H_v$ | $\text{J/kg}$ | latent heat of vaporization |
606+
| $D_v$ | $\text{m}^2/\text{s}$ | vapor diffusivity |
607+
| $Δq_{v,\text{sat}}$ | $\text{kg/kg}$ | difference in saturation vapor mixing ratio between the surface of the ice particle and the surrounding air |
608+
| $L_f$ | $\text{J/kg}$ | latent heat of fusion |
609+
| $C_p$ | $\text{J/(kg K)}$ | specific heat capacity of water |
610+
611+
For collisions with ice particles of size $D_i$, the freezing (riming) rate $∂_t\mathcal{M}_\text{frz}$ [$\text{kg/s}$] is thus
612+
```math
613+
∂_t\mathcal{M}_\text{frz} = \min \left( ∂_t\mathcal{M}_\text{col}, ∂_t\mathcal{M}_\text{max} \right),
614+
\quad \textsf{where} \quad
615+
∂_t\mathcal{M}_\text{col} = ∂_t\mathcal{M}_\text{col,c} + ∂_t\mathcal{M}_\text{col,r},
616+
```
617+
where $∂_t\mathcal{M}_{col}$ is the total mass of the collected liquid droplets. The fraction of the collected liquid that freezes is given by
618+
```math
619+
f_\text{frz}(D_i)
620+
= \frac{∂_t \mathcal{M}_\text{frz}(D_i)}{∂_t \mathcal{M}_\text{col}(D_i)}
621+
```
622+
623+
!!! note "Bulk freezing rates"
624+
The bulk freezing rates are obtained by integrating the per-particle freezing rates over the ice particle size distribution,
625+
```math
626+
\textcolor{brown}{\text{QCCOL}}
627+
= ∫_0^∞ ∂_t\mathcal{M}_{\text{col},c}(D_i) f_\text{frz}(D_i) N'_i(D_i) \mathrm{d}D_i,
628+
\quad
629+
\textcolor{brown}{\text{QRCOL}}
630+
= ∫_0^∞ ∂_t\mathcal{M}_{\text{col},r}(D_i) f_\text{frz}(D_i) N'_i(D_i) \mathrm{d}D_i,
631+
```
632+
which can be combined to give the bulk liquid freezing rate
633+
```math
634+
\textcolor{brown}{\text{QCCOL} + \text{QRCOL}}
635+
= ∫_0^∞ ∂_t \mathcal{M}_\text{frz}(D_i) N'_i(D_i) \mathrm{d}D_i
636+
= ∫_0^∞ ∂_t \mathcal{M}_\text{col}(D_i) f_\text{frz}(D_i) N'_i(D_i) \mathrm{d}D_i
637+
```
638+
639+
Any excess mass $∂_t\mathcal{M}_{shd}$ [$\text{kg/s}$] is shed as rain, where
640+
```math
641+
∂_t\mathcal{M}_\text{shd} = ∂_t\mathcal{M}_{col} - ∂_t\mathcal{M}_\text{frz}
642+
```
643+
We assume that all shed droplets are shed at some fixed diameter $D_\text{shd}$.
644+
Following [MorrisonMilbrandt2015](@cite), we set $D_\text{shd} = 1 \text{ mm}$.
645+
This implies that the number of shed droplets is
646+
```math
647+
∂_t\mathcal{N}_\text{shd} = \frac{∂_t\mathcal{M}_\text{shd}}{m(D_\text{shd})}
648+
```
649+
650+
!!! note "Bulk source to rain mass and number"
651+
The corresponding bulk source to rain mass and number is
652+
```math
653+
\begin{align*}
654+
\textcolor{brown}{\text{QRSHD}}
655+
&= ∫_0^∞ ∂_t\mathcal{M}_\text{shd}(D_i) N'_i(D_i) \mathrm{d}D_i
656+
= ∫_0^∞ (1 - f_{\text{frz}}(D_i)) ∂_t\mathcal{M}_\text{col}(D_i) N'_i(D_i) \mathrm{d}D_i \\
657+
&= ∫_0^∞ ∂_t\mathcal{M}_\text{col}(D_i) N'_i(D_i) \mathrm{d}D_i - (\textcolor{brown}{\text{QCCOL} + \text{QRCOL}})
658+
\\
659+
\textcolor{brown}{\text{NRSHD}} &= ∫_0^∞ ∂_t\mathcal{N}_\text{shd} N'_i(D_i) \mathrm{d}D_i
660+
= ∫_0^∞ \frac{∂_t\mathcal{M}_\text{shd}}{m(D_\text{shd})} N'_i(D_i) \mathrm{d}D_i
661+
= \frac{\textcolor{brown}{\text{QRSHD}}}{m(D_\text{shd})},
662+
\end{align*}
663+
```
664+
665+
Finally, how does this affect the rime volume?
666+
667+
Collisions with rain, cloud, and shedded wet growth all contribute to the rime volume.
668+
669+
The change in rime volume is given by the change in mass divided by some rime density.
670+
671+
!!! note "Local rime density parameterization"
672+
Following [MorrisonMilbrandt2015](@cite), we apply the Cober & List (1993) [CoberList1993](@cite) parameterization to calculate the local rime density for each collision pair.
673+
For an ice particle of size $D_i$ colliding with a liquid droplet of size $D_l$ (where $l ∈ \{c, r\}$ for cloud or rain), the rime density is given by
674+
```math
675+
\begin{align}
676+
V_\text{term}(D_i, D_l) &= |v_i(D_i) - v_l(D_l)|,
677+
\\
678+
R_i(D_i, D_l) &= \frac{D_l \, V_\text{term}(D_i, D_l)}{2 |T - T_\text{freeze}|},
679+
\\
680+
ρ'_{rim}(R_i) &=
681+
\begin{cases}
682+
a_\text{CL} + b_\text{CL} R_i + c_\text{CL} R_i^2 & \text{if }\, 1 ≤ R_i ≤ 8,
683+
\\
684+
(1 - \frac{R_i-8}{12-8}) ρ'_{rim}(8) + \frac{R_i-8}{12-8} ρ^* & \text{if }\, 8 < R_i ≤ 12,
685+
\end{cases}
686+
\end{align}
687+
```
688+
where $ρ^* = 900$ kg/m³ is the density of solid bulk ice, and $a_\text{CL}, b_\text{CL}, c_\text{CL}$ are the coefficients for the Cober & List (1993) parameterization [CoberList1993](@cite).
689+
The $R_i$ quantity is limited to the range $1 ≤ R_i ≤ 12$.
690+
Their values are $a_\text{CL} = 51$, $b_\text{CL} = 114$, and $c_\text{CL} = -5.5$.
691+
692+
The change in rime volume at some diameter $D_i$ is then given by
693+
```math
694+
∂_t\mathcal{B}_{\text{rim},l}(D_i)
695+
= ∫_0^∞ \frac{∂_t\mathcal{V}_l(D_i, D_l) m_l(D_l)}{ρ'_{rim}(R_i(D_i, D_l))} N_l(D_l) \mathrm{d}D_l
696+
```
697+
The bulk rate of rime volume change is limited by the fraction of mass that freezes, and is then given by
698+
```math
699+
\textcolor{brown}{\text{BCCOL} + \text{BRCOL}}
700+
= \sum_{l ∈ \{c, r\}} ∫_0^∞ ∂_t\mathcal{B}_{\text{rim},l}(D_i) f_\text{frz}(D_i) N'_i(D_i) \mathrm{d}D_i
701+
```
702+
703+
In the wet growth regime, the rime compacts, rapidly relaxing towards the solid bulk ice density $ρ^* = 900$ kg/m³.
704+
Because the P3 scheme only tracks the bulk rime volume $B_\text{rim}$, the densification rate is applied to the bulk rime volume.
705+
As a measure of the "intensity" of wet growth, we use the fraction of the total liquid collection that occurs by particles in wet growth regime,
706+
```math
707+
\begin{align*}
708+
f_\text{wet}
709+
&= \frac{
710+
∫_0^∞ \mathbb{1}_\text{wet}(D_i) ⋅ ∂_t \mathcal{M}_{col}(D_i) N'(D_i) \mathrm{d}D_i
711+
}{∫_0^∞ ∂_t \mathcal{M}_\text{col}(D_i) N'(D_i) \mathrm{d}D_i
712+
}, \\
713+
&= \frac{
714+
∫_0^∞ \mathbb{1}_\text{wet}(D_i) ⋅ ∂_t \mathcal{M}_{col}(D_i) N'(D_i) \mathrm{d}D_i
715+
}{\textcolor{brown}{\text{QCCOL} + \text{QRCOL} + \text{QRSHD}}
716+
},
717+
\end{align*}
718+
```
719+
The bulk rate of rime volume change is then given by
720+
```math
721+
\begin{align*}
722+
\textcolor{brown}{\text{QIWET}}
723+
&= f_\text{wet} ⋅ \frac{1}{τ_\text{wet}} (L_\text{ice} - L_\text{rim})
724+
= f_\text{wet} ⋅ \frac{1}{τ_\text{wet}} L_\text{ice} (1 - F_{rim})
725+
\\
726+
\textcolor{brown}{\text{BIWET}}
727+
&= f_\text{wet} ⋅ \frac{1}{τ_\text{wet}} \left(\frac{L_\text{ice}}{ρ^*} - B_\text{rim}\right) \\
728+
\end{align*}
729+
```
730+
731+
### Heterogeneous Freezing
436732

437733
Immersion freezing is parameterized based on water activity and follows the ABIFM
438734
parameterization from [KnopfAlpert2013](@cite).
@@ -460,7 +756,7 @@ include("plots/P3ImmersionFreezing.jl")
460756
```
461757
![](P3_het_ice_nucleation.svg)
462758

463-
## Melting
759+
### Melting
464760

465761
Melting rate is derived in the same way as in the
466762
[1-moment scheme](https://clima.github.io/CloudMicrophysics.jl/dev/Microphysics1M/#Snow-melt).

‎src/P3.jl

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -11,6 +11,7 @@ import SpecialFunctions as SF
1111
import QuadGK as QGK
1212
import RootSolvers as RS
1313
import LogExpFunctions
14+
import StaticArrays as SA
1415

1516
import ClimaParams as CP
1617
import CloudMicrophysics.Parameters as CMP

0 commit comments

Comments
 (0)