Add getter for observed variables

docs/src/tutorials/acausal_components.md

@@ -320,6 +320,7 @@ sol = solve(prob)
 plot(sol)
 ```
 
+By default, this plots only the unknown variables that had to be solved for.
 However, what if we wanted to plot the timeseries of a different variable? Do
 not worry, that information was not thrown away! Instead, transformations
 like `structural_simplify` simply change unknown variables into observables which are
@@ -346,3 +347,9 @@ or we can plot the timeseries of the resistor's voltage:
 ```@example acausal
 plot(sol, idxs = [rc_model.resistor.v])
 ```
+
+Although it may be more confusing than helpful here, we can of course also plot all unknown and observed variables together:
+
+```@example acausal
+plot(sol, idxs = [unknowns(rc_model); observables(rc_model)])
+```

src/ModelingToolkit.jl

@@ -85,8 +85,8 @@ import Symbolics: rename, get_variables!, _solve, hessian_sparsity,
                   substituter, scalarize, getparent, hasderiv, hasdiff
 
 import DiffEqBase: @add_kwonly
-export independent_variables, unknowns, parameters, full_parameters, continuous_events,
-       discrete_events
+export independent_variables, unknowns, observables, parameters, full_parameters,
+       continuous_events, discrete_events
 @reexport using Symbolics
 @reexport using UnPack
 RuntimeGeneratedFunctions.init(@__MODULE__)

src/systems/abstractsystem.jl

@@ -430,7 +430,7 @@ function has_observed_with_lhs(sys, sym)
     if has_index_cache(sys) && (ic = get_index_cache(sys)) !== nothing
         return haskey(ic.observed_syms_to_timeseries, sym)
     else
-        return any(isequal(sym), [eq.lhs for eq in observed(sys)])
+        return any(isequal(sym), observables(sys))
     end
 end
 
@@ -489,7 +489,7 @@ function _all_ts_idxs!(ts_idxs, ::ScalarSymbolic, sys, sym::Symbol)
     if has_index_cache(sys) && (ic = get_index_cache(sys)) !== nothing
         return _all_ts_idxs!(ts_idxs, sys, ic.symbol_to_variable[sym])
     elseif is_variable(sys, sym) || is_independent_variable(sys, sym) ||
-           any(isequal(sym), [getname(eq.lhs) for eq in observed(sys)])
+           any(isequal(sym), getname.(observables(sys)))
         push!(ts_idxs, ContinuousTimeseries())
     elseif is_timeseries_parameter(sys, sym)
         push!(ts_idxs, timeseries_parameter_index(sys, s).timeseries_idx)
@@ -579,7 +579,7 @@ SymbolicIndexingInterface.constant_structure(::AbstractSystem) = true
 
 function SymbolicIndexingInterface.all_variable_symbols(sys::AbstractSystem)
     syms = variable_symbols(sys)
-    obs = getproperty.(observed(sys), :lhs)
+    obs = observables(sys)
     return isempty(obs) ? syms : vcat(syms, obs)
 end
 
@@ -1411,6 +1411,7 @@ _nonum(@nospecialize x) = x isa Num ? x.val : x
 $(TYPEDSIGNATURES)
 
 Get the unknown variables of the system `sys` and its subsystems.
+These must be explicitly solved for, unlike `observables(sys)`.
 
 See also [`ModelingToolkit.get_unknowns`](@ref).
 """
@@ -1677,6 +1678,14 @@ function controls(sys::AbstractSystem)
     isempty(systems) ? ctrls : [ctrls; reduce(vcat, namespace_controls.(systems))]
 end
 
+"""
+$(TYPEDSIGNATURES)
+
+Get the observed equations of the system `sys` and its subsystems.
+These can be expressed in terms of `unknowns(sys)`, and do not have to be explicitly solved for.
+
+See also [`observables`](@ref) and [`ModelingToolkit.get_observed()`](@ref).
+"""
 function observed(sys::AbstractSystem)
     obs = get_observed(sys)
     systems = get_systems(sys)
@@ -1686,6 +1695,19 @@ function observed(sys::AbstractSystem)
          init = Equation[])]
 end
 
+"""
+$(TYPEDSIGNATURES)
+
+Get the observed variables of the system `sys` and its subsystems.
+These can be expressed in terms of `unknowns(sys)`, and do not have to be explicitly solved for.
+It is equivalent to all left hand sides of `observed(sys)`.
+
+See also [`observed`](@ref).
+"""
+function observables(sys::AbstractSystem)
+    return map(eq -> eq.lhs, observed(sys))
+end
+
 Base.@deprecate default_u0(x) defaults(x) false
 Base.@deprecate default_p(x) defaults(x) false
 

src/systems/diffeqs/abstractodesystem.jl

@@ -1496,7 +1496,7 @@ function InitializationProblem{iip, specialize}(sys::AbstractSystem,
         @warn errmsg
     end
 
-    uninit = setdiff(unknowns(sys), [unknowns(isys); getfield.(observed(isys), :lhs)])
+    uninit = setdiff(unknowns(sys), [unknowns(isys); observables(isys)])
 
     # TODO: throw on uninitialized arrays
     filter!(x -> !(x isa Symbolics.Arr), uninit)

src/systems/diffeqs/odesystem.jl

@@ -47,7 +47,7 @@ struct ODESystem <: AbstractODESystem
     var_to_name::Any
     """Control parameters (some subset of `ps`)."""
     ctrls::Vector
-    """Observed variables."""
+    """Observed equations."""
     observed::Vector{Equation}
     """System of constraints that must be satisfied by the solution to the system."""
     constraintsystem::Union{Nothing, ConstraintsSystem}
@@ -532,7 +532,7 @@ function build_explicit_observed_function(sys, ts;
 
     vs = ModelingToolkit.vars(ts; op)
     namespace_subs = Dict()
-    ns_map = Dict{Any, Any}(renamespace(sys, eq.lhs) => eq.lhs for eq in observed(sys))
+    ns_map = Dict{Any, Any}(renamespace(sys, obs) => obs for obs in observables(sys))
     for sym in unknowns(sys)
         ns_map[renamespace(sys, sym)] = sym
         if iscall(sym) && operation(sym) === getindex

src/systems/diffeqs/sdesystem.jl

@@ -48,7 +48,7 @@ struct SDESystem <: AbstractODESystem
     var_to_name::Any
     """Control parameters (some subset of `ps`)."""
     ctrls::Vector
-    """Observed variables."""
+    """Observed equations."""
     observed::Vector{Equation}
     """
     Time-derivative matrix. Note: this field will not be defined until

src/systems/jumps/jumpsystem.jl

@@ -59,7 +59,7 @@ struct JumpSystem{U <: ArrayPartition} <: AbstractTimeDependentSystem
     ps::Vector
     """Array variables."""
     var_to_name::Any
-    """Observed variables."""
+    """Observed equations."""
     observed::Vector{Equation}
     """The name of the system."""
     name::Symbol

src/systems/nonlinear/nonlinearsystem.jl

@@ -32,7 +32,7 @@ struct NonlinearSystem <: AbstractTimeIndependentSystem
     ps::Vector
     """Array variables."""
     var_to_name::Any
-    """Observed variables."""
+    """Observed equations."""
     observed::Vector{Equation}
     """
     Jacobian matrix. Note: this field will not be defined until

src/systems/optimization/constraints_system.jl

@@ -33,7 +33,7 @@ struct ConstraintsSystem <: AbstractTimeIndependentSystem
     ps::Vector
     """Array variables."""
     var_to_name::Any
-    """Observed variables."""
+    """Observed equations."""
     observed::Vector{Equation}
     """
     Jacobian matrix. Note: this field will not be defined until

src/systems/optimization/optimizationsystem.jl

@@ -31,7 +31,7 @@ struct OptimizationSystem <: AbstractOptimizationSystem
     ps::Vector
     """Array variables."""
     var_to_name::Any
-    """Observed variables."""
+    """Observed equations."""
     observed::Vector{Equation}
     """List of constraint equations of the system."""
     constraints::Vector{Union{Equation, Inequality}}

src/systems/problem_utils.jl

@@ -1016,7 +1016,7 @@ function get_u0_p(sys,
         u0map = Dict(u0map)
     end
     if u0map isa Dict
-        allobs = Set(getproperty.(observed(sys), :lhs))
+        allobs = Set(observables(sys))
         if any(in(allobs), keys(u0map))
             u0s_in_obs = filter(in(allobs), keys(u0map))
             @warn "Observed variables cannot be assigned initial values. Initial values for $u0s_in_obs will be ignored."

src/systems/systemstructure.jl

@@ -721,7 +721,7 @@ function _structural_simplify!(state::TearingState, io; simplify = false,
         sys = ModelingToolkit.tearing(
             sys, state; simplify, mm, check_consistency, kwargs...)
     end
-    fullunknowns = [map(eq -> eq.lhs, observed(sys)); unknowns(sys)]
+    fullunknowns = [observables(sys); unknowns(sys)]
     @set! sys.observed = ModelingToolkit.topsort_equations(observed(sys), fullunknowns)
 
     ModelingToolkit.invalidate_cache!(sys), input_idxs

test/serialization.jl

@@ -31,7 +31,7 @@ sys = include_string(@__MODULE__, str)
 
 # check answer
 ss = structural_simplify(rc_model)
-all_obs = [o.lhs for o in observed(ss)]
+all_obs = observables(ss)
 prob = ODEProblem(ss, [capacitor.v => 0.0], (0, 0.1))
 sol = solve(prob, ImplicitEuler())
 

test/structural_transformation/utils.jl

@@ -51,8 +51,8 @@ end
         [D(x) ~ z[1] + z[2] + foo(z)[1], y[1] ~ 2t, y[2] ~ 3t, z ~ foo(y)], t)
     @test length(equations(sys)) == 1
     @test length(observed(sys)) == 7
-    @test any(eq -> isequal(eq.lhs, y), observed(sys))
-    @test any(eq -> isequal(eq.lhs, z), observed(sys))
+    @test any(obs -> isequal(obs, y), observables(sys))
+    @test any(obs -> isequal(obs, z), observables(sys))
     prob = ODEProblem(sys, [x => 1.0], (0.0, 1.0), [foo => _tmp_fn])
     @test_nowarn prob.f(prob.u0, prob.p, 0.0)