Merge pull request #2862 from hersle/check_iv_is_par

Create independent variables with @independent_variables

docs/src/basics/FAQ.md

@@ -186,7 +186,7 @@ p, replace, alias = SciMLStructures.canonicalize(Tunable(), prob.p)
 
 This error can come up after running `structural_simplify` on a system that generates dummy derivatives (i.e. variables with `ˍt`).  For example, here even though all the variables are defined with initial values, the `ODEProblem` generation will throw an error that defaults are missing from the variable map.
 
-```
+```julia
 using ModelingToolkit
 using ModelingToolkit: t_nounits as t, D_nounits as D
 
@@ -197,13 +197,13 @@ eqs = [x1 + x2 + 1 ~ 0
        2 * D(D(x1)) + D(D(x2)) + D(D(x3)) + D(x4) + 4 ~ 0]
 @named sys = ODESystem(eqs, t)
 sys = structural_simplify(sys)
-prob = ODEProblem(sys, [], (0,1))
+prob = ODEProblem(sys, [], (0, 1))
 ```
 
 We can solve this problem by using the `missing_variable_defaults()` function
 
-```
-prob = ODEProblem(sys, ModelingToolkit.missing_variable_defaults(sys), (0,1))
+```julia
+prob = ODEProblem(sys, ModelingToolkit.missing_variable_defaults(sys), (0, 1))
 ```
 
 This function provides 0 for the default values, which is a safe assumption for dummy derivatives of most models.  However, the 2nd argument allows for a different default value or values to be used if needed.
@@ -221,12 +221,26 @@ julia> ModelingToolkit.missing_variable_defaults(sys, [1,2,3])
 Use the `u0_constructor` keyword argument to map an array to the desired
 container type. For example:
 
-```
+```julia
 using ModelingToolkit, StaticArrays
 using ModelingToolkit: t_nounits as t, D_nounits as D
 
-sts = @variables x1(t)=0.0
+sts = @variables x1(t) = 0.0
 eqs = [D(x1) ~ 1.1 * x1]
 @mtkbuild sys = ODESystem(eqs, t)
-prob = ODEProblem{false}(sys, [], (0,1); u0_constructor = x->SVector(x...))
+prob = ODEProblem{false}(sys, [], (0, 1); u0_constructor = x -> SVector(x...))
+```
+
+## Using a custom independent variable
+
+When possible, we recommend `using ModelingToolkit: t_nounits as t, D_nounits as D` as the independent variable and its derivative.
+However, if you want to use your own, you can do so:
+
+```julia
+using ModelingToolkit
+
+@independent_variables x
+D = Differential(x)
+@variables y(x)
+@named sys = ODESystem([D(y) ~ x], x)
 ```

docs/src/basics/Validation.md

@@ -8,15 +8,19 @@ Units may be assigned with the following syntax.
 
 ```@example validation
 using ModelingToolkit, DynamicQuantities
-@variables t [unit = u"s"] x(t) [unit = u"m"] g(t) w(t) [unit = "Hz"]
+@independent_variables t [unit = u"s"]
+@variables x(t) [unit = u"m"] g(t) w(t) [unit = u"Hz"]
 
-@variables(t, [unit = u"s"], x(t), [unit = u"m"], g(t), w(t), [unit = "Hz"])
+@parameters(t, [unit = u"s"])
+@variables(x(t), [unit = u"m"], g(t), w(t), [unit = u"Hz"])
 
+@parameters begin
+    t, [unit = u"s"]
+end
 @variables(begin
-    t, [unit = u"s"],
     x(t), [unit = u"m"],
     g(t),
-    w(t), [unit = "Hz"]
+    w(t), [unit = u"Hz"]
 end)
 
 # Simultaneously set default value (use plain numbers, not quantities)
@@ -46,10 +50,11 @@ Example usage below. Note that `ModelingToolkit` does not force unit conversions
 
 ```@example validation
 using ModelingToolkit, DynamicQuantities
+@independent_variables t [unit = u"ms"]
 @parameters τ [unit = u"ms"]
-@variables t [unit = u"ms"] E(t) [unit = u"kJ"] P(t) [unit = u"MW"]
+@variables E(t) [unit = u"kJ"] P(t) [unit = u"MW"]
 D = Differential(t)
-eqs = eqs = [D(E) ~ P - E / τ,
+eqs = [D(E) ~ P - E / τ,
     0 ~ P]
 ModelingToolkit.validate(eqs)
 ```
@@ -70,10 +75,11 @@ An example of an inconsistent system: at present, `ModelingToolkit` requires tha
 
 ```@example validation
 using ModelingToolkit, DynamicQuantities
+@independent_variables t [unit = u"ms"]
 @parameters τ [unit = u"ms"]
-@variables t [unit = u"ms"] E(t) [unit = u"J"] P(t) [unit = u"MW"]
+@variables E(t) [unit = u"J"] P(t) [unit = u"MW"]
 D = Differential(t)
-eqs = eqs = [D(E) ~ P - E / τ,
+eqs = [D(E) ~ P - E / τ,
     0 ~ P]
 ModelingToolkit.validate(eqs) #Returns false while displaying a warning message
 ```
@@ -115,7 +121,8 @@ In order for a function to work correctly during both validation & execution, th
 
 ```julia
 using ModelingToolkit, DynamicQuantities
-@variables t [unit = u"ms"] E(t) [unit = u"J"] P(t) [unit = u"MW"]
+@independent_variables t [unit = u"ms"]
+@variables E(t) [unit = u"J"] P(t) [unit = u"MW"]
 D = Differential(t)
 eqs = [D(E) ~ P - E / 1u"ms"]
 ModelingToolkit.validate(eqs) #Returns false while displaying a warning message
@@ -129,8 +136,9 @@ Instead, they should be parameterized:
 
 ```@example validation3
 using ModelingToolkit, DynamicQuantities
+@independent_variables t [unit = u"ms"]
 @parameters τ [unit = u"ms"]
-@variables t [unit = u"ms"] E(t) [unit = u"kJ"] P(t) [unit = u"MW"]
+@variables E(t) [unit = u"kJ"] P(t) [unit = u"MW"]
 D = Differential(t)
 eqs = [D(E) ~ P - E / τ]
 ModelingToolkit.validate(eqs) #Returns true

docs/src/tutorials/SampledData.md

@@ -97,7 +97,8 @@ H(z) = \dfrac{b_2 z^2 + b_1 z + b_0}{a_2 z^2 + a_1 z + a_0}
 may thus be modeled as
 
 ```julia
-@variables t y(t) [description = "Output"] u(t) [description = "Input"]
+t = ModelingToolkit.t_nounits
+@variables y(t) [description = "Output"] u(t) [description = "Input"]
 k = ShiftIndex(Clock(t, dt))
 eqs = [
     a2 * y(k) + a1 * y(k - 1) + a0 * y(k - 2) ~ b2 * u(k) + b1 * u(k - 1) + b0 * u(k - 2)

src/ModelingToolkit.jl

@@ -127,6 +127,7 @@ using .BipartiteGraphs
 
 include("variables.jl")
 include("parameters.jl")
+include("independent_variables.jl")
 include("constants.jl")
 
 include("utils.jl")
@@ -183,13 +184,13 @@ for S in subtypes(ModelingToolkit.AbstractSystem)
 end
 
 const t_nounits = let
-    only(@parameters t)
+    only(@independent_variables t)
 end
 const t_unitful = let
-    only(@parameters t [unit = Unitful.u"s"])
+    only(@independent_variables t [unit = Unitful.u"s"])
 end
 const t = let
-    only(@parameters t [unit = DQ.u"s"])
+    only(@independent_variables t [unit = DQ.u"s"])
 end
 
 const D_nounits = Differential(t_nounits)
@@ -262,7 +263,7 @@ export generate_initializesystem
 export alg_equations, diff_equations, has_alg_equations, has_diff_equations
 export get_alg_eqs, get_diff_eqs, has_alg_eqs, has_diff_eqs
 
-export @variables, @parameters, @constants, @brownian
+export @variables, @parameters, @independent_variables, @constants, @brownian
 export @named, @nonamespace, @namespace, extend, compose, complete
 export debug_system
 

src/discretedomain.jl

@@ -98,7 +98,7 @@ $(FIELDS)
 ```jldoctest
 julia> using Symbolics
 
-julia> @variables t;
+julia> t = ModelingToolkit.t_nounits
 
 julia> Δ = Sample(t, 0.01)
 (::Sample) (generic function with 2 methods)
@@ -166,7 +166,9 @@ The `ShiftIndex` operator allows you to index a signal and obtain a shifted disc
 # Examples
 
 ```
-julia> @variables t x(t);
+julia> t = ModelingToolkit.t_nounits;
+
+julia> @variables x(t);
 
 julia> k = ShiftIndex(t, 0.1);
 

src/independent_variables.jl

@@ -0,0 +1,11 @@
+"""
+    @independent_variables t₁ t₂ ...
+
+Define one or more independent variables. For example:
+
+    @independent_variables t
+    @variables x(t)
+"""
+macro independent_variables(ts...)
+    :(@parameters $(ts...)) |> esc # TODO: treat independent variables separately from variables and parameters
+end

src/systems/abstractsystem.jl

@@ -2237,11 +2237,10 @@ This example builds the following feedback interconnection and linearizes it fro
 
 ```julia
 using ModelingToolkit
-@variables t
+using ModelingToolkit: t_nounits as t, D_nounits as D
 function plant(; name)
     @variables x(t) = 1
     @variables u(t)=0 y(t)=0
-    D = Differential(t)
     eqs = [D(x) ~ -x + u
            y ~ x]
     ODESystem(eqs, t; name = name)
@@ -2250,7 +2249,6 @@ end
 function ref_filt(; name)
     @variables x(t)=0 y(t)=0
     @variables u(t)=0 [input = true]
-    D = Differential(t)
     eqs = [D(x) ~ -2 * x + u
            y ~ x]
     ODESystem(eqs, t, name = name)

src/systems/alias_elimination.jl

@@ -388,8 +388,10 @@ Use Kahn's algorithm to topologically sort observed equations.
 
 Example:
 ```julia
-julia> @variables t x(t) y(t) z(t) k(t)
-(t, x(t), y(t), z(t), k(t))
+julia> t = ModelingToolkit.t_nounits
+
+julia> @variables x(t) y(t) z(t) k(t)
+(x(t), y(t), z(t), k(t))
 
 julia> eqs = [
            x ~ y + z

src/systems/dependency_graphs.jl

@@ -15,8 +15,9 @@ Example:
 
 ```julia
 using ModelingToolkit
+using ModelingToolkit: t_nounits as t
 @parameters β γ κ η
-@variables t S(t) I(t) R(t)
+@variables S(t) I(t) R(t)
 
 rate₁ = β * S * I
 rate₂ = γ * I + t

src/systems/diffeqs/basic_transformations.jl

@@ -16,7 +16,8 @@ Example:
 ```julia
 using ModelingToolkit, OrdinaryDiffEq, Test
 
-@parameters t α β γ δ
+@independent_variables t
+@parameters α β γ δ
 @variables x(t) y(t)
 D = Differential(t)
 

src/systems/diffeqs/modelingtoolkitize.jl

@@ -242,7 +242,7 @@ Generate `SDESystem`, dependent variables, and parameters from an `SDEProblem`.
 function modelingtoolkitize(prob::DiffEqBase.SDEProblem; kwargs...)
     prob.f isa DiffEqBase.AbstractParameterizedFunction &&
         return (prob.f.sys, prob.f.sys.unknowns, prob.f.sys.ps)
-    @parameters t
+    @independent_variables t
     p = prob.p
     has_p = !(p isa Union{DiffEqBase.NullParameters, Nothing})
 

src/systems/diffeqs/odesystem.jl

@@ -10,10 +10,10 @@ $(FIELDS)
 
 ```julia
 using ModelingToolkit
+using ModelingToolkit: t_nounits as t, D_nounits as D
 
 @parameters σ ρ β
-@variables t x(t) y(t) z(t)
-D = Differential(t)
+@variables x(t) y(t) z(t)
 
 eqs = [D(x) ~ σ*(y-x),
        D(y) ~ x*(ρ-z)-y,
@@ -184,6 +184,7 @@ struct ODESystem <: AbstractODESystem
             discrete_subsystems = nothing, solved_unknowns = nothing,
             split_idxs = nothing, parent = nothing; checks::Union{Bool, Int} = true)
         if checks == true || (checks & CheckComponents) > 0
+            check_independent_variables([iv])
             check_variables(dvs, iv)
             check_parameters(ps, iv)
             check_equations(deqs, iv)

src/systems/diffeqs/sdesystem.jl

@@ -10,10 +10,10 @@ $(FIELDS)
 
 ```julia
 using ModelingToolkit
+using ModelingToolkit: t_nounits as t, D_nounits as D
 
 @parameters σ ρ β
-@variables t x(t) y(t) z(t)
-D = Differential(t)
+@variables x(t) y(t) z(t)
 
 eqs = [D(x) ~ σ*(y-x),
        D(y) ~ x*(ρ-z)-y,
@@ -137,6 +137,7 @@ struct SDESystem <: AbstractODESystem
             complete = false, index_cache = nothing, parent = nothing;
             checks::Union{Bool, Int} = true)
         if checks == true || (checks & CheckComponents) > 0
+            check_independent_variables([iv])
             check_variables(dvs, iv)
             check_parameters(ps, iv)
             check_equations(deqs, iv)
@@ -320,10 +321,10 @@ experiments. Springer Science & Business Media.
 
 ```julia
 using ModelingToolkit
+using ModelingToolkit: t_nounits as t, D_nounits as D
 
 @parameters α β
-@variables t x(t) y(t) z(t)
-D = Differential(t)
+@variables x(t) y(t) z(t)
 
 eqs = [D(x) ~ α*x]
 noiseeqs = [β*x]

src/systems/discrete_system/discrete_system.jl

@@ -101,6 +101,7 @@ struct DiscreteSystem <: AbstractTimeDependentSystem
             complete = false, index_cache = nothing, parent = nothing;
             checks::Union{Bool, Int} = true)
         if checks == true || (checks & CheckComponents) > 0
+            check_independent_variables([iv])
             check_variables(dvs, iv)
             check_parameters(ps, iv)
         end

src/systems/jumps/jumpsystem.jl

@@ -34,9 +34,10 @@ $(FIELDS)
 
 ```julia
 using ModelingToolkit, JumpProcesses
+using ModelingToolkit: t_nounits as t
 
 @parameters β γ
-@variables t S(t) I(t) R(t)
+@variables S(t) I(t) R(t)
 rate₁   = β*S*I
 affect₁ = [S ~ S - 1, I ~ I + 1]
 rate₂   = γ*I
@@ -118,6 +119,7 @@ struct JumpSystem{U <: ArrayPartition} <: AbstractTimeDependentSystem
             complete = false, index_cache = nothing;
             checks::Union{Bool, Int} = true) where {U <: ArrayPartition}
         if checks == true || (checks & CheckComponents) > 0
+            check_independent_variables([iv])
             check_variables(unknowns, iv)
             check_parameters(ps, iv)
         end

src/systems/pde/pdesystem.jl

@@ -11,8 +11,8 @@ $(FIELDS)
 ```julia
 using ModelingToolkit
 
-@parameters x
-@variables t u(..)
+@parameters x t
+@variables u(..)
 Dxx = Differential(x)^2
 Dtt = Differential(t)^2
 Dt = Differential(t)

src/utils.jl

@@ -113,6 +113,13 @@ const CheckAll = 1 << 0
 const CheckComponents = 1 << 1
 const CheckUnits = 1 << 2
 
+function check_independent_variables(ivs)
+    for iv in ivs
+        isparameter(iv) ||
+            @warn "Independent variable $iv should be defined with @independent_variables $iv."
+    end
+end
+
 function check_parameters(ps, iv)
     for p in ps
         isequal(iv, p) &&
@@ -349,7 +356,8 @@ Return a `Set` containing all variables in `x` that appear in
 Example:
 
 ```
-@variables t u(t) y(t)
+t = ModelingToolkit.t_nounits
+@variables u(t) y(t)
 D  = Differential(t)
 v  = ModelingToolkit.vars(D(y) ~ u)
 v == Set([D(y), u])
@@ -421,7 +429,8 @@ collect_differential_variables(sys) = collect_operator_variables(sys, Differenti
 Return  a `Set` with all applied operators in `x`, example:
 
 ```
-@variables t u(t) y(t)
+@independent_variables t
+@variables u(t) y(t)
 D = Differential(t)
 eq = D(y) ~ u
 ModelingToolkit.collect_applied_operators(eq, Differential) == Set([D(y)])