refactor(Common,P3): consolidate the Chen 2022 velocity curves#762
Draft
haakon-e wants to merge 18 commits into
Draft
refactor(Common,P3): consolidate the Chen 2022 velocity curves#762haakon-e wants to merge 18 commits into
haakon-e wants to merge 18 commits into
Conversation
Member
Author
|
This change is part of the following stack: Change managed by git-spice. |
haakon-e
changed the title
Codecov Report❌ Patch coverage is Additional details and impacted files@@ Coverage Diff @@
## he/ice-selfcol-triangle #762 +/- ##
===========================================================
+ Coverage 92.81% 92.83% +0.01%
===========================================================
Files 56 57 +1
Lines 2952 2957 +5
===========================================================
+ Hits 2740 2745 +5
Misses 212 212
🚀 New features to boost your workflow:
|
haakon-e
force-pushed
the
he/ice-selfcol-triangle
branch
from
July 10, 2026 23:15
8ef7fe6 to
bbb3da8
Compare
haakon-e
force-pushed
the
he/velocity-unify
branch
from
July 10, 2026 23:15
698c7ee to
b9b0cd9
Compare
…unds Store a size-distribution quadrature rule on `P3IceParams` (Gauss-Legendre by default; set the order via `quadrature_order`, or pass a full rule via `quad`). `quad` is required by `integrate` and the P3 process functions; `build_quadrature` is removed. Integration bounds place the Chen 2022 small/large-ice velocity breakpoint on a subinterval boundary (`velocity_integral_bounds`), so a fixed-order rule stays accurate. This applies to the velocity-weighted integrals and to `ice_melt`, whose ventilation factor carries the regime break. The weighted-velocity integrals run the integral unconditionally and select the degenerate-state result branchlessly; `compute_max_freeze_rate` selects its limit values branchlessly for the same reason. Ice self-collection integrates over the upper triangle `D₁ ≤ D₂`, which puts the `D₁ = D₂` derivative discontinuity of the relative fall speed on a subinterval boundary and halves the integrand evaluations.
…ion benchmark Un-gate the 2M+P3 GPU kernel test and the P3 collision benchmark from VERSION >= v"1.12", set their allocation and memory budgets to zero, and recalibrate three time budgets from observed CI timings.
Check the triangular ice self-collection integral against the full-square double integral with the pair-counting factor, replacing the reference TODO. Cover number_tendency_from_mass_limits with x_min = 0, where the lower-bound relaxation never fires.
Assert zero allocations, type stability, and a time budget on the production entry point, not only its collision sub-integral.
Place the known derivative discontinuities of the liquid-ice collision integrands on subinterval boundaries: the ice terminal-velocity regime break at D_cutoff in the outer (ice) integral, and the ice-liquid fall-speed crossing in the inner (liquid) integrals. The rain N and M components already split at the crossing through the closed form; the cloud path and the rain B_rim quadrature did not. The crossover diameter is computed once per outer node and shared between the closed-form N and M components and the B_rim bounds.
The upper integration bound is the (1 - p) quantile of the size distribution, so the last subinterval can span most of the log(1/p) decay lengths of the exponential tail; in large-mean-size states it covers several decades and dominates the low-order quadrature error (ice self-collection error at Gauss-Legendre order 12 is 2.4% in a hail-core state). A breakpoint at three decay lengths, clamped into the integration interval, keeps each subinterval resolvable at low order (the same state and order reach 7e-5). Update the collision smoke references: with the crossing and decay-scale breakpoints the quadrature resolves the previously unresolved integrand structure, which shifts the converged values. The wet-growth indicator component moves the most; its integrand is a step function in the ice diameter, and its quadrature error is dominated by the step at any order.
The freeze/shed partition min() and the wet-growth indicator change branch at the diameter where the collected liquid mass rate balances the Musil freeze limit; interior to a subinterval, the indicator step limits quadrature convergence of the shedding and wet-densification components at ~1e-2 relative error regardless of order. wet_growth_onset_diameter locates up to two balance points (the wet-growth window closes again at large sizes where the freeze limit outgrows collection) from a closed-form balance: the cloud collection term neglects the droplet fall speed, making it a polynomial moment of the cloud size distribution, and the rain term is the closed-form mass component. The onset placement is within 0.7% of the exact balance across tested regimes; the step components' error at order 12 drops from 1.5e-3 to 4e-5 (shedding) and from 1e-1 to 7e-3 (wet indicator).
With the crossing, decay-scale, and wet-growth-onset breakpoints, the three-level error study (17 regime states, GL(128) reference) splits the bulk-tendency error at Gauss-Legendre order 6 into two families. The transport components (sedimentation velocities, melt) converge to at most 0.45%. Every error above that sits in a component proportional to the collision efficiency (collision sources, ice self-collection, rime), with a worst case of 1.4%; these rates carry the E = 1 assumption, whose parameterization uncertainty exceeds the quadrature error by an order of magnitude. Order 5 exceeds 2% in the transport family and is rejected. The quadrature-order sweep locks orders 6 through 12 in as a regression test against a GL(200) reference, with graupel/hail-core states in the state set.
VolumetricCollisionRate carries the ice and liquid terminal-velocity closures as fields, so its consumers query them instead of receiving the velocities separately: the fall-speed crossing bounds, the closed-form rain integrals (whose velocity-curve coefficients now come from the closure rather than a second coefficient evaluation), and the wet-growth onset balance. Integrands without the velocity structure integrate over unmodified bounds. velocity_breakpoints(v_term) returns the diameters where a terminal-velocity closure changes functional form; velocity_integral_bounds takes the closure and no longer reaches into the Chen 2022 cutoff field. The wet-growth onset bisection runs through the RootSolvers API in log diameter. The freeze-limit and weighted-velocity denominators divide directly; the final selection discards the non-finite branch, so the intermediate clipping is unnecessary. The liquid collision integrand names its number, mass, and rime-volume rates. The documentation plots use the package-default quadrature rule.
The study behind the default quadrature order lives in test/p3_quadrature_error_study.jl: run_quadrature_error_study measures each quadrature-consuming P3 quantity against a high-order reference across regime states, attributed to transport, collision-efficiency, and bulk levels. The order-sweep regression test includes it for the shared column-state set, and the developers guide documents when to rerun it.
Replace ChebyshevGauss with Gauss-Legendre to match the scheme default. Scheme-path calls use the default GL(6); the standalone reference integrals use the smallest GL order meeting their existing tolerance. Raise the rain ND tolerance to 3e-6 to reflect the 2p truncation of the integration bounds. Drop the redundant Gauss-Legendre convergence testset already covered by the dedicated quadrature testset.
…ax-freeze clamp Add a docstring for `collision_cross_section_ice_ice` matching the sibling `collision_cross_section_ice_liquid` convention, and list the symbol in the P3Scheme API-docs block. The `compute_max_freeze_rate` closure returns `floatmax(FT)` when the Musil (1970) denominator is non-positive. ForwardDiff propagation was verified to be AD-safe: the `ifelse` discards the non-finite branch (finite derivative in Float32/Float64), `floatmax` and `zero` resolve for the `Dual` type, and both branches infer to the same concrete type. No numerical change.
haakon-e
force-pushed
the
he/ice-selfcol-triangle
branch
from
July 11, 2026 00:58
bbb3da8 to
8a2faf4
Compare
haakon-e
force-pushed
the
he/velocity-unify
branch
from
July 11, 2026 00:58
b9b0cd9 to
dfc5e53
Compare
Store terms::NTuple{N, NTuple{3, FT}} of (a, b, c) per summand of
v(D) = sum(a D^b exp(-c D)) instead of parallel ai/bi/ci tuples. Add a
Chen2022_vel_coeffs method on the curve returning the (ai, bi, ci)
tuples; the closed-form rain consumers read coefficients through it.
…nctor P3IceParticleVelocityFunctor holds the aspect-ratio treatment as a field, set from state.params.aspect_ratio at construction, instead of reading it through the state at call time.
Move Chen2022_vel_coeffs, Chen2022VelocityCurve, particle_terminal_velocity, and the Chen 2022 pdf helpers from Common.jl to TerminalVelocityCurves.jl, included inside the Common module. No code change.
Change the Chen 2022 terminal-velocity coefficient representation from
parallel tuples `(ai, bi, ci)` to per-term triples
`terms::NTuple{N, NTuple{3, FT}}` of `(aₖ, bₖ, cₖ)`, aligning the return
of `Chen2022_vel_coeffs` with the `terms` already stored on
`Chen2022VelocityCurve` and with its evaluation `∑ₖ aₖ D^bₖ exp(-cₖ D)`.
- `Chen2022_vel_coeffs(coeffs, ρₐ[, ρᵢ])` returns `terms`; the same unit
conversions are applied per term, so values are unchanged.
- `Chen2022VelocityCurve` builds `terms` directly; the parallel-array
constructor is replaced by one taking `terms` and promoting to a common
`FT`. `Chen2022_vel_coeffs(v::Chen2022VelocityCurve)` returns `v.terms`.
- `closed_rain_inner_NM` takes `terms` and loops over the triples instead
of indexing parallel arrays; call sites pass `terms` and drop the
`SVector` wrapping.
- Update the Microphysics1M/2M consumers to sum over `terms`, preserving
the exact summation order (per-term `ϕ^κ` factors kept inside the sum).
Representation change only: `Chen2022_vel_coeffs`, curve evaluation, and a
representative `closed_rain_inner_NM` call are byte-identical before and
after in Float64 and Float32.
haakon-e
force-pushed
the
he/velocity-unify
branch
from
July 11, 2026 02:32
dfc5e53 to
c89f4e8
Compare
haakon-e
force-pushed
the
he/ice-selfcol-triangle
branch
from
July 11, 2026 03:28
dc1db96 to
355b541
Compare
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
1 participant
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Purpose
group Chen2022VelocityCurve coefficients by term
Content
Store
terms::NTuple{N, NTuple{3, FT}}of(a, b, c)per summand ofv(D) = sum(a D^b exp(-c D))instead of parallelai/bi/cituples. Add aChen2022_vel_coeffsmethod on the curve returning the(ai, bi, ci)tuples; the closed-form rain consumers read coefficients through it.P3IceParticleVelocityFunctorholds the aspect-ratio treatment as a field, set fromstate.params.aspect_ratioat construction, instead of reading it through the state at call time.Move Chen2022_vel_coeffs, Chen2022VelocityCurve, particle_terminal_velocity, and the Chen 2022 pdf helpers from Common.jl to TerminalVelocityCurves.jl, included inside the Common module. No code change.