From 180d62d6f7ed2266df5e844f99ac943727644e93 Mon Sep 17 00:00:00 2001 From: Brian Groenke <brian.groenke@awi.de> Date: Tue, 26 Nov 2024 00:40:02 +0100 Subject: [PATCH] Fix failing tests and remove Presets submodule --- README.md | 6 +++--- docs/src/manual/architecture.md | 2 +- docs/src/manual/overview.md | 2 +- docs/src/quickstart.md | 8 +++---- examples/cglite_parameter_ensembles.jl | 2 +- examples/heat_freeW_bucketW_samoylov.jl | 4 ++-- examples/heat_freeW_lake_lite_implicit.jl | 2 +- examples/heat_freeW_lite_implicit.jl | 2 +- examples/heat_freeW_samoylov.jl | 9 +++----- .../heat_freeW_seb_snow_bucketW_samoylov.jl | 6 +++--- examples/heat_freeW_snow_samoylov.jl | 4 ++-- examples/heat_sfcc_constantbc.jl | 4 ++-- examples/heat_sfcc_richardseq_samoylov.jl | 2 +- examples/heat_sfcc_salt_constantbc.jl | 2 +- examples/heat_sfcc_samoylov.jl | 9 ++++---- examples/heat_simple_autodiff_grad.jl | 6 +++--- src/CryoGrid.jl | 8 +++++-- src/IO/forcings/forcings_json.jl | 2 +- src/Physics/Soils/para/simple.jl | 4 ++-- src/Physics/Surface/SEB/seb.jl | 2 +- src/Physics/Surface/SEB/seb_state.jl | 6 +++--- src/{Presets/Presets.jl => models.jl} | 21 ------------------- test/CryoGridSBIExt/runtests.jl | 5 +++-- test/IO/forcing_tests.jl | 21 ++++++++++--------- test/Physics/Heat/heat_conduction_tests.jl | 17 ++++++++------- test/Physics/Hydrology/water_balance_tests.jl | 12 +++++------ test/Physics/Salt/runtests.jl | 2 +- test/Physics/Surface/seb_tests.jl | 12 +++++------ test/Physics/Surface/swb_tests.jl | 11 +++------- test/Solvers/cglite_implicit_tests.jl | 5 ++--- test/benchmarks.jl | 4 ++-- test/test_problems.jl | 6 +++--- 32 files changed, 94 insertions(+), 114 deletions(-) rename src/{Presets/Presets.jl => models.jl} (86%) diff --git a/README.md b/README.md index 7d41667a..2127659d 100755 --- a/README.md +++ b/README.md @@ -51,12 +51,12 @@ using Plots: plot # The forcing file will be automatically downloaded to the input/ folder if not already present. forcings = loadforcings(CryoGrid.Forcings.Samoylov_ERA_obs_fitted_1979_2014_spinup_extended_2044); # get preset soil and initial temperature profile for Samoylov -soilprofile, tempprofile = CryoGrid.Presets.SamoylovDefault +soilprofile, tempprofile = CryoGrid.SamoylovDefault initT = initializer(:T, tempprofile) # choose grid with 5cm spacing -grid = CryoGrid.Presets.DefaultGrid_5cm +grid = CryoGrid.DefaultGrid_5cm # build a default Tile with heat conduction on the given soil profile -tile = CryoGrid.Presets.SoilHeatTile( +tile = CryoGrid.SoilHeatTile( TemperatureBC(forcings.Tair), GeothermalHeatFlux(0.053u"W/m^2"), soilprofile, diff --git a/docs/src/manual/architecture.md b/docs/src/manual/architecture.md index eb0720e2..5cadab06 100644 --- a/docs/src/manual/architecture.md +++ b/docs/src/manual/architecture.md @@ -65,7 +65,7 @@ using CryoGrid using CryoGrid.Diagnostics soil = Ground() -grid = CryoGrid.Presets.DefaultGrid_5cm +grid = CryoGrid.DefaultGrid_5cm state = Diagnostics.build_dummy_state(grid, soil) @which CryoGrid.computediagnostic!(soil, state) diff --git a/docs/src/manual/overview.md b/docs/src/manual/overview.md index 03357710..414f950a 100644 --- a/docs/src/manual/overview.md +++ b/docs/src/manual/overview.md @@ -18,7 +18,7 @@ strat = Stratigraphy( 0.0u"m" => Ground(soilprofile, HeatBalance(:H; freezecurve=DallAmico())), 1000.0u"m" => Bottom(GeothermalHeatFlux(0.053u"J/s/m^2")) ); -grid = CryoGrid.Presets.DefaultGrid_5cm +grid = CryoGrid.DefaultGrid_5cm # define initial conditions for temperature using a given profile; # The default initializer linearly interpolates between profile points. initT = initializer(:T, tempprofile) diff --git a/docs/src/quickstart.md b/docs/src/quickstart.md index 2d887e00..75c28cb8 100644 --- a/docs/src/quickstart.md +++ b/docs/src/quickstart.md @@ -1,6 +1,6 @@ ## Quick start -After [installing](installation.md) `CryoGrid.jl`, you can get started right away with a simple soil heat model. The [`Presets`](@ref) module (aliased `CryoGrid.Presets`) provides pre-specified model configurations that can be obtained with a single function call. It is also possible to provide a custom soil and temperature profile to `SoilHeatTile`; here the values given in `Presets.SamoylovDefault` is used. +After [installing](installation.md) `CryoGrid.jl`, you can get started right away with a simple soil heat model. The [`Presets`](@ref) module (aliased `CryoGrid.Presets`) provides pre-specified model configurations that can be obtained with a single function call. It is also possible to provide a custom soil and temperature profile to `SoilHeatTile`; here the values given in `SamoylovDefault` is used. ```julia using CryoGrid @@ -10,12 +10,12 @@ using Plots # The forcing file will be automatically downloaded to the input/ folder if not already present. forcings = loadforcings(CryoGrid.Forcings.Samoylov_ERA_obs_fitted_1979_2014_spinup_extended_2044, :Tair => u"°C"); # get preset soil and initial temperature profile for Samoylov -soilprofile, tempprofile = CryoGrid.Presets.SamoylovDefault +soilprofile, tempprofile = CryoGrid.SamoylovDefault initT = initializer(:T, tempprofile) # choose grid with 5cm spacing -grid = CryoGrid.Presets.DefaultGrid_5cm +grid = CryoGrid.DefaultGrid_5cm # build Tile from the given soil and temperature profiles -tile = CryoGrid.Presets.SoilHeatTile(TemperatureBC(forcings.Tair), GeothermalHeatFlux(0.053u"W/m^2"), soilprofile, initT, grid=grid) +tile = CryoGrid.SoilHeatTile(TemperatureBC(forcings.Tair), GeothermalHeatFlux(0.053u"W/m^2"), soilprofile, initT, grid=grid) # define time span (1 year) tspan = (DateTime(2010,11,30),DateTime(2011,11,30)) u0, du0 = initialcondition!(tile, tspan) diff --git a/examples/cglite_parameter_ensembles.jl b/examples/cglite_parameter_ensembles.jl index 91bb16c8..4c6c40e5 100644 --- a/examples/cglite_parameter_ensembles.jl +++ b/examples/cglite_parameter_ensembles.jl @@ -42,7 +42,7 @@ strat = Stratigraphy( soil_layers, z_bot => Bottom(GeothermalHeatFlux(0.053u"W/m^2")) ); -modelgrid = CryoGrid.Presets.DefaultGrid_2cm +modelgrid = CryoGrid.DefaultGrid_2cm tile = Tile(strat, modelgrid, ssinit); # Since the solver can take daily timesteps, we can easily specify longer simulation time spans at minimal cost. # Here we specify a time span of 10 years. diff --git a/examples/heat_freeW_bucketW_samoylov.jl b/examples/heat_freeW_bucketW_samoylov.jl index 96d67625..355b734e 100644 --- a/examples/heat_freeW_bucketW_samoylov.jl +++ b/examples/heat_freeW_bucketW_samoylov.jl @@ -6,7 +6,7 @@ # Frist, load forcings and define boundary conditions. using CryoGrid forcings = loadforcings(CryoGrid.Forcings.Samoylov_ERA_obs_fitted_1979_2014_spinup_extended_2044); -_, tempprofile = CryoGrid.Presets.SamoylovDefault; +_, tempprofile = CryoGrid.SamoylovDefault; initT = initializer(:T, tempprofile) initsat = initializer(:sat, (l,state) -> state.sat .= l.para.sat); upperbc = WaterHeatBC( @@ -23,7 +23,7 @@ strat = @Stratigraphy( 2.0u"m" => Ground(SimpleSoil(por=0.10,sat=1.0,org=0.0), heat=HeatBalance(), water=WaterBalance(BucketScheme())), 1000.0u"m" => Bottom(GeothermalHeatFlux(0.053u"W/m^2")), ); -modelgrid = CryoGrid.Presets.DefaultGrid_2cm +modelgrid = CryoGrid.DefaultGrid_2cm tile = Tile(strat, modelgrid, forcings, initT, initsat); # Now we set up the problem and solve using the integrator interface. diff --git a/examples/heat_freeW_lake_lite_implicit.jl b/examples/heat_freeW_lake_lite_implicit.jl index c0690472..dc4d2da8 100644 --- a/examples/heat_freeW_lake_lite_implicit.jl +++ b/examples/heat_freeW_lake_lite_implicit.jl @@ -20,7 +20,7 @@ tempprofile_linear = TemperatureProfile( 10.0u"m" => -10.0u"°C", 1000.0u"m" => 10.2u"°C" ) -modelgrid = Grid(vcat(-1.0u"m":0.02u"m":-0.02u"m", CryoGrid.Presets.DefaultGrid_2cm)) +modelgrid = Grid(vcat(-1.0u"m":0.02u"m":-0.02u"m", CryoGrid.DefaultGrid_2cm)) z_top = -1.0u"m" z_sub = keys(soilprofile) z_bot = modelgrid[end] diff --git a/examples/heat_freeW_lite_implicit.jl b/examples/heat_freeW_lite_implicit.jl index eb5b9931..aa88e557 100644 --- a/examples/heat_freeW_lite_implicit.jl +++ b/examples/heat_freeW_lite_implicit.jl @@ -36,7 +36,7 @@ strat = @Stratigraphy( soil_layers..., z_bot => Bottom(GeothermalHeatFlux(0.053u"W/m^2")) ); -modelgrid = Grid(vcat(z_top:0.02u"m":-0.02u"m", CryoGrid.Presets.DefaultGrid_2cm)) +modelgrid = Grid(vcat(z_top:0.02u"m":-0.02u"m", CryoGrid.DefaultGrid_2cm)) tile = Tile(strat, modelgrid, forcings, ssinit); # Since the solver can take daily timesteps, we can easily specify longer simulation time spans at minimal cost. diff --git a/examples/heat_freeW_samoylov.jl b/examples/heat_freeW_samoylov.jl index e2eb192b..a54e0806 100644 --- a/examples/heat_freeW_samoylov.jl +++ b/examples/heat_freeW_samoylov.jl @@ -23,14 +23,14 @@ soilprofile = SoilProfile( ); # We construct a state variable initializer for temperature `T` from the temperature profile preset for Samoylov. -initT = initializer(:T, CryoGrid.Presets.SamoylovDefault.tempprofile); +initT = initializer(:T, CryoGrid.SamoylovDefault.tempprofile); # We choose the default grid discretization with 5 cm spacing at the surface. grid = CryoGrid.DefaultGrid_5cm; # Now we construct the Tile using the built-in model configuration `SoilHeatTile` which defines a # standalone soil straigraphy with only heat conduction and no water flow. -tile = CryoGrid.Presets.SoilHeatTile( +tile = CryoGrid.SoilHeatTile( :H, TemperatureBC(Input(:Tair), NFactor(nf=Param(0.6), nt=Param(0.9))), GeothermalHeatFlux(0.053u"W/m^2"), @@ -59,7 +59,4 @@ out = CryoGridOutput(sol) import Plots zs = [1,10,20,30,50,100,200,500,1000]u"cm" cg = Plots.cgrad(:copper,rev=true); -Plots.plot(out.T[Z(Near(zs))], color=cg[LinRange(0.0,1.0,length(zs))]', ylabel="Temperature", leg=false, size=(800,500), dpi=150) - -using BenchmarkTools -@profview @btime $tile($du0, $u0, $prob.p, $prob.tspan[1]) +Plots.plot(out.T[Z(Near(zs))], color=cg[LinRange(0.0,1.0,length(zs))]', ylabel="Temperature", leg=false, size=(800,500), dpi=150) \ No newline at end of file diff --git a/examples/heat_freeW_seb_snow_bucketW_samoylov.jl b/examples/heat_freeW_seb_snow_bucketW_samoylov.jl index 4248173e..2f47c9d5 100644 --- a/examples/heat_freeW_seb_snow_bucketW_samoylov.jl +++ b/examples/heat_freeW_seb_snow_bucketW_samoylov.jl @@ -8,7 +8,7 @@ using CryoGrid using OrdinaryDiffEq # First, load the forcings and construct the Tile. -modelgrid = CryoGrid.Presets.DefaultGrid_2cm; +modelgrid = CryoGrid.DefaultGrid_2cm; soilprofile = SoilProfile( 0.0u"m" => SimpleSoil(por=0.80,sat=1.0,org=0.75), 0.1u"m" => SimpleSoil(por=0.80,sat=1.0,org=0.25), @@ -17,9 +17,9 @@ soilprofile = SoilProfile( 10.0u"m" => SimpleSoil(por=0.30,sat=1.0,org=0.0), ); ## mid-winter temperature profile -tempprofile = CryoGrid.Presets.SamoylovDefault.tempprofile +tempprofile = CryoGrid.SamoylovDefault.tempprofile forcings = loadforcings(CryoGrid.Forcings.Samoylov_ERA_obs_fitted_1979_2014_spinup_extended_2044); -tempprofile = CryoGrid.Presets.SamoylovDefault.tempprofile +tempprofile = CryoGrid.SamoylovDefault.tempprofile initT = initializer(:T, tempprofile) seb = SurfaceEnergyBalance() swb = SurfaceWaterBalance() diff --git a/examples/heat_freeW_snow_samoylov.jl b/examples/heat_freeW_snow_samoylov.jl index 2a423222..2bdc5d3e 100644 --- a/examples/heat_freeW_snow_samoylov.jl +++ b/examples/heat_freeW_snow_samoylov.jl @@ -15,7 +15,7 @@ soilprofile = SoilProfile( 3.0u"m" => SimpleSoil(por=0.50,sat=1.0,org=0.0), 10.0u"m" => SimpleSoil(por=0.30,sat=1.0,org=0.0), ); -initT = initializer(:T, CryoGrid.Presets.SamoylovDefault.tempprofile) +initT = initializer(:T, CryoGrid.SamoylovDefault.tempprofile) initsat = initializer(:sat, 1.0) z_top = -2.0u"m" z_sub = keys(soilprofile) @@ -39,7 +39,7 @@ strat = @Stratigraphy( ground_layers..., z_bot => Bottom(GeothermalHeatFlux(0.053u"J/s/m^2")) ); -modelgrid = CryoGrid.Presets.DefaultGrid_5cm +modelgrid = CryoGrid.DefaultGrid_5cm tile = Tile(strat, modelgrid, forcings, initT, initsat) # define time span, 2 years + 3 months tspan = (DateTime(2010,9,30), DateTime(2012,9,30)) diff --git a/examples/heat_sfcc_constantbc.jl b/examples/heat_sfcc_constantbc.jl index 11cee116..6cc3f12c 100644 --- a/examples/heat_sfcc_constantbc.jl +++ b/examples/heat_sfcc_constantbc.jl @@ -8,7 +8,7 @@ using CryoGrid # Select default grid and initial temperature profile. -grid = CryoGrid.Presets.DefaultGrid_10cm +grid = CryoGrid.DefaultGrid_10cm tempprofile = TemperatureProfile( 0.0u"m" => -20.0u"°C", 1000.0u"m" => 10.0u"°C", @@ -31,7 +31,7 @@ Plots.plot(-2.0u"°C":0.01u"K":0.0u"°C", sfcc) # specified. heatop = Heat.Diffusion1D(:H) initT = initializer(:T, tempprofile) -tile = CryoGrid.Presets.SoilHeatTile( +tile = CryoGrid.SoilHeatTile( heatop, ## 10 W/m^2 in and out, i.e. net zero flux ConstantBC(HeatBalance, CryoGrid.Neumann, 10.0u"W/m^2"), diff --git a/examples/heat_sfcc_richardseq_samoylov.jl b/examples/heat_sfcc_richardseq_samoylov.jl index 36aada34..877d07f8 100644 --- a/examples/heat_sfcc_richardseq_samoylov.jl +++ b/examples/heat_sfcc_richardseq_samoylov.jl @@ -38,7 +38,7 @@ strat = @Stratigraphy( 2.0u"m" => Ground(SimpleSoil(por=0.10,sat=1.0,org=0.0,freezecurve=sfcc); heat, water), 1000.0u"m" => Bottom(GeothermalHeatFlux(0.053u"W/m^2")) ); -grid = CryoGrid.Presets.DefaultGrid_2cm +grid = CryoGrid.DefaultGrid_2cm tile = Tile(strat, grid, forcings, initT); u0, du0 = @time initialcondition!(tile, tspan) prob = CryoGridProblem(tile, u0, tspan, saveat=3*3600, savevars=(:T,:θw,:θwi,:kw)); diff --git a/examples/heat_sfcc_salt_constantbc.jl b/examples/heat_sfcc_salt_constantbc.jl index bfc22d37..e686b013 100644 --- a/examples/heat_sfcc_salt_constantbc.jl +++ b/examples/heat_sfcc_salt_constantbc.jl @@ -14,7 +14,7 @@ strat = @Stratigraphy( 0.0u"m" => SalineGround(SimpleSoil(freezecurve=sfcc), heat=HeatBalance(:T), salt=SaltMassBalance()), 1000.0u"m" => Bottom(GeothermalHeatFlux(0.053u"W/m^2")) ); -modelgrid = CryoGrid.Presets.DefaultGrid_10cm +modelgrid = CryoGrid.DefaultGrid_10cm tile = Tile(strat, modelgrid, initT, initsalt, initpor) tspan = (DateTime(1990,1,1),DateTime(2000,12,31)) u0, du0 = initialcondition!(tile, tspan) diff --git a/examples/heat_sfcc_samoylov.jl b/examples/heat_sfcc_samoylov.jl index 3b6b8119..add6edc0 100644 --- a/examples/heat_sfcc_samoylov.jl +++ b/examples/heat_sfcc_samoylov.jl @@ -5,14 +5,15 @@ using CryoGrid -# First we set up the soil heat model: +# First we set up the soil heat model. Note that the default soil profile +# for Samoylov already includes the appropriate freeze curves. forcings = loadforcings(CryoGrid.Forcings.Samoylov_ERA5_fitted_daily_1979_2020); -grid = CryoGrid.Presets.DefaultGrid_5cm -soilprofile, tempprofile = CryoGrid.Presets.SamoylovDefault +grid = CryoGrid.DefaultGrid_5cm +soilprofile, tempprofile = CryoGrid.SamoylovDefault initT = initializer(:T, tempprofile) upperbc = TemperatureBC(Input(:Tair), NFactor(nf=0.6, nt=0.9)) lowerbc = GeothermalHeatFlux(0.053u"W/m^2") -tile = CryoGrid.Presets.SoilHeatTile(upperbc, lowerbc, soilprofile, forcings, initT; grid=grid) +tile = CryoGrid.SoilHeatTile(upperbc, lowerbc, soilprofile, forcings, initT; grid=grid) tspan = (DateTime(2010,10,30),DateTime(2011,10,30)) u0, du0 = @time initialcondition!(tile, tspan) prob = CryoGridProblem(tile, u0, tspan, saveat=3*3600.0, savevars=(:T,)); diff --git a/examples/heat_simple_autodiff_grad.jl b/examples/heat_simple_autodiff_grad.jl index 5b88521e..99f74eed 100644 --- a/examples/heat_simple_autodiff_grad.jl +++ b/examples/heat_simple_autodiff_grad.jl @@ -7,10 +7,10 @@ # Set up forcings and boundary conditions similarly to other examples: using CryoGrid forcings = loadforcings(CryoGrid.Forcings.Samoylov_ERA_obs_fitted_1979_2014_spinup_extended_2044); -soilprofile, tempprofile = CryoGrid.Presets.SamoylovDefault -grid = CryoGrid.Presets.DefaultGrid_5cm +soilprofile, tempprofile = CryoGrid.SamoylovDefault +grid = CryoGrid.DefaultGrid_5cm initT = initializer(:T, tempprofile) -tile = CryoGrid.Presets.SoilHeatTile( +tile = CryoGrid.SoilHeatTile( :T, TemperatureBC(Input(:Tair), NFactor(nf=Param(0.5), nt=Param(0.9))), GeothermalHeatFlux(0.053u"W/m^2"), diff --git a/src/CryoGrid.jl b/src/CryoGrid.jl index bad4afaa..1f13ec07 100755 --- a/src/CryoGrid.jl +++ b/src/CryoGrid.jl @@ -36,6 +36,8 @@ import Interpolations @reexport using Unitful @reexport using UnPack +using SciMLBase: intervals # resolve name collision + export Interpolations # Common types and methods @@ -79,6 +81,7 @@ include("IO/InputOutput.jl") include("Tiles/Tiles.jl") @reexport using .Tiles +using .Tiles: layers # resolve name collision parameters = Tiles.parameters export ConstantBC, PeriodicBC, ConstantValue, PeriodicValue, ConstantFlux, PeriodicFlux @@ -99,8 +102,9 @@ include("problem.jl") export CGEuler, CryoGridIntegrator, CryoGridSolution include("Solvers/Solvers.jl") -# Presets -include("Presets/Presets.jl") +# Built-in model definitions +export SoilHeatTile, SamoylovDefault +include("models.jl") using .InputOutput: Resource diff --git a/src/IO/forcings/forcings_json.jl b/src/IO/forcings/forcings_json.jl index 3c56b0d8..0a19ff87 100644 --- a/src/IO/forcings/forcings_json.jl +++ b/src/IO/forcings/forcings_json.jl @@ -41,7 +41,7 @@ function loadforcings(format::ForcingFormatJSON{1}, filename::String, kwargs...) end return Interpolated1D(DimStack((; forcings...)); kwargs...) end -function loadforcings(format::ForcingFormatJSON{2}, filename::String) +function loadforcings(format::ForcingFormatJSON{2}, filename::String; kwargs...) dict = open(filename, "r") do file; JSON3.read(file) end # convert JSON3 dict for data field to Julia dict data = Dict(dict[:data]...) diff --git a/src/Physics/Soils/para/simple.jl b/src/Physics/Soils/para/simple.jl index 213b902b..13c7b825 100644 --- a/src/Physics/Soils/para/simple.jl +++ b/src/Physics/Soils/para/simple.jl @@ -35,12 +35,12 @@ SoilThermalProperties( kh_i = ThermalProperties().kh_i, kh_a = ThermalProperties().kh_a, kh_o=Param(0.25, units=u"W/m/K", domain=0..Inf), # organic [Hillel (1982)] - kh_m=Param(3.8, untis=u"W/m/K", domain=0..Inf), # mineral [Hillel (1982)] + kh_m=Param(3.8, units=u"W/m/K", domain=0..Inf), # mineral [Hillel (1982)] ch_w = ThermalProperties().ch_w, ch_i = ThermalProperties().ch_i, ch_a = ThermalProperties().ch_a, ch_o=Param(2.5e6, units=u"J/K/m^3", domain=0..Inf), # heat capacity organic - ch_m=Param(2.0e6, untis=u"J/K/m^3", domain=0..Inf), # heat capacity mineral + ch_m=Param(2.0e6, units=u"J/K/m^3", domain=0..Inf), # heat capacity mineral kwargs..., ) = ThermalProperties(; kh_w, kh_i, kh_a, kh_m, kh_o, ch_w, ch_i, ch_a, ch_m, ch_o, kwargs...) diff --git a/src/Physics/Surface/SEB/seb.jl b/src/Physics/Surface/SEB/seb.jl index 4a98aa7b..b5245cdb 100644 --- a/src/Physics/Surface/SEB/seb.jl +++ b/src/Physics/Surface/SEB/seb.jl @@ -96,7 +96,7 @@ function SurfaceEnergyBalance(; stabfun::StabilityFunctions = HøgstrømSHEBA(), ) inputs = (; Tair, Lin, Sin, pr, qh, wind) - SurfaceEnergyBalance(; inputs, para, solscheme, stabfun) + SurfaceEnergyBalance(inputs, para, solscheme, stabfun) end """ diff --git a/src/Physics/Surface/SEB/seb_state.jl b/src/Physics/Surface/SEB/seb_state.jl index d2703a48..97217968 100644 --- a/src/Physics/Surface/SEB/seb_state.jl +++ b/src/Physics/Surface/SEB/seb_state.jl @@ -1,10 +1,10 @@ """ - SEBInputs{TT1,TT2,TR,TP,TQ,TW,TZ,Tsurf} + SEBInputs{TT1,TT2,TR,TP,TQ,TW,Tsurf} Non-prognostic input variables to the SEB that are assumed to be given. """ -Base.@kwdef struct SEBInputs{TT1,TT2,TR,TP,TQ,TW,TZ,Tsurf} - Ts::TT1 = nothing +Base.@kwdef struct SEBInputs{TT1,TT2,TR,TP,TQ,TW,Tsurf} + Ts::TT1 = 0.0u"°C" Tair::TT2 = 0.0u"°C" Lin::TR = 0.0u"W/m^2" Sin::TR = 0.0u"W/m^2" diff --git a/src/Presets/Presets.jl b/src/models.jl similarity index 86% rename from src/Presets/Presets.jl rename to src/models.jl index 2978ec7e..ce5ebbce 100644 --- a/src/Presets/Presets.jl +++ b/src/models.jl @@ -1,22 +1,3 @@ -""" -Pre-built CryoGrid configurations for rapid prototyping. -""" -module Presets - -using CryoGrid -using CryoGrid.InputOutput: Resource -using CryoGrid.Numerics -using CryoGrid.Utils - -# physics modules -using CryoGrid.Heat -using CryoGrid.Hydrology -using CryoGrid.Soils - -using Statistics - -export SoilHeatTile, SamoylovDefault - """ SoilHeatTile([heatop=:H], upperbc::BoundaryProcess, soilprofile::Profile, init::CryoGrid.Initializer...; grid::Grid=DefaultGrid_10cm, tile_kwargs...) where {F<:FreezeCurve} @@ -61,5 +42,3 @@ const SamoylovDefault = ( 1000.0u"m" => 10.2u"°C", ) ) - -end \ No newline at end of file diff --git a/test/CryoGridSBIExt/runtests.jl b/test/CryoGridSBIExt/runtests.jl index 5ddb3e0b..aa2f09df 100644 --- a/test/CryoGridSBIExt/runtests.jl +++ b/test/CryoGridSBIExt/runtests.jl @@ -15,14 +15,15 @@ function cryogrid_test_problem(tspan) ); heatop = Heat.Diffusion1D(:H) initT = initializer(:T, tempprofile) - tile = CryoGrid.Presets.SoilHeatTile( + tile = CryoGrid.SoilHeatTile( heatop, ## 10 W/m^2 in and out, i.e. net zero flux PeriodicBC(HeatBalance, CryoGrid.Dirichlet, 365.0u"d", 15.0u"°C", 0.0, -5.0u"°C"), ConstantBC(HeatBalance, CryoGrid.Neumann, -0.1u"W/m^2"), soilprofile, + inputs(), initT; - grid=CryoGrid.Presets.DefaultGrid_10cm, + grid=CryoGrid.DefaultGrid_10cm, ); u0, _ = initialcondition!(tile, tspan); prob = CryoGridProblem(tile, u0, tspan) diff --git a/test/IO/forcing_tests.jl b/test/IO/forcing_tests.jl index 47be0699..60f2dfda 100644 --- a/test/IO/forcing_tests.jl +++ b/test/IO/forcing_tests.jl @@ -2,23 +2,24 @@ using CryoGrid using Dates using Test, BenchmarkTools -@testset "InterpolatedForcing" begin +@testset "Interpolated1D" begin ts = DateTime(2020,1,1):Day(1):DateTime(2020,3,1) ys = randn(length(ts)) - forcing = CryoGrid.InterpolatedForcing(ts, ys, :test) - @test all(forcing.interpolant.knots[1] .≈ convert_t.(ts)) - @test forcing(Dates.datetime2epochms(ts[1])/1000.0) ≈ ys[1] - @test forcing(Dates.datetime2epochms(ts[end])/1000.0) ≈ ys[end] - @test_throws BoundsError forcing(-1.0) - @test_throws BoundsError forcing(Dates.datetime2epochms(ts[end])/1000.0+1.0) + arr = DimStack((test=DimArray(ys, (Ti(ts),)),)) + interp = InputOutput.Interpolated1D(arr) + f = interp.test + @test f(Dates.datetime2epochms(ts[1])/1000.0) ≈ ys[1] + @test f(Dates.datetime2epochms(ts[end])/1000.0) ≈ ys[end] + @test_throws BoundsError f(-1.0) + @test_throws BoundsError f(Dates.datetime2epochms(ts[end])/1000.0+1.0) t1,y1 = ts[1], ys[1] t2,y2 = ts[2], ys[2] - @test forcing((Dates.datetime2epochms(t1) + Dates.datetime2epochms(t2))/2000.0) ≈ (y1+y2)/2 + @test f((Dates.datetime2epochms(t1) + Dates.datetime2epochms(t2))/2000.0) ≈ (y1+y2)/2 t = Dates.datetime2epochms(t1)/1000.0 - benchres = @benchmark $forcing($t) + benchres = @benchmark $f($t) @test benchres.allocs == 0 out = zeros(100) queries = t .+ (1:100); - benchres = @benchmark $out .= $forcing.($queries) + benchres = @benchmark $out .= $f.($queries) @test benchres.allocs == 0 end \ No newline at end of file diff --git a/test/Physics/Heat/heat_conduction_tests.jl b/test/Physics/Heat/heat_conduction_tests.jl index 692cf824..78ed4831 100644 --- a/test/Physics/Heat/heat_conduction_tests.jl +++ b/test/Physics/Heat/heat_conduction_tests.jl @@ -79,19 +79,22 @@ using Test end @testset "Boundary conditions" begin @testset "n-factors" begin - ts = DateTime(2010,1,1):Hour(1):DateTime(2010,1,1,4) - forcing = InterpolatedForcing(ts, [1.0,0.5,-0.5,-1.0,0.1]u"°C", :Tair) - tgrad = TemperatureBC(forcing, NFactor(nf=0.5, nt=1.0)) + T_pos = TemperatureBC(10.0, NFactor(nf=0.5, nt=0.75)) + T_neg = TemperatureBC(-10.0, NFactor(nf=0.5, nt=0.75)) heat = HeatBalance() sub = TestGroundLayer(heat) zerobc = ConstantBC(HeatBalance, CryoGrid.Dirichlet, 0.0u"°C") function f1(t) state = (T_ub=[Inf], nfactor=[Inf], t=t) - computediagnostic!(Top(zerobc), tgrad, state) - return boundaryvalue(tgrad, state) + computediagnostic!(Top(T_pos), T_pos, state) + T_ub = boundaryvalue(T_pos, state) + computediagnostic!(Top(T_neg), T_neg, state) + T_lb = boundaryvalue(T_neg, state) + return T_ub, T_lb end - Tres = f1.(Dates.datetime2epochms.(ts)./1000.0) - @test all(Tres .≈ [1.0,0.5,-0.25,-0.5,0.1]) + T1, T2 = f1(0.0) + @test T1 == 7.5 + @test T2 == -5.0 end end @testset "Fourier solution" begin diff --git a/test/Physics/Hydrology/water_balance_tests.jl b/test/Physics/Hydrology/water_balance_tests.jl index ac819ccd..b4f08730 100644 --- a/test/Physics/Hydrology/water_balance_tests.jl +++ b/test/Physics/Hydrology/water_balance_tests.jl @@ -6,9 +6,9 @@ using CryoGrid.Utils using Test @testset "Bucket scheme" begin - testgrid = CryoGrid.Presets.DefaultGrid_2cm + testgrid = CryoGrid.DefaultGrid_2cm @testset "variables" begin - sub = TestGroundLayer(WaterBalance(BucketScheme())) + sub = TestGroundLayer(WaterBalance(BucketScheme())) |> stripparams vars = CryoGrid.variables(sub) prognostic_vars = filter(CryoGrid.isprognostic, vars) diagnostic_vars = filter(CryoGrid.isdiagnostic, vars) @@ -23,7 +23,7 @@ using Test @test isa(proc, WaterBalance) end @testset "computediagnostic!" begin - sub = TestGroundLayer(WaterBalance(BucketScheme())) + sub = TestGroundLayer(WaterBalance(BucketScheme())) |> stripparams state = Diagnostics.build_dummy_state(testgrid, sub, with_units=true) # initialize fully saturated state.sat .= 1.0 @@ -33,7 +33,7 @@ using Test @test allfinite(state.θw) end @testset "computeprognostic!" begin - sub = TestGroundLayer(WaterBalance(BucketScheme())) + sub = TestGroundLayer(WaterBalance(BucketScheme())) |> stripparams state = Diagnostics.build_dummy_state(testgrid, sub, with_units=true) # initialize fully saturated state.sat .= 1.0 @@ -53,8 +53,8 @@ using Test @test all(state.jw[2:end-1] .> zero(eltype(state.jw))) end @testset "interact!" begin - sub1 = TestGroundLayer(WaterBalance(BucketScheme())) - sub2 = TestGroundLayer(WaterBalance(BucketScheme())) + sub1 = TestGroundLayer(WaterBalance(BucketScheme())) |> stripparams + sub2 = TestGroundLayer(WaterBalance(BucketScheme())) |> stripparams state1 = Diagnostics.build_dummy_state(testgrid[0.0u"m"..10.0u"m"], sub1, with_units=true) state2 = Diagnostics.build_dummy_state(testgrid[10.0u"m"..1000.0u"m"], sub2, with_units=true) # case 1: both saturated diff --git a/test/Physics/Salt/runtests.jl b/test/Physics/Salt/runtests.jl index d9dc4845..70f5154a 100644 --- a/test/Physics/Salt/runtests.jl +++ b/test/Physics/Salt/runtests.jl @@ -5,7 +5,7 @@ using CryoGrid.Utils using Test @testset "Salt" begin - testgrid = CryoGrid.Presets.DefaultGrid_2cm + testgrid = CryoGrid.DefaultGrid_2cm @testset "variables" begin soil = SalineGround() vars = CryoGrid.variables(soil) diff --git a/test/Physics/Surface/seb_tests.jl b/test/Physics/Surface/seb_tests.jl index c7762670..d4d96bee 100644 --- a/test/Physics/Surface/seb_tests.jl +++ b/test/Physics/Surface/seb_tests.jl @@ -9,8 +9,8 @@ Surface.surfaceproperties(seb::SurfaceEnergyBalance, sub::TestGroundLayer) = seb Sin = upreferred(200.0u"W/m^2") forcings = (; Tair, pr, qh, wind, Lin, Sin) @testset "Iterative" begin - seb = SurfaceEnergyBalance( - forcings; + seb = SurfaceEnergyBalance(; + forcings..., solscheme = Surface.Iterative(), ) # TODO: would be good to have SEB work correctly with units to verify correctness; @@ -29,8 +29,8 @@ Surface.surfaceproperties(seb::SurfaceEnergyBalance, sub::TestGroundLayer) = seb @test seb_output.Qg > zero(seb_output.Qg) end @testset "Analytical" begin - seb = SurfaceEnergyBalance( - forcings; + seb = SurfaceEnergyBalance(; + forcings..., solscheme = Surface.Analytical(), ) seb = pstrip(seb) @@ -46,8 +46,8 @@ Surface.surfaceproperties(seb::SurfaceEnergyBalance, sub::TestGroundLayer) = seb @test seb_output.Qg > zero(seb_output.Qg) end @testset "Numerical" begin - seb = SurfaceEnergyBalance( - forcings; + seb = SurfaceEnergyBalance(; + forcings..., solscheme=Surface.Numerical(), ) seb = pstrip(seb) diff --git a/test/Physics/Surface/swb_tests.jl b/test/Physics/Surface/swb_tests.jl index 6a524984..96cf4dbd 100644 --- a/test/Physics/Surface/swb_tests.jl +++ b/test/Physics/Surface/swb_tests.jl @@ -1,13 +1,8 @@ @testset "Surface Water Balance" begin - testgrid = CryoGrid.Presets.DefaultGrid_2cm + testgrid = CryoGrid.DefaultGrid_2cm @testset "interact!" begin - sub = TestGroundLayer(WaterBalance(BucketScheme())) - top = Top( - SurfaceWaterBalance( - ConstantForcing(1e-6u"m/s", :rainfall), - ConstantForcing(0.0u"m/s", :snowfall), - ) - ) + sub = stripparams(TestGroundLayer(WaterBalance(BucketScheme()))) + top = Top(SurfaceWaterBalance(1e-6u"m/s", 0.0u"m/s")) stop = Diagnostics.build_dummy_state(testgrid, top, with_units=true) ssub = Diagnostics.build_dummy_state(testgrid, sub, with_units=true) # initialize fully saturated diff --git a/test/Solvers/cglite_implicit_tests.jl b/test/Solvers/cglite_implicit_tests.jl index bc91216c..5c437e86 100644 --- a/test/Solvers/cglite_implicit_tests.jl +++ b/test/Solvers/cglite_implicit_tests.jl @@ -27,7 +27,7 @@ end α = upreferred(strat.soil.para.heat.kh_m) / upreferred(strat.soil.para.heat.ch_m) T_analytic = Heat.heat_conduction_linear_periodic_ub(T₀, A, P, ustrip(α)) initT = initializer(:T, (layer, state) -> state.T .= T_analytic.(cells(state.grid), 0.0)) - modelgrid = CryoGrid.Presets.DefaultGrid_2cm + modelgrid = CryoGrid.DefaultGrid_2cm # modelgrid = Grid(z_top:0.02u"m":z_bot) # modelgrid = Grid(vcat(0.0:0.02:1.0, 1.05:0.05:5.0, 5.1:0.1:10.0)*u"m") tile = Tile(strat, modelgrid, initT) @@ -58,7 +58,7 @@ end z_bot => Bottom(ConstantFlux(HeatBalance, 0.0)) ); initT = initializer(:T, -1.0) - modelgrid = CryoGrid.Presets.DefaultGrid_2cm + modelgrid = CryoGrid.DefaultGrid_2cm tile = Tile(strat, modelgrid, initT) # define time span, 5 years tspan = (0.0, 5*365*24*3600.0) @@ -68,7 +68,6 @@ end sol = solve(prob, LiteImplicitEuler(), dt=24*3600) out = CryoGridOutput(sol) - kh_w, kh_i, kh_a, kh_m, kh_o = Heat.thermalconductivities(strat.soil) θm = mineral(strat.soil) θo = organic(strat.soil) θ_s = (θw=0.0, θi=porosity(strat.soil), θa=0.0, θm=θm, θo=θo) diff --git a/test/benchmarks.jl b/test/benchmarks.jl index 87c01f4f..19e67e87 100644 --- a/test/benchmarks.jl +++ b/test/benchmarks.jl @@ -15,7 +15,7 @@ function run_solver_benchmark(solvers, grids, dts, tspan, upperbc; # high accuracy reference solution, 1 minute time steps println("Computing reference solution for grid spacing $min_dx with forward Euler @ $dt_ref s time steps ...") # basic 1-layer heat conduction model (defaults to free water freezing scheme) - model = Presets.SoilHeat(:H, upperbc, profile; grid=grid, freezecurve=freezecurve) + model = SoilHeat(:H, upperbc, profile; grid=grid, freezecurve=freezecurve) p = copy(model.pproto) p .= params prob = CryoGridProblem(model,tspan,p) @@ -55,7 +55,7 @@ forcings = loadforcings(filename); # use air temperature as upper boundary forcing tair = TemperatureBC(forcings.Tair) solvers = [Euler, DP5, ROCK2, ROCK4, Trapezoid, ROS3P] -grids = [Presets.DefaultGrid_2cm, Presets.DefaultGrid_5cm, Presets.DefaultGrid_10cm, Presets.DefaultGrid_20cm] +grids = [DefaultGrid_2cm, DefaultGrid_5cm, DefaultGrid_10cm, DefaultGrid_20cm] dts = [2*60.0, 10*60.0, 30*60.0, 3600.0] tspan = (DateTime(2010,10,30),DateTime(2011,10,30)) results_freeW_tair = run_solver_benchmark(solvers, grids, dts, tspan, tair) diff --git a/test/test_problems.jl b/test/test_problems.jl index d12036f2..93c81569 100644 --- a/test/test_problems.jl +++ b/test/test_problems.jl @@ -13,9 +13,9 @@ function test_heat_conduction_freeW_periodic_bc() initT = initializer(:T, tempprofile) # Define the simulation time span. tspan = (0.0,2*24*3600.0) - discretization = CryoGrid.Presets.DefaultGrid_10cm + discretization = CryoGrid.DefaultGrid_10cm Qbot = -0.01u"W/m^2" - tile = CryoGrid.Presets.SoilHeatTile( + tile = CryoGrid.SoilHeatTile( heatop, PeriodicBC(HeatBalance, CryoGrid.Neumann, 24*3600.0, 20.0, 0.0, -5.0), ConstantBC(HeatBalance, CryoGrid.Neumann, Qbot), @@ -32,7 +32,7 @@ function test_water_flow_bucket_scheme(topflux=0.0u"m/s") initsat = initializer(:sat, 0.5) # Define the simulation time span. tspan = (0.0, 24*3600.0) - discretization = CryoGrid.Presets.DefaultGrid_10cm + discretization = CryoGrid.DefaultGrid_10cm strat = Stratigraphy( 0.0u"m" => Top(ConstantBC(WaterBalance, CryoGrid.Neumann, topflux),), 0.0u"m" => Ground(SimpleSoil(por=0.5, sat=0.5), heat=nothing, water=WaterBalance(BucketScheme())), -- GitLab