add support for symbolic parameters

Author: pepijndevos

src/systems/abstractsystem.jl

@@ -181,6 +181,15 @@ function Base.getproperty(sys::AbstractSystem, name::Symbol)
     throw(error("Variable $name does not exist"))
 end
 
+function Base.setproperty!(sys::AbstractSystem, prop::Symbol, val)
+    param = Sym{Parameter{Real}}(prop)
+    if param in parameters(sys)
+        sys.default_p[param] = val
+    else
+        setfield!(sys, prop, val)
+    end
+end
+
 function renamespace(namespace, x)
     if x isa Num
         renamespace(namespace, value(x))
@@ -203,7 +212,7 @@ end
 
 function namespace_default_p(sys)
     d_p = default_p(sys)
-    Dict(parameters(sys, k) => d_p[k] for k in keys(d_p))
+    Dict(parameters(sys, k) => namespace_expr(d_p[k], nameof(sys), []) for k in keys(d_p))
 end
 
 function namespace_equations(sys::AbstractSystem)

src/variables.jl

@@ -213,6 +213,10 @@ applicable.
 function varmap_to_vars(varmap::Dict, varlist; defaults=Dict())
     varmap = merge(defaults, varmap) # prefers the `varmap`
     varmap = Dict(value(k)=>value(varmap[k]) for k in keys(varmap))
+    # resolve symbolic parameter expressions
+    for (p, v) in pairs(varmap)
+        varmap[p] = fixpoint_sub(v, varmap)
+    end
     T′ = eltype(values(varmap))
     T = Base.isconcretetype(T′) ? T′ : Base.promote_typeof(values(varmap)...)
     out = Vector{T}(undef, length(varlist))

test/runtests.jl

@@ -1,5 +1,6 @@
 using SafeTestsets, Test
 
+@safetestset "Symbolic parameters test" begin include("symbolic_parameters.jl") end
 @safetestset "Parsing Test" begin include("variable_parsing.jl") end
 @safetestset "Differentiation Test" begin include("derivatives.jl") end
 @safetestset "Simplify Test" begin include("simplify.jl") end

test/symbolic_parameters.jl

@@ -0,0 +1,35 @@
+using ModelingToolkit
+using NonlinearSolve
+using Test
+
+@variables x y z
+@parameters σ ρ β
+
+eqs = [0 ~ σ*(y-x),
+       0 ~ x*(ρ-z)-y,
+       0 ~ x*y - β*z]
+
+par = [
+    σ => 1,
+    ρ => 0.1+σ,
+    β => ρ*1.1
+]
+ns = NonlinearSystem(eqs, [x,y,z],[σ,ρ,β], name=:ns, default_p=par)
+ModelingToolkit.default_p(ns)
+resolved = ModelingToolkit.varmap_to_vars(Dict(), parameters(ns), defaults=ModelingToolkit.default_p(ns))
+@test resolved == [1, 0.1+1, (0.1+1)*1.1]
+
+prob = NonlinearProblem(ns, ones(3), Pair[])
+@show sol = solve(prob,NewtonRaphson())
+
+@variables a
+@parameters b
+top = NonlinearSystem([0 ~ -a + ns.x+1], [a], [b], systems=[ns], name=:top)
+flat = flatten(top)
+flat.b = ns.σ*0.5
+
+flatres = ModelingToolkit.varmap_to_vars(Dict(), parameters(flat), defaults=ModelingToolkit.default_p(flat))
+@test flatres == [0.5, 1, 0.1+1, (0.1+1)*1.1]
+
+prob = NonlinearProblem(flat, ones(4), Pair[])
+@show sol = solve(prob,NewtonRaphson())