very much work in progress

Author: trontrytel

src/P3Scheme.jl

@@ -467,20 +467,26 @@ function terminal_velocity(
     bₛ(b) = b + μ
     cₛ(c) = c + λ
 
+    aₛ_m(a) = aₛ(a) * FT(π) / 6 * p3.ρ_i
+    bₛ_m(b) = bₛ(b) + 3
 
-    h₁(a) = a * N_0 * (16 * p3.ρ_i^2 * p3.γ^3 / (9 * FT(π) * α_va^2))^κ
-    h₂(b) = b + μ + κ * (3 * p3.σ - 2 * p3.β_va)
-    h₃(c) = c + λ
-    h₄(a) = a * N_0 * (p3.ρ_i / th.ρ_g)^(2 * κ)
+    spheres(a, b, c) = (aₛ(a),  bₛ(b),  cₛ(c))
+    spheres_m(a, b, c) = (aₛ_m(a), bₛ_m(b), cₛ(c))
 
-    spheres_n(a, b, c) = (aₛ(a), bₛ(b), cₛ(c))
-    spheres_m(a, b, c) = (a * N_0 * FT(π) / 6 * p3.ρ_i, b + μ + 3, h₃(c))
+    aₙₛ(a) = aₛ(a) * (16 * p3.ρ_i^2 * p3.γ^3 / (9 * FT(π) * α_va^2))^κ
+    bₙₛ(a) = bₛ(b) + κ * (3 * p3.σ - 2 * p3.β_va)
 
-    non_spheres_n(a, b, c) = (h₁(a),        h₂(b),           h₃(c))
-    non_spheres_m(a, b, c) = (h₁(a) * α_va, h₂(b) + p3.β_va, h₃(c))
+    aₙₛ_m = aₙₛ(a) * α_va
+    bₙₛ_m = bₙₛ(b) + p3.β_va
 
-    rimed_spheres_n(a, b, c) = (h₄(a),                       h₂(b),            h₃(c))
-    rimed_spheres_m(a, b, c) = (h₄(a) * FT(π) / 6 * th.ρ_g , h₂(b) + p3.β_va , h₃(c))
+    non_spheres(a, b, c) = (aₙₛ(a), bₙₛ(b), cₛ(c))
+    non_spheres_m(a, b, c) = (aₙₛ_m(a), bₙₛ_m(b), cₛ(c))
+
+    aᵣₛ(a) = aₛ(a) * (p3.ρ_i / th.ρ_g)^(2 * κ)
+    aᵣₛ_m(a) = aᵣₛ(a) * FT(π) / 6 * th.ρ_g
+
+    rimed_spheres_n(a, b, c) = (aᵣₛ(a), bₙₛ(b), cₛ(c))
+    rimed_spheres_m(a, b, c) = (aᵣₛ_m(a), bₙₛ_m(b), cₛ(c))
 
     #TODO - double check there is some equation typo here
     v_n_D_cr(D, a, b, c) =
@@ -494,11 +500,11 @@ function terminal_velocity(
         if F_r == 0
             # S mall particles
             # TODO - figure out the splatting vs indexes 1,2
-            v_m += ∫ᵧ(FT(0), D_th, spheres_m(aₛ[i], bₛ[i], cₛ[i])...)
-            v_n += integrate_to_gamma(FT(0), D_th, spheres_n(aₛ[i], bₛ[i], cₛ[i])...)
+            v_m += ∫_Γ(FT(0), D_th, spheres_m(aₛ[i], bₛ[i], cₛ[i])...)
+            v_n += ∫_Γ(FT(0), D_th, spheres(aₛ[i], bₛ[i], cₛ[i])...)
 
-            v_m += integrate_to_gamma(D_th, cutoff, non_spheres_m(aₛ[i], bₛ[i], cₛ[i])...)
-            v_n += integrate_to_gamma(D_th, cutoff, non_spheres_n(aₛ[i], bₛ[i], cₛ[i])...)
+            v_m += ∫_Γ(D_th, cutoff, non_spheres_m(aₛ[i], bₛ[i], cₛ[i])...)
+            v_n += ∫_Γ(D_th, cutoff, non_spheres_n(aₛ[i], bₛ[i], cₛ[i])...)
 
             # Get velocity coefficients for large particles
             v_m += integrate_to_gamma(cutoff, Inf, non_spheres_m(aₗ[i], bₗ[i], cₗ[i])...)