more tests

CliMA · Feb 6, 2024 · 520be90 · 520be90
1 parent cbbb4a1
commit 520be90
Show file tree

Hide file tree

Showing 2 changed files with 23 additions and 13 deletions.
diff --git a/src/P3Scheme.jl b/src/P3Scheme.jl
@@ -228,7 +228,7 @@ Returns the shape parameter μ for a given λ value
 Eq. 3 in Morrison and Milbrandt (2015).
 """
 function DSD_μ(p3::PSP3, λ::FT) where {FT}
-    # @assert λ > FT(0)
+    #@assert λ > FT(0)
     return min(p3.μ_max, max(FT(0), p3.a * λ^p3.b - p3.c))
 end
 
@@ -306,7 +306,8 @@ function q_gamma(
 
     D_th = D_th_helper(p3)
     λ = exp(log_λ)
-    # println(" log λ = ", log_λ, ", λ = ", λ)
+    #= println(" log λ = ", log_λ, ", λ = ", λ)
+    println("F_r = ", F_r) =#
     N_0 = DSD_N₀(p3, N, λ)
 
     return ifelse(
@@ -343,6 +344,9 @@ function distribution_parameter_solver(
 
     N̂ = FT(1e20)
 
+    min = 2.0 * 1e4 
+    max = 2.5 * 1e4
+
     # Get the thresholds for different particles regimes
     th = thresholds(p3, ρ_r, F_r)
 
@@ -353,7 +357,7 @@ function distribution_parameter_solver(
     sol =
         RS.find_zero(
             shape_problem,
-            RS.SecantMethod(FT(log(4 * 1e4)), FT(log(5 * 1e4))),
+            RS.SecantMethod(FT(log(min)), FT(log(max))),
             RS.CompactSolution(),
             RS.RelativeSolutionTolerance(eps(FT)),
             10,

diff --git a/test/p3_LambdaTesting.jl b/test/p3_LambdaTesting.jl
@@ -265,8 +265,10 @@ function buckets(λ_min::FT, λ_max::FT, F_r_min::FT, F_r_max::FT, p3::PSP3, ρ_
                     E[i, j] = FT(1)
                 elseif er <= FT(1) 
                     E[i, j] = FT(2)
-                else 
+                elseif er <= FT(1e10) 
                     E[i, j] = FT(3)
+                else 
+                    E[i, j] = FT(4)
                 end
             else 
                 if er <= min_error
@@ -275,8 +277,10 @@ function buckets(λ_min::FT, λ_max::FT, F_r_min::FT, F_r_max::FT, p3::PSP3, ρ_
                     E[i, j] = FT(1.5)
                 elseif er <= FT(1) 
                     E[i, j] = FT(2.5)
-                else 
+                elseif er <= FT(1e10) 
                     E[i, j] = FT(3.5)
+                else 
+                    E[i, j] = FT(4.5)
                 end
             end
         end
@@ -297,7 +301,7 @@ function plot_buckets(λ_min::FT, λ_max::FT, F_r_min::FT, F_r_max::FT, ρ_r_min
         (λs, F_rs, E) = buckets(λ_min, λ_max, F_r_min, F_r_max, p3, ρ_r, N, λSteps, F_rSteps)
 
         Plt.heatmap!(λs, F_rs, E)
-        Plt.Colorbar(f[i, 2], limits = (0, 3.5), colormap = :viridis, flipaxis = false)
+        Plt.Colorbar(f[i, 2], limits = (0, 4.5), colormap = :viridis, flipaxis = false)
 
         # Plt.Colorbar(f[1, 2], label = "error", limits = FT.(range_val))
     end
@@ -320,16 +324,18 @@ F_r_max = FT(1-eps(FT))
 
 #constant_Fr_ρr(FT(0.8), FT(400), FT(0), FT(2000))
 
-#plot_all_gaps(λ_min, λ_max, F_r_min, F_r_max, ρ_r_min, ρ_r_max, p3, 200, 200, 10)
 
 #plot_relerrors(λ_min, λ_max, p3, 200, 200, 1)
 
-F = FT(0.34) 
-ρ = FT(650)
-λ = FT(35000)
-N = FT(1e8) 
+F = FT(0.5) 
+ρ = FT(100)
+λ = FT(26002)
+N = FT(1e12) 
 
-#diff = λ_diff(F, ρ, N, λ, p3)
+diff = λ_diff(F, ρ, N, λ, p3)
+println(diff) 
 
 #println(λ_diff(FT(0.5), FT(400), FT(1e8), FT(15000), p3))
-plot_buckets(λ_min, λ_max, F_r_min, F_r_max, ρ_r_min, ρ_r_max, p3, 250, 250, 9)
+# plot_all_gaps(λ_min, λ_max, F_r_min, F_r_max, ρ_r_min, ρ_r_max, p3, 250, 250, 9)
+
+plot_buckets(λ_min, λ_max, F_r_min, F_r_max, ρ_r_min, ρ_r_max, p3, 500, 500, 9)