test surface_temperature

jl-pkgs · Oct 13, 2024 · 92ffd3c · 92ffd3c
1 parent 5b50d7b
commit 92ffd3c
Show file tree

Hide file tree

Showing 3 changed files with 56 additions and 6 deletions.
diff --git a/src/surface_temperature.jl b/src/surface_temperature.jl
@@ -1,5 +1,8 @@
+"""
+- cp: capacity_heat, [J kg-1 K-1]
+"""
 function surface_temperature_jl(T_air::FT, rh_air::FT, depth_snow::FT, depth_water::FT,
-  capacity_heat_soil1::FT, capacity_heat_soil0::FT, Gheat_g::FT,
+  cp_soil1::FT, cp_soil0::FT, Gheat_g::FT,
   depth_soil1::FT, ρ_snow::FT, tempL_u::FT, Rn_g::FT,
   E_soil::FT, E_water_g::FT, E_snow_g::FT, λ_soil1::FT,
   perc_snow_g::FT, G_soil1::FT,
@@ -30,7 +33,7 @@ function surface_temperature_jl(T_air::FT, rh_air::FT, depth_snow::FT, depth_wat
   G::FT = 0.0
 
   if depth_snow <= 0.02
-    ttt = capacity_heat_soil1 * 0.02 / length_step
+    ttt = cp_soil1 * 0.02 / length_step
     T_ground = (T_ground_last * ttt * ra_g * depth_soil1 + Gg * ra_g * depth_soil1 + ρₐ * cp * T_air * depth_soil1 + ra_g * λ_soil1 * T_soil1_last) /
                (ρₐ * cp * depth_soil1 + ra_g * λ_soil1 + ttt * ra_g * depth_soil1)
     T_ground = clamp(T_ground, T_ground_last - 25, T_ground_last + 25)
@@ -45,7 +48,7 @@ function surface_temperature_jl(T_air::FT, rh_air::FT, depth_snow::FT, depth_wat
     G = clamp(G, -100.0, 100.0)
 
   elseif depth_snow > 0.02 && depth_snow <= 0.05
-    ttt = capacity_heat_soil1 * 0.02 / length_step  # for soil fraction part
+    ttt = cp_soil1 * 0.02 / length_step  # for soil fraction part
 
     T_soil0 = (T_soil0_last * ttt * ra_g * depth_soil1 + Gg * ra_g * depth_soil1 + ρₐ * cp * T_air * depth_soil1 + 2 * ra_g * λ_soil1 * T_soil1_last) /
               (ρₐ * cp * depth_soil1 + 2 * ra_g * λ_soil1 + ttt * ra_g * depth_soil1)
@@ -58,8 +61,8 @@ function surface_temperature_jl(T_air::FT, rh_air::FT, depth_snow::FT, depth_wat
 
     T_snow = clamp(T_snow, T_air - 25.0, T_air + 25.0)
 
-    ttt = (λ_soil1 * T_soil1_last / depth_soil1 + T_snow * κ_snow + 0.02 * capacity_heat_soil1 / length_step * T_any0_last) /
-          (λ_soil1 / depth_soil1 + κ_snow / depth_snow + 0.02 * capacity_heat_soil1 / length_step)
+    ttt = (λ_soil1 * T_soil1_last / depth_soil1 + T_snow * κ_snow + 0.02 * cp_soil1 / length_step * T_any0_last) /
+          (λ_soil1 / depth_soil1 + κ_snow / depth_snow + 0.02 * cp_soil1 / length_step)
     T_any0 = T_soil0 * (1 - perc_snow_g) + ttt * perc_snow_g
 
     G_snow = κ_snow / (depth_snow + 0.5 * depth_soil1) * (T_snow - T_soil1_last)
@@ -98,7 +101,7 @@ function surface_temperature_jl(T_air::FT, rh_air::FT, depth_snow::FT, depth_wat
     G_snow2 = (T_snow2_last - T_any0_last) / ((0.5 * (depth_snow - 0.04) / κ_snow) + (0.02 / λ_soil1))
     T_snow2 = T_snow2_last + ((G_snow1 - G_snow2) / (cp_ice * ρ_snow * (depth_snow - 0.04))) * length_step
 
-    T_any0 = T_any0_last + ((G_snow2 - G_soil1) / (capacity_heat_soil0 * 0.02)) * length_step
+    T_any0 = T_any0_last + ((G_snow2 - G_soil1) / (cp_soil0 * 0.02)) * length_step
     T_soil0 = T_any0
 
     if T_snow > 0.0 && T_snow_last <= 0.0 && depth_snow > 0.0

diff --git a/test/modules/modules.jl b/test/modules/modules.jl
@@ -1,3 +1,4 @@
+include("test-surface_temperature.jl")
 include("test-radiation.jl")
 include("test-snowpack.jl")
 include("test-photosynthesis.jl")

diff --git a/test/modules/test-surface_temperature.jl b/test/modules/test-surface_temperature.jl
@@ -0,0 +1,46 @@
+@testset "surface_temperature_jl tests" begin
+  T_air = 2.0  # Air temperature in Celsius
+  rh_air = 0.5  # Relative humidity
+  depth_snow = 0.03  # Depth of snow in meters
+  depth_water = 0.1  # Depth of water in meters
+  cp_soil1 = 1000.0  # Heat capacity of soil layer 1
+  cp_soil0 = 800.0  # Heat capacity of soil layer 0
+  Gheat_g = 0.1  # Ground heat flux
+  depth_soil1 = 0.5  # Depth of soil layer 1 in meters
+  ρ_snow = 300.0  # Density of snow in kg/m^3
+  tempL_u = 0.0  # Temperature of the lower layer in Celsius
+  Rn_g = 200.0  # Net radiation on the ground
+  E_soil = 0.1  # Evaporation from soil
+  E_water_g = 0.05  # Evaporation from water on the ground
+  E_snow_g = 0.02  # Evaporation from snow on the ground
+  λ_soil1 = 0.5  # Thermal conductivity of soil layer 1
+  perc_snow_g = 0.3  # Percentage of snow cover on the ground
+  G_soil1 = 0.2  # Soil heat flux in layer 1
+  T_ground_last = 0.001 # Last ground temperature in Celsius
+  T_soil1_last = 0.001 # Last soil layer 1 temperature in Celsius
+  T_any0_last = 0.001 # Last temperature of any layer 0 in Celsius
+  T_soil0_last = 0.001 # Last soil layer 0 temperature in Celsius
+  T_snow_last = 0.001 # Last snow temperature in Celsius
+  T_snow1_last = 0.001 # Last snow layer 1 temperature in Celsius
+  T_snow2_last = 0.001 # Last snow layer 2 temperature in Celsius
+
+  # T_ground, T_any0, T_soil0, T_snow, T_snow1, T_snow2, G
+  r_jl = surface_temperature_jl(
+    T_air, rh_air,
+    depth_snow, depth_water,
+    cp_soil1, cp_soil0, Gheat_g,
+    depth_soil1, ρ_snow, tempL_u, Rn_g, E_soil, E_water_g, E_snow_g, λ_soil1, perc_snow_g, G_soil1,
+    T_ground_last, T_soil1_last, T_any0_last, T_soil0_last, T_snow_last, T_snow1_last, T_snow2_last
+  )
+
+  r_c = clang.surface_temperature_c(
+    T_air, rh_air,
+    depth_snow, depth_water,
+    cp_soil1, cp_soil0, Gheat_g,
+    depth_soil1, ρ_snow, tempL_u, Rn_g, E_soil, E_water_g, E_snow_g, λ_soil1, perc_snow_g, G_soil1,
+    T_ground_last, T_soil1_last, T_any0_last, T_soil0_last, T_snow_last, T_snow1_last, T_snow2_last
+  )
+  # println(r_jl)
+  # println(r_c)
+  @test all(isapprox.(r_jl, r_c, rtol=1e-8))
+end