Skip to content

Commit

Permalink
force bucket and richards model with ERA5
Browse files Browse the repository at this point in the history
  • Loading branch information
kmdeck committed Apr 30, 2024
1 parent 7b0c273 commit d611ad6
Show file tree
Hide file tree
Showing 6 changed files with 253 additions and 77 deletions.
4 changes: 0 additions & 4 deletions .buildkite/pipeline.yml
Original file line number Diff line number Diff line change
Expand Up @@ -49,10 +49,6 @@ steps:
command: "julia --color=yes --project=experiments experiments/standalone/Soil/evaporation.jl"
artifact_paths: "experiments/standalone/Soil/evaporation*png"

- label: "vaira_test"
command: "julia --color=yes --project=experiments experiments/integrated/fluxnet/run_fluxnet.jl US-Var"
artifact_paths: "experiments/integrated/fluxnet/US-Var/out/*png"

- label: "ozark_test"
command: "julia --color=yes --project=experiments experiments/integrated/fluxnet/run_fluxnet.jl US-MOz"
artifact_paths: "experiments/integrated/fluxnet/US-MOz/out/*png"
Expand Down
15 changes: 15 additions & 0 deletions experiments/Artifacts.toml
Original file line number Diff line number Diff line change
@@ -0,0 +1,15 @@
[lwp_pallardy_etal2018]
git-tree-sha1 = "29ed37c288601ab353f4f5dc13cf652ec4aec305"

["processed_topographic_index 2.5 deg"]
git-tree-sha1 = "31d4110d008fcfb1ecfd77eb9675c983b650eedf"

[era5_land_forcing_data2021]
git-tree-sha1 = "ec424296df6b60cfe273ac8f981701fbbed0bd8a"

[soil_params_Gupta2020_2022]
git-tree-sha1 = "ffc4c3fa76ea632ba46ae5d58671d4161124af29"

[[soil_params_Gupta2020_2022.download]]
sha256 = "636db48ecb387a9aa59afbff8af86f0ccc092690c688f0cbb52041c133361e22"
url = "https://caltech.box.com/shared/static/f2y23qx0lggjskftzgh7ht7fsbh36gmm.gz"
6 changes: 3 additions & 3 deletions experiments/Manifest.toml
Original file line number Diff line number Diff line change
Expand Up @@ -336,10 +336,10 @@ uuid = "595c0a79-7f3d-439a-bc5a-b232dc3bde79"
version = "0.7.19"

[[deps.ClimaUtilities]]
deps = ["Artifacts", "CFTime", "Dates"]
git-tree-sha1 = "d1b30a07829248b325ffffbb739d0f6fe79dfe7c"
deps = ["Artifacts", "Dates"]
git-tree-sha1 = "617ef474d78e2fd713626e22db8e9fb4881e90ee"
uuid = "b3f4f4ca-9299-4f7f-bd9b-81e1242a7513"
version = "0.1.3"
version = "0.1.5"
weakdeps = ["Adapt", "ClimaComms", "ClimaCore", "ClimaCoreTempestRemap", "Interpolations", "NCDatasets"]

[deps.ClimaUtilities.extensions]
Expand Down
171 changes: 118 additions & 53 deletions experiments/standalone/Bucket/global_bucket_staticmap.jl
Original file line number Diff line number Diff line change
Expand Up @@ -20,9 +20,10 @@ using CairoMakie
using Dates
using DelimitedFiles
using Statistics

using ClimaUtilities.ClimaArtifacts
import Interpolations
import ClimaUtilities.TimeVaryingInputs: TimeVaryingInput

import ClimaUtilities.ClimaArtifacts: @clima_artifact
import ClimaTimeSteppers as CTS
import NCDatasets
using ClimaCore
Expand Down Expand Up @@ -55,7 +56,14 @@ function compute_extrema(v)
return (minimum(mins), maximum(maxes))
end

anim_plots = false
function compute_clims(v)
means = [mean(u) for u in v]
sigmas = [std(u) for u in v]
return (minimum(means) - maximum(sigmas), maximum(means) + maximum(sigmas))
end

anim_plots = true
regridder_type = :InterpolationsRegridder
FT = Float64;
context = ClimaComms.context()
earth_param_set = LP.LandParameters(FT);
Expand All @@ -70,11 +78,11 @@ soil_depth = FT(3.5);
bucket_domain = ClimaLand.Domains.SphericalShell(;
radius = FT(6.3781e6),
depth = soil_depth,
nelements = (50, 10),
nelements = (30, 10),
npolynomial = 1,
dz_tuple = FT.((1.0, 0.05)),
);
ref_time = DateTime(2005);
ref_time = DateTime(2021);

# Set up parameters
σS_c = FT(0.2);
Expand All @@ -85,43 +93,90 @@ z_0b = FT(1e-3);
ρc_soil = FT(2e6);
τc = FT(3600);
t0 = 0.0;
tf = 7 * 86400;
Δt = 3600.0;
tf = 14 * 86400;
Δt = 3600.0 / 3;

# Construct albedo parameter object using static map
# Use separate regridding directory for CPU and GPU runs to avoid race condition
device_suffix =
typeof(ClimaComms.context().device) <: ClimaComms.CPUSingleThreaded ?
"cpu" : "gpu"
regrid_dir = joinpath(
pkgdir(ClimaLand),
"experiments/standalone/Bucket/$device_suffix/regrid-static/",
)
!ispath(regrid_dir) && mkpath(regrid_dir)
t_start = t0
surface_space = bucket_domain.space.surface
α_snow = FT(0.8)
albedo = PrescribedBaregroundAlbedo{FT}(α_snow, surface_space);

bucket_parameters = BucketModelParameters(FT; albedo, z_0m, z_0b, τc);

# Forcing data
era5_artifact_path = @clima_artifact("era5_land_forcing_data2021", context)
# Below, the preprocess_func argument is used to
# 1. Convert precipitation to be negative (as it is downwards)
# 2. Convert accumulations over an hour to a rate per second
# Precipitation:
precip = (t) -> 0;
snow_precip = (t) -> -5e-7 * (t < 1 * 86400);
# Diurnal temperature variations:
T_atmos = (t) -> 275.0 + 5.0 * sin(2.0 * π * t / 86400 - π / 2);
# Constant otherwise:
u_atmos = (t) -> 3.0;
q_atmos = (t) -> 0.001;
h_atmos = FT(2);
P_atmos = (t) -> 101325;
precip = TimeVaryingInput(
joinpath(era5_artifact_path, "era5_2021_0.9x1.25_clima.nc"),
"rf",
surface_space;
reference_date = ref_time,
t_start,
regridder_type,
file_reader_kwargs = (; preprocess_func = (data) -> -data / 3600,),
)

snow_precip = TimeVaryingInput(
joinpath(era5_artifact_path, "era5_2021_0.9x1.25.nc"),
"sf",
surface_space;
reference_date = ref_time,
t_start,
regridder_type,
file_reader_kwargs = (; preprocess_func = (data) -> -data / 3600,),
)

u_atmos = TimeVaryingInput(
joinpath(era5_artifact_path, "era5_2021_0.9x1.25_clima.nc"),
"ws",
surface_space;
reference_date = ref_time,
t_start,
regridder_type,
)
q_atmos = TimeVaryingInput(
joinpath(era5_artifact_path, "era5_2021_0.9x1.25_clima.nc"),
"q",
surface_space;
reference_date = ref_time,
t_start,
regridder_type,
)
P_atmos = TimeVaryingInput(
joinpath(era5_artifact_path, "era5_2021_0.9x1.25.nc"),
"sp",
surface_space;
reference_date = ref_time,
t_start,
regridder_type,
)

T_atmos = TimeVaryingInput(
joinpath(era5_artifact_path, "era5_2021_0.9x1.25.nc"),
"t2m",
surface_space;
reference_date = ref_time,
t_start,
regridder_type,
)
h_atmos = FT(10);


bucket_atmos = PrescribedAtmosphere(
TimeVaryingInput(precip),
TimeVaryingInput(snow_precip),
TimeVaryingInput(T_atmos),
TimeVaryingInput(u_atmos),
TimeVaryingInput(q_atmos),
TimeVaryingInput(P_atmos),
precip,
snow_precip,
T_atmos,
u_atmos,
q_atmos,
P_atmos,
ref_time,
h_atmos,
);
Expand All @@ -130,14 +185,26 @@ bucket_atmos = PrescribedAtmosphere(
# peak at local noon, and a prescribed downwelling LW radiative
# flux, assuming the air temperature is on average 275 degrees
# K with a diurnal amplitude of 5 degrees K:
SW_d = (t) -> max(1361 * sin(2π * t / 86400 - π / 2), 0.0);
LW_d = (t) -> 5.67e-8 * (275.0 + 5.0 * sin(2.0 * π * t / 86400 - π / 2))^4;
bucket_rad = PrescribedRadiativeFluxes(
FT,
TimeVaryingInput(SW_d),
TimeVaryingInput(LW_d),
ref_time,
);
SW_d = TimeVaryingInput(
joinpath(era5_artifact_path, "era5_2021_0.9x1.25.nc"),
"ssrd",
surface_space;
reference_date = ref_time,
t_start,
regridder_type,
file_reader_kwargs = (; preprocess_func = (data) -> data / 3600,),
)
LW_d = TimeVaryingInput(
joinpath(era5_artifact_path, "era5_2021_0.9x1.25.nc"),
"strd",
surface_space;
reference_date = ref_time,
t_start,
regridder_type,
file_reader_kwargs = (; preprocess_func = (data) -> data / 3600,),
)

bucket_rad = PrescribedRadiativeFluxes(FT, SW_d, LW_d, ref_time);

model = BucketModel(
parameters = bucket_parameters,
Expand Down Expand Up @@ -214,18 +281,22 @@ F_sfc = [
) for k in 1:length(sol.t)
];

sw_forcing = [
Remapping.interpolate(remapper, saved_values.saveval[k].drivers.SW_d)
for k in 1:length(sol.t)
];

# save prognostic state to CSV - for comparison between GPU and CPU output
open(joinpath(outdir, "tf_state_$device_suffix.txt"), "w") do io
writedlm(io, hcat(T_sfc[end][:], W[end][:], Ws[end][:], σS[end][:]), ',')
end;
# animation settings
nframes = length(W)
framerate = 2
nframes = length(T_sfc) # hourly data
fig_ts = Figure(size = (1000, 1000))
for (i, (field_ts, field_name)) in enumerate(
zip(
[W, σS, T_sfc, evaporation, F_sfc],
["W", "σS", "T_sfc", "evaporation", "F_sfc"],
[W, σS, T_sfc, evaporation, F_sfc, sw_forcing],
["W", "σS", "T_sfc", "evaporation", "F_sfc", "SW forcing"],
),
)
if anim_plots
Expand All @@ -236,13 +307,18 @@ for (i, (field_ts, field_name)) in enumerate(
ylabel = "Latitude",
title = field_name,
)
clims = compute_extrema(field_ts)
clims = compute_clims(field_ts)
CairoMakie.Colorbar(fig[:, end + 1], colorrange = clims)
outfile = joinpath(
outdir,
string("anim_$(device_suffix)_", field_name, ".mp4"),
)
record(fig, outfile, 1:nframes; framerate = framerate) do frame
record(
fig,
outfile,
(nframes - 7 * 24):2:nframes;
framerate = 3,
) do frame
CairoMakie.heatmap!(
longpts,
latpts,
Expand All @@ -251,18 +327,7 @@ for (i, (field_ts, field_name)) in enumerate(
)
end
end
# Plot the timeseries of the mean value as well.
xlabel = i == 5 ? "Time (days)" : ""
ax2 = Axis(
fig_ts[i, 1],
xlabel = xlabel,
ylabel = field_name,
title = "Global bucket with static map albedo",
)
CairoMakie.lines!(ax2, sol.t ./ 3600 ./ 24, [mean(x) for x in field_ts])

end
outfile = joinpath(outdir, string("ts_$device_suffix.png"))
CairoMakie.save(outfile, fig_ts)

# delete regrid directory
rm(regrid_dir, recursive = true)
10 changes: 10 additions & 0 deletions experiments/standalone/Soil/Artifacts.toml
Original file line number Diff line number Diff line change
@@ -0,0 +1,10 @@
["processed_topographic_index 2.5 deg"]
git-tree-sha1 = "31d4110d008fcfb1ecfd77eb9675c983b650eedf"

[soil_params_Gupta2020_2022]
git-tree-sha1 = "ffc4c3fa76ea632ba46ae5d58671d4161124af29"

[[soil_params_Gupta2020_2022.download]]
sha256 = "636db48ecb387a9aa59afbff8af86f0ccc092690c688f0cbb52041c133361e22"
url = "https://caltech.box.com/shared/static/f2y23qx0lggjskftzgh7ht7fsbh36gmm.gz"

Loading

0 comments on commit d611ad6

Please sign in to comment.