diff --git a/Project.toml b/Project.toml index 8a45647e3add22cf76a4b8433ab1ce1c8a0bdba6..e6a046ddd90b042c35a158e52bba7b77d6590a36 100755 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "CryoGrid" uuid = "a535b82e-5f3d-4d97-8b0b-d6483f5bebd5" authors = ["Brian Groenke <brian.groenke@awi.de>", "Jan Nitzbon <jan.nitzbon@awi.de>", "Moritz Langer <moritz.langer@awi.de>"] -version = "0.23.1" +version = "0.23.2" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" diff --git a/examples/heat_sfcc_samoylov_implicitT.jl b/examples/heat_sfcc_samoylov_implicitT.jl index 469d71c3b4afe8299e5ff1562e06e2851924d87c..6c355aa3fbb75ba0aabe9ec298a45d22e3bfce29 100644 --- a/examples/heat_sfcc_samoylov_implicitT.jl +++ b/examples/heat_sfcc_samoylov_implicitT.jl @@ -38,7 +38,7 @@ tile = Tile(strat, modelgrid, forcings, ssinit); tspan = (DateTime(2000,10,1), DateTime(2002,10,1)) u0, du0 = initialcondition!(tile, tspan); prob = CryoGridProblem(tile, u0, tspan, saveat=24*3600, savevars=(:T,:θw,:∂H∂T), step_limiter=nothing) -integrator = init(prob, ImplicitEuler(autodiff=false)) +integrator = init(prob, ImplicitEuler()) # test one step @time step!(integrator) last_t = integrator.t diff --git a/src/Solvers/LiteImplicit/cglite_types.jl b/src/Solvers/LiteImplicit/cglite_types.jl index 47c993f91d840839ecd1e45eca97049199fb59df..11d3253f61d7d3306a5b073c7bea643ef39f1bde 100644 --- a/src/Solvers/LiteImplicit/cglite_types.jl +++ b/src/Solvers/LiteImplicit/cglite_types.jl @@ -22,7 +22,16 @@ struct LiteImplicitEulerCache{Tu,TA} <: SciMLBase.DECache D::TA end -function DiffEqBase.__init(prob::CryoGridProblem, alg::LiteImplicitEuler, args...; dt=24*3600.0, saveat=dt, kwargs...) +function DiffEqBase.__init( + prob::CryoGridProblem, + alg::LiteImplicitEuler, + args...; + dt=24*3600.0, + dtmax=dt, + dtmin=1.0, + saveat=dt, + kwargs... +) tile = Tile(prob.f) grid = tile.grid u0 = copy(prob.u0) @@ -54,6 +63,6 @@ function DiffEqBase.__init(prob::CryoGridProblem, alg::LiteImplicitEuler, args.. similar(u0, eltype(u0), length(prob.u0.H)), ) p = isnothing(prob.p) ? prob.p : collect(prob.p) - opts = CryoGridIntegratorOptions(; dtmax=dt, saveat=CryoGrid.expandtstep(saveat, prob.tspan), kwargs...) + opts = CryoGridIntegratorOptions(; saveat=CryoGrid.expandtstep(saveat, prob.tspan), dtmax, dtmin, kwargs...) return CryoGridIntegrator(alg, cache, opts, sol, copy(u0), p, t0, convert(eltype(prob.tspan), dt), 1, 1) end diff --git a/src/Solvers/basic_solvers.jl b/src/Solvers/basic_solvers.jl index 5e9afe2d97de1aafaaf3df2672a5d26fccf7dd92..c03a4b99f970f39324d6195d003be0c7f7af9429 100644 --- a/src/Solvers/basic_solvers.jl +++ b/src/Solvers/basic_solvers.jl @@ -74,6 +74,6 @@ function perform_step!(integrator::CryoGridIntegrator{CGEuler}) # set next dt # compute maximum timestep dtmax = CryoGrid.timestep(tile, du, u, p, t₀) - integrator.dt = dtmax + integrator.dt = max(min(dtmax, integrator.opts.dtmax), integrator.opts.dtmin) return nothing end diff --git a/src/Solvers/integrator.jl b/src/Solvers/integrator.jl index 53ed2732c90f03400e633e5a5d4d73c9f2965111..195bab2b42834c3a48ab95f8c80e9901803820a1 100755 --- a/src/Solvers/integrator.jl +++ b/src/Solvers/integrator.jl @@ -52,7 +52,7 @@ function DiffEqBase.sensitivity_solution(sol::CryoGridSolution, u, t) end Base.@kwdef mutable struct CryoGridIntegratorOptions{Tt,Tsaveat,Ttstops} - dtmax::Tt = 24*3600.0 + dtmax::Tt = 3600.0 dtmin::Tt = one(typeof(dtmax)) saveat::Tsaveat = typeof(dtmax)[] tstops::Ttstops = SortedSet{typeof(dtmin)}() diff --git a/src/problem.jl b/src/problem.jl index 6c05881f386d286a86159e3d4c921dab5e90d1f9..db3ebf1952e624598185d28800507f4d219be102 100644 --- a/src/problem.jl +++ b/src/problem.jl @@ -47,7 +47,8 @@ Constructor for `CryoGridProblem` that automatically generates all necessary cal function CryoGridProblem( tile::Tile, u0::ComponentVector, - tspan::NTuple{2,Float64}; + tspan::NTuple{2,Float64}, + p::Union{Nothing,AbstractVector}=nothing; diagnostic_stepsize=3600.0, saveat=3600.0, savevars=(), @@ -68,7 +69,7 @@ function CryoGridProblem( # strip all "fixed" parameters tile = stripparams(FixedParam, tile) # retrieve variable parameters - p = length(ModelParameters.params(tile)) > 0 ? parameters(tile) : nothing + tilepara = length(ModelParameters.params(tile)) > 0 ? parameters(tile) : nothing du0 = zero(u0) # remove units tile = stripunits(tile) @@ -95,6 +96,8 @@ function CryoGridProblem( tile.data.outputs = savevals # build mass matrix mass_matrix = Numerics.build_mass_matrix(tile.state) + # get params + p = isnothing(p) && !isnothing(tilepara) ? ustrip.(vec(tilepara)) : p func = odefunction(tile, u0, p, tspan; mass_matrix, specialization, function_kwargs...) return CryoGridProblem{true}(func, u0, tspan, p, callbacks, saveat, getsavestate, isoutofdomain, prob_kwargs) end