Give the reduced system its own time span

Project.toml

@@ -17,7 +17,7 @@ TSVD = "9449cd9e-2762-5aa3-a617-5413e99d722e"
 [compat]
 DocStringExtensions = "0.9.5"
 LinearAlgebra = "1"
-ModelingToolkit = "11.39.1"
+ModelingToolkit = "11.42"
 PrecompileTools = "1.2.1"
 RandomizedLinAlg = "0.1"
 SciMLBase = "3.48"

README.md

@@ -68,7 +68,7 @@ sol = solve(full_prob)
 # apply POD-DEIM to obtain the reduced-order model
 pod_dim = deim_dim = 5
 deim_sys = deim(full_prob, sol, pod_dim; deim_dim = deim_dim)
-deim_prob = DAEProblem(deim_sys, nothing, full_prob.tspan; build_initializeprob = false)
+deim_prob = DAEProblem(deim_sys, nothing; build_initializeprob = false)
 deim_sol = solve(deim_prob)
 
 # retrieve the approximate solution of the original full-order model

docs/src/tutorials/deim_FitzHugh-Nagumo.md

@@ -138,7 +138,7 @@ terms. This can be done by simply calling [`deim`](@ref).
 using ModelOrderReduction
 pod_dim = deim_dim = 5
 deim_sys = deim(full_prob, sol, pod_dim; deim_dim = deim_dim)
-deim_prob = DAEProblem(deim_sys, nothing, full_prob.tspan; build_initializeprob = false)
+deim_prob = DAEProblem(deim_sys, nothing; build_initializeprob = false)
 deim_sol = solve(deim_prob)
 nₜ_deim = length(deim_sol[t])
 sol_deim_x = deim_sol[x]

src/deim.jl

@@ -258,7 +258,7 @@ the first snapshot column, whose time is `first(times)` when `times` is given.
 """
 function _reduced_system(
         fom::FullOrderModel, snapshot::AbstractMatrix, V::AbstractMatrix, U::AbstractMatrix,
-        name::Symbol; times = nothing, kwargs...
+        name::Symbol; times = nothing, tspan = nothing, kwargs...
     )
     iv = fom.iv
     D = Differential(iv)
@@ -305,7 +305,7 @@ function _reduced_system(
     reduced = System(
         [D(reduced_state) ~ rhs], iv, Symbolics.unwrap.(Symbolics.scalarize(reduced_state)),
         [fom.parameters; array_parameters; first(lift_value)];
-        name, observed = [reconstruction; preserved], initial_conditions
+        name, observed = [reconstruction; preserved], initial_conditions, tspan
     )
     for key in (ModelingToolkit.ProblemTypeCtx, ModelingToolkit.MiscSystemData)
         SymbolicUtils.hasmetadata(fom.system, key) || continue
@@ -319,7 +319,7 @@ end
 function _deim(
         sys::System, snapshot::AbstractMatrix, pod_dim::Integer, deim_dim::Integer,
         name::Symbol; snapshot_times = nothing, training_parameters = Dict{Any, Any}(),
-        kwargs...
+        tspan = nothing, kwargs...
     )
     rows = length(ModelingToolkit.unknowns(sys))
     size(snapshot, 1) == rows || throw(
@@ -339,7 +339,8 @@ function _deim(
     )
     nonlinear_basis = _pod_basis(nonlinear_snapshot, deim_dim)
     return _reduced_system(
-        fom, snapshot, state_basis, nonlinear_basis, name; times = snapshot_times, kwargs...
+        fom, snapshot, state_basis, nonlinear_basis, name;
+        times = snapshot_times, tspan, kwargs...
     )
 end
 
@@ -379,6 +380,10 @@ Nonlinear terms are evaluated at the numeric parameter defaults of `sys`. Use th
 problem/solution method for problem-specific parameter values. The reduced state and its
 derivative at the first snapshot column are stored as initial conditions.
 
+The reduced system inherits the time span of `sys`, so a problem can be built from it
+without repeating one. Passing a time span explicitly still overrides it, and a source
+system without one produces a reduced system without one.
+
 Construct a `DAEProblem` with `build_initializeprob = false` from the returned system
 to keep array code generation, or call `ModelingToolkit.mtkcompile` on it and construct
 an `ODEProblem` to generate scalar code for the reduced equation. Reconstruct fields with
@@ -415,7 +420,10 @@ function deim(
         deim_dim::Integer = pod_dim, name::Symbol = Symbol(nameof(sys), :_deim),
         snapshot_times = nothing, kwargs...
     )::System
-    return _deim(sys, snapshot, pod_dim, deim_dim, name; snapshot_times, kwargs...)
+    return _deim(
+        sys, snapshot, pod_dim, deim_dim, name;
+        snapshot_times, tspan = ModelingToolkit.get_tspan(sys), kwargs...
+    )
 end
 
 const SymbolicProblem = Union{SciMLBase.AbstractODEProblem, SciMLBase.AbstractDAEProblem}
@@ -437,6 +445,8 @@ Reduce the symbolic system behind `prob` with POD-DEIM, using one of its saved s
 `DAEProblem` returned by MethodOfLines. The saved states are the snapshot columns, the
 saved times supply the independent variable for time-dependent terms, and the parameter
 values of `prob` are used for training and become the defaults of the reduced system.
+The reduced system also inherits the time span of `prob`, which is the interval it was
+trained on, so a problem can be built from it without repeating one.
 `sol` may be the `SciMLBase.PDETimeSeriesSolution` returned by MethodOfLines or its
 underlying `SciMLBase.AbstractODESolution`. See the system method for the reduction
 itself and for how to construct problems from the returned system.
@@ -466,9 +476,7 @@ itself and for how to construct problems from the returned system.
 full_problem = discretize(pde_system, discretization; fallback = false)
 full_solution = solve(full_problem)
 reduced_system = deim(full_problem, full_solution, 4)
-reduced_problem = DAEProblem(
-    reduced_system, nothing, full_problem.tspan; build_initializeprob = false
-)
+reduced_problem = DAEProblem(reduced_system, nothing; build_initializeprob = false)
 ```
 """
 function deim(
@@ -494,7 +502,7 @@ function deim(
     reduced_name = isnothing(name) ? Symbol(nameof(sys), :_deim) : name
     return _deim(
         sys, snapshot, pod_dim, deim_dim, reduced_name;
-        snapshot_times = sol.t, training_parameters, kwargs...
+        snapshot_times = sol.t, training_parameters, tspan = prob.tspan, kwargs...
     )
 end
 

test/Project.toml

@@ -13,7 +13,7 @@ Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 SymbolicIndexingInterface = "0.3"
 LinearAlgebra = "1"
 MethodOfLines = "1"
-ModelingToolkit = "11"
+ModelingToolkit = "11.42"
 OrdinaryDiffEq = "7"
 SafeTestsets = "0.1"
 SciMLBase = "3"

test/deim.jl

@@ -366,3 +366,44 @@ end
     )
     @test norm(reconstruction - truth) / norm(truth) < 1.0e-6
 end
+
+@testset "the reduced system carries a tspan" begin
+    @variables p(t)[1:6]
+    p0 = [0.4, 0.3, 0.2, 0.1, -0.2, 0.5]
+    pstates = ModelingToolkit.Symbolics.unwrap.(ModelingToolkit.Symbolics.scalarize(p))
+    span = (0.0, 0.5)
+    @named source = System(
+        [Dt(p) ~ -p - p .^ 3], t, pstates, [];
+        initial_conditions = [p => p0, Dt(p) => -p0 - p0 .^ 3],
+    )
+    source_problem = DAEProblem(complete(source), nothing, span; build_initializeprob = false)
+    source_solution = solve(source_problem; saveat = 0.05, abstol = 1.0e-10, reltol = 1.0e-10)
+    @test successful_retcode(source_solution)
+
+    # The problem method trains on `prob`, so its interval is the reduced model's.
+    from_problem = deim(source_problem, source_solution, 3)
+    @test ModelingToolkit.get_tspan(from_problem) == span
+    reduced_problem = DAEProblem(from_problem, nothing; build_initializeprob = false)
+    @test reduced_problem.tspan == span
+    reduced_solution = solve(reduced_problem; saveat = 0.05, abstol = 1.0e-10, reltol = 1.0e-10)
+    @test successful_retcode(reduced_solution)
+    truth = reduce(hcat, source_solution[p])
+    reconstruction = reduce(
+        hcat, [
+            SII.observed(from_problem, p)(state, reduced_problem.p, time)
+                for (state, time) in zip(reduced_solution.u, reduced_solution.t)
+        ]
+    )
+    @test norm(reconstruction - truth) / norm(truth) < 1.0e-4
+
+    # The system method inherits whatever the source system stores.
+    snapshot = [sin(0.3 * i * j) + cos(0.2 * i * (j + 1)) for i in 1:6, j in 1:8]
+    @mtkcompile with_span = System([Dt(p) ~ -p - p .^ 3], t; tspan = span)
+    @test ModelingToolkit.get_tspan(deim(with_span, snapshot, 2)) == span
+
+    # A source system without one must not acquire a fabricated interval.
+    @mtkcompile without_span = System([Dt(p) ~ -p - p .^ 3], t)
+    bare = deim(without_span, snapshot, 2)
+    @test ModelingToolkit.get_tspan(bare) === nothing
+    @test DAEProblem(bare, nothing, span; build_initializeprob = false).tspan == span
+end