diff --git a/docs/src/API/variables.md b/docs/src/API/variables.md index 04d85e06b9..45b24100e1 100644 --- a/docs/src/API/variables.md +++ b/docs/src/API/variables.md @@ -293,7 +293,7 @@ For systems that contain parameters with metadata like described above, have som In the example below, we define a system with tunable parameters and extract bounds vectors ```@example metadata -@variables x(t)=0 u(t)=0 [input = true] y(t)=0 [output = true] +@variables x(t)=0 u(t)=0 [input=true] y(t)=0 [output=true] @parameters T [tunable = true, bounds = (0, Inf)] @parameters k [tunable = true, bounds = (0, Inf)] eqs = [D(x) ~ (-x + k * u) / T # A first-order system with time constant T and gain k diff --git a/docs/src/basics/Events.md b/docs/src/basics/Events.md index 89f874b08b..e38c4a63d6 100644 --- a/docs/src/basics/Events.md +++ b/docs/src/basics/Events.md @@ -92,8 +92,8 @@ The basic purely symbolic continuous event interface to encode *one* continuous event is ```julia -AbstractSystem(eqs, ...; continuous_events::Vector{Equation}) -AbstractSystem(eqs, ...; continuous_events::Pair{Vector{Equation}, Vector{Equation}}) +AbstractSystem(eqs, _...; continuous_events::Vector{Equation}) +AbstractSystem(eqs, _...; continuous_events::Pair{Vector{Equation}, Vector{Equation}}) ``` In the former, equations that evaluate to 0 will represent conditions that should @@ -272,7 +272,7 @@ In addition to continuous events, discrete events are also supported. The general interface to represent a collection of discrete events is ```julia -AbstractSystem(eqs, ...; discrete_events = [condition1 => affect1, condition2 => affect2]) +AbstractSystem(eqs, _...; discrete_events = [condition1 => affect1, condition2 => affect2]) ``` where conditions are symbolic expressions that should evaluate to `true` when an @@ -497,7 +497,8 @@ so far we aren't using anything that's not possible with the implicit interface. You can also write ```julia -[temp ~ furnace_off_threshold] => ModelingToolkit.ImperativeAffect(modified = (; +[temp ~ + furnace_off_threshold] => ModelingToolkit.ImperativeAffect(modified = (; furnace_on)) do x, o, i, c @set! x.furnace_on = false end diff --git a/docs/src/basics/FAQ.md b/docs/src/basics/FAQ.md index e3b12b46ab..f6fdd81210 100644 --- a/docs/src/basics/FAQ.md +++ b/docs/src/basics/FAQ.md @@ -63,7 +63,7 @@ The same principle applies to any parameter type that is not `Float64`. @parameters p1::Int # integer-valued @parameters p2::Bool # boolean-valued @parameters p3::MyCustomStructType # non-numeric -@parameters p4::ComponentArray{...} # non-standard array +@parameters p4::ComponentArray{_...} # non-standard array ``` ## Getting the index for a symbol diff --git a/docs/src/basics/MTKLanguage.md b/docs/src/basics/MTKLanguage.md index 09fa6d2daa..b15eec126f 100644 --- a/docs/src/basics/MTKLanguage.md +++ b/docs/src/basics/MTKLanguage.md @@ -381,7 +381,8 @@ Refer the following example for different ways to define symbolic arrays. @parameters begin p1[1:4] p2[1:N] - p3[1:N, 1:M] = 10, + p3[1:N, + 1:M] = 10, [description = "A multi-dimensional array of arbitrary length with description"] (p4[1:N, 1:M] = 10), [description = "An alternate syntax for p3 to match the syntax of vanilla parameters macro"] diff --git a/docs/src/basics/Validation.md b/docs/src/basics/Validation.md index 6e17beeded..3f36a06e5e 100644 --- a/docs/src/basics/Validation.md +++ b/docs/src/basics/Validation.md @@ -108,7 +108,7 @@ function ModelingToolkit.get_unit(op::typeof(dummycomplex), args) end sts = @variables a(t)=0 [unit = u"cm"] -ps = @parameters s=-1 [unit = u"cm"] c=c [unit = u"cm"] +ps = @parameters s=-1 [unit=u"cm"] c=c [unit=u"cm"] eqs = [D(a) ~ dummycomplex(c, s);] sys = System( eqs, t, [sts...;], [ps...;], name = :sys, checks = ~ModelingToolkit.CheckUnits) diff --git a/docs/src/examples/sparse_jacobians.md b/docs/src/examples/sparse_jacobians.md index 79a93c3471..a87f824d8d 100644 --- a/docs/src/examples/sparse_jacobians.md +++ b/docs/src/examples/sparse_jacobians.md @@ -23,15 +23,20 @@ function brusselator_2d_loop(du, u, p, t) @inbounds for I in CartesianIndices((N, N)) i, j = Tuple(I) x, y = xyd_brusselator[I[1]], xyd_brusselator[I[2]] - ip1, im1, jp1, jm1 = limit(i + 1, N), limit(i - 1, N), limit(j + 1, N), + ip1, im1, jp1, + jm1 = limit(i + 1, N), limit(i - 1, N), limit(j + 1, N), limit(j - 1, N) - du[i, j, 1] = alpha * (u[im1, j, 1] + u[ip1, j, 1] + u[i, jp1, 1] + u[i, jm1, 1] - - 4u[i, j, 1]) + - B + u[i, j, 1]^2 * u[i, j, 2] - (A + 1) * u[i, j, 1] + - brusselator_f(x, y, t) - du[i, j, 2] = alpha * (u[im1, j, 2] + u[ip1, j, 2] + u[i, jp1, 2] + u[i, jm1, 2] - - 4u[i, j, 2]) + - A * u[i, j, 1] - u[i, j, 1]^2 * u[i, j, 2] + du[i, + j, + 1] = alpha * (u[im1, j, 1] + u[ip1, j, 1] + u[i, jp1, 1] + u[i, jm1, 1] - + 4u[i, j, 1]) + + B + u[i, j, 1]^2 * u[i, j, 2] - (A + 1) * u[i, j, 1] + + brusselator_f(x, y, t) + du[i, + j, + 2] = alpha * (u[im1, j, 2] + u[ip1, j, 2] + u[i, jp1, 2] + u[i, jm1, 2] - + 4u[i, j, 2]) + + A * u[i, j, 1] - u[i, j, 1]^2 * u[i, j, 2] end end p = (3.4, 1.0, 10.0, step(xyd_brusselator)) diff --git a/docs/src/examples/tearing_parallelism.md b/docs/src/examples/tearing_parallelism.md index f123f8b7b3..924102eff0 100644 --- a/docs/src/examples/tearing_parallelism.md +++ b/docs/src/examples/tearing_parallelism.md @@ -15,7 +15,7 @@ using ModelingToolkit: t_nounits as t, D_nounits as D # Basic electric components @connector function Pin(; name) - @variables v(t)=1.0 i(t)=1.0 [connect = Flow] + @variables v(t)=1.0 i(t)=1.0 [connect=Flow] System(Equation[], t, [v, i], [], name = name) end @@ -36,7 +36,7 @@ function ConstantVoltage(; name, V = 1.0) end @connector function HeatPort(; name) - @variables T(t)=293.15 Q_flow(t)=0.0 [connect = Flow] + @variables T(t)=293.15 Q_flow(t)=0.0 [connect=Flow] System(Equation[], t, [T, Q_flow], [], name = name) end diff --git a/docs/src/tutorials/disturbance_modeling.md b/docs/src/tutorials/disturbance_modeling.md index 341077b76b..0178873e38 100644 --- a/docs/src/tutorials/disturbance_modeling.md +++ b/docs/src/tutorials/disturbance_modeling.md @@ -188,7 +188,9 @@ disturbance_inputs = [ssys.d1, ssys.d2] P = ssys.system_model outputs = [P.inertia1.phi, P.inertia2.phi, P.inertia1.w, P.inertia2.w] -(f_oop, f_ip), x_sym, p_sym, io_sys = ModelingToolkit.generate_control_function( +(f_oop, f_ip), x_sym, +p_sym, +io_sys = ModelingToolkit.generate_control_function( model_with_disturbance, inputs, disturbance_inputs; disturbance_argument = true) g = ModelingToolkit.build_explicit_observed_function( diff --git a/docs/src/tutorials/linear_analysis.md b/docs/src/tutorials/linear_analysis.md index 250dbaa6a7..057dd03910 100644 --- a/docs/src/tutorials/linear_analysis.md +++ b/docs/src/tutorials/linear_analysis.md @@ -39,7 +39,7 @@ This is signified by the name being the middle argument to `connect`. Of the above mentioned functions, all except for [`open_loop`](@ref) return the output of [`ModelingToolkit.linearize`](@ref), which is ```julia -matrices, simplified_sys = linearize(...) +matrices, simplified_sys = linearize(_...) # matrices = (; A, B, C, D) ``` diff --git a/ext/MTKCasADiDynamicOptExt.jl b/ext/MTKCasADiDynamicOptExt.jl index 73e4a77ed4..8cefb7b0e4 100644 --- a/ext/MTKCasADiDynamicOptExt.jl +++ b/ext/MTKCasADiDynamicOptExt.jl @@ -17,7 +17,9 @@ struct MXLinearInterpolation t::Vector{Float64} dt::Float64 end -Base.getindex(m::MXLinearInterpolation, i...) = length(i) == length(size(m.u)) ? m.u[i...] : m.u[i..., :] +function Base.getindex(m::MXLinearInterpolation, i...) + length(i) == length(size(m.u)) ? m.u[i...] : m.u[i..., :] +end mutable struct CasADiModel model::Opti @@ -55,7 +57,7 @@ function (M::MXLinearInterpolation)(τ) (i > length(M.t) || i < 1) && error("Cannot extrapolate past the tspan.") colons = ntuple(_ -> (:), length(size(M.u)) - 1) if i < length(M.t) - M.u[colons..., i] + Δ*(M.u[colons..., i+1] - M.u[colons..., i]) + M.u[colons..., i] + Δ * (M.u[colons..., i + 1] - M.u[colons..., i]) else M.u[colons..., i] end @@ -65,7 +67,9 @@ function MTK.CasADiDynamicOptProblem(sys::System, op, tspan; dt = nothing, steps = nothing, guesses = Dict(), kwargs...) - prob, _ = MTK.process_DynamicOptProblem(CasADiDynamicOptProblem, CasADiModel, sys, op, tspan; dt, steps, guesses, kwargs...) + prob, + _ = MTK.process_DynamicOptProblem( + CasADiDynamicOptProblem, CasADiModel, sys, op, tspan; dt, steps, guesses, kwargs...) prob end @@ -127,10 +131,10 @@ function MTK.lowered_integral(model::CasADiModel, expr, lo, hi) for (i, t) in enumerate(model.U.t) if lo < t < hi Δt = min(dt, t - lo) - total += (0.5*Δt*(expr[i] + expr[i-1])) + total += (0.5 * Δt * (expr[i] + expr[i - 1])) elseif t >= hi && (t - dt < hi) Δt = hi - t + dt - total += (0.5*Δt*(expr[i] + expr[i-1])) + total += (0.5 * Δt * (expr[i] + expr[i - 1])) end end model.tₛ * total @@ -186,9 +190,13 @@ struct CasADiCollocation <: AbstractCollocation tableau::DiffEqBase.ODERKTableau end -MTK.CasADiCollocation(solver, tableau = MTK.constructDefault()) = CasADiCollocation(solver, tableau) +function MTK.CasADiCollocation(solver, tableau = MTK.constructDefault()) + CasADiCollocation(solver, tableau) +end -function MTK.prepare_and_optimize!(prob::CasADiDynamicOptProblem, solver::CasADiCollocation; verbose = false, solver_options = Dict(), plugin_options = Dict(), kwargs...) +function MTK.prepare_and_optimize!( + prob::CasADiDynamicOptProblem, solver::CasADiCollocation; verbose = false, + solver_options = Dict(), plugin_options = Dict(), kwargs...) solver_opti = add_solve_constraints!(prob, solver.tableau) verbose || (solver_options["print_level"] = 0) solver!(solver_opti, "$(solver.solver)", plugin_options, solver_options) @@ -211,7 +219,7 @@ end function MTK.get_V_values(model::CasADiModel) value_getter = MTK.successful_solve(model) ? CasADi.debug_value : CasADi.value (nu, nt) = size(model.V.u) - if nu*nt != 0 + if nu * nt != 0 V_vals = value_getter(model.solver_opti, model.V.u) size(V_vals, 2) == 1 && (V_vals = V_vals') V_vals = [[V_vals[i, j] for i in 1:nu] for j in 1:nt] @@ -224,9 +232,11 @@ function MTK.get_t_values(model::CasADiModel) value_getter = MTK.successful_solve(model) ? CasADi.debug_value : CasADi.value ts = value_getter(model.solver_opti, model.tₛ) .* model.U.t end -MTK.objective_value(model::CasADiModel) = CasADi.pyconvert(Float64, model.solver_opti.py.value(model.solver_opti.py.f)) +function MTK.objective_value(model::CasADiModel) + CasADi.pyconvert(Float64, model.solver_opti.py.value(model.solver_opti.py.f)) +end -function MTK.successful_solve(m::CasADiModel) +function MTK.successful_solve(m::CasADiModel) isnothing(m.solver_opti) && return false retcode = CasADi.return_status(m.solver_opti) retcode == "Solve_Succeeded" || retcode == "Solved_To_Acceptable_Level" diff --git a/ext/MTKInfiniteOptExt.jl b/ext/MTKInfiniteOptExt.jl index 8150cd06f7..364a972725 100644 --- a/ext/MTKInfiniteOptExt.jl +++ b/ext/MTKInfiniteOptExt.jl @@ -48,26 +48,32 @@ struct InfiniteOptDynamicOptProblem{uType, tType, isinplace, P, F, K} <: end MTK.generate_internal_model(m::Type{InfiniteOptModel}) = InfiniteModel() -MTK.generate_time_variable!(m::InfiniteModel, tspan, tsteps) = @infinite_parameter(m, t in [tspan[1], tspan[2]], num_supports = length(tsteps)) -MTK.generate_state_variable!(m::InfiniteModel, u0::Vector, ns, ts) = @variable(m, U[i = 1:ns], Infinite(m[:t]), start=u0[i]) -MTK.generate_input_variable!(m::InfiniteModel, c0, nc, ts) = @variable(m, V[i = 1:nc], Infinite(m[:t]), start=c0[i]) +function MTK.generate_time_variable!(m::InfiniteModel, tspan, tsteps) + @infinite_parameter(m, t in [tspan[1], tspan[2]], num_supports=length(tsteps)) +end +function MTK.generate_state_variable!(m::InfiniteModel, u0::Vector, ns, ts) + @variable(m, U[i = 1:ns], Infinite(m[:t]), start=u0[i]) +end +function MTK.generate_input_variable!(m::InfiniteModel, c0, nc, ts) + @variable(m, V[i = 1:nc], Infinite(m[:t]), start=c0[i]) +end function MTK.generate_timescale!(m::InfiniteModel, guess, is_free_t) - @variable(m, tₛ ≥ 0, start = guess) + @variable(m, tₛ≥0, start=guess) if !is_free_t - fix(tₛ, 1, force=true) + fix(tₛ, 1, force = true) set_start_value(tₛ, 1) end tₛ end -function MTK.add_constraint!(m::InfiniteOptModel, expr::Union{Equation, Inequality}) +function MTK.add_constraint!(m::InfiniteOptModel, expr::Union{Equation, Inequality}) if expr isa Equation - @constraint(m.model, expr.lhs - expr.rhs == 0) + @constraint(m.model, expr.lhs - expr.rhs==0) elseif expr.relational_op === Symbolics.geq - @constraint(m.model, expr.lhs - expr.rhs ≥ 0) + @constraint(m.model, expr.lhs - expr.rhs≥0) else - @constraint(m.model, expr.lhs - expr.rhs ≤ 0) + @constraint(m.model, expr.lhs - expr.rhs≤0) end end MTK.set_objective!(m::InfiniteOptModel, expr) = @objective(m.model, Min, expr) @@ -76,7 +82,9 @@ function MTK.JuMPDynamicOptProblem(sys::System, op, tspan; dt = nothing, steps = nothing, guesses = Dict(), kwargs...) - prob, _ = MTK.process_DynamicOptProblem(JuMPDynamicOptProblem, InfiniteOptModel, sys, op, tspan; dt, steps, guesses, kwargs...) + prob, + _ = MTK.process_DynamicOptProblem(JuMPDynamicOptProblem, InfiniteOptModel, sys, + op, tspan; dt, steps, guesses, kwargs...) prob end @@ -84,12 +92,17 @@ function MTK.InfiniteOptDynamicOptProblem(sys::System, op, tspan; dt = nothing, steps = nothing, guesses = Dict(), kwargs...) - prob, pmap = MTK.process_DynamicOptProblem(InfiniteOptDynamicOptProblem, InfiniteOptModel, sys, op, tspan; dt, steps, guesses, kwargs...) + prob, + pmap = MTK.process_DynamicOptProblem( + InfiniteOptDynamicOptProblem, InfiniteOptModel, + sys, op, tspan; dt, steps, guesses, kwargs...) MTK.add_equational_constraints!(prob.wrapped_model, sys, pmap, tspan) prob end -MTK.lowered_integral(model::InfiniteOptModel, expr, lo, hi) = model.tₛ * InfiniteOpt.∫(expr, model.model[:t], lo, hi) +function MTK.lowered_integral(model::InfiniteOptModel, expr, lo, hi) + model.tₛ * InfiniteOpt.∫(expr, model.model[:t], lo, hi) +end MTK.lowered_derivative(model::InfiniteOptModel, i) = ∂(model.U[i], model.model[:t]) function MTK.process_integral_bounds(model::InfiniteOptModel, integral_span, tspan) @@ -125,7 +138,7 @@ function add_solve_constraints!(prob::JuMPDynamicOptProblem, tableau) nᵥ = length(V) if MTK.is_explicit(tableau) K = Any[] - for τ in tsteps[1:end-1] + for τ in tsteps[1:(end - 1)] for (i, h) in enumerate(c) ΔU = sum([A[i, j] * K[j] for j in 1:(i - 1)], init = zeros(nᵤ)) Uₙ = [U[i](τ) + ΔU[i] * dt for i in 1:nᵤ] @@ -142,14 +155,15 @@ function add_solve_constraints!(prob::JuMPDynamicOptProblem, tableau) K = @variable(model, K[1:length(α), 1:nᵤ], Infinite(model[:t])) ΔUs = A * K ΔU_tot = dt * (K' * α) - for τ in tsteps[1:end-1] + for τ in tsteps[1:(end - 1)] for (i, h) in enumerate(c) ΔU = @view ΔUs[i, :] Uₙ = U + ΔU * dt @constraint(model, [j = 1:nᵤ], K[i, j]==(tₛ * f(Uₙ, V, p, τ + h * dt)[j]), DomainRestrictions(t => τ), base_name="solve_K$i($τ)") end - @constraint(model, [n = 1:nᵤ], U[n](τ) + ΔU_tot[n]==U[n](min(τ + dt, tsteps[end])), + @constraint(model, + [n = 1:nᵤ], U[n](τ) + ΔU_tot[n]==U[n](min(τ + dt, tsteps[end])), DomainRestrictions(t => τ), base_name="solve_U($τ)") end end @@ -159,15 +173,21 @@ struct JuMPCollocation <: AbstractCollocation solver::Any tableau::DiffEqBase.ODERKTableau end -MTK.JuMPCollocation(solver, tableau = MTK.constructDefault()) = JuMPCollocation(solver, tableau) +function MTK.JuMPCollocation(solver, tableau = MTK.constructDefault()) + JuMPCollocation(solver, tableau) +end struct InfiniteOptCollocation <: AbstractCollocation solver::Any derivative_method::InfiniteOpt.AbstractDerivativeMethod end -MTK.InfiniteOptCollocation(solver, derivative_method = InfiniteOpt.FiniteDifference(InfiniteOpt.Backward())) = InfiniteOptCollocation(solver, derivative_method) +function MTK.InfiniteOptCollocation( + solver, derivative_method = InfiniteOpt.FiniteDifference(InfiniteOpt.Backward())) + InfiniteOptCollocation(solver, derivative_method) +end -function MTK.prepare_and_optimize!(prob::JuMPDynamicOptProblem, solver::JuMPCollocation; verbose = false, kwargs...) +function MTK.prepare_and_optimize!( + prob::JuMPDynamicOptProblem, solver::JuMPCollocation; verbose = false, kwargs...) model = prob.wrapped_model.model verbose || set_silent(model) # Unregister current solver constraints @@ -190,7 +210,8 @@ function MTK.prepare_and_optimize!(prob::JuMPDynamicOptProblem, solver::JuMPColl model end -function MTK.prepare_and_optimize!(prob::InfiniteOptDynamicOptProblem, solver::InfiniteOptCollocation; verbose = false, kwargs...) +function MTK.prepare_and_optimize!(prob::InfiniteOptDynamicOptProblem, + solver::InfiniteOptCollocation; verbose = false, kwargs...) model = prob.wrapped_model.model verbose || set_silent(model) set_derivative_method(model[:t], solver.derivative_method) @@ -223,8 +244,8 @@ function MTK.successful_solve(model::InfiniteModel) error("Model not solvable; please report this to github.com/SciML/ModelingToolkit.jl with a MWE.") pstatus === FEASIBLE_POINT && - (tstatus === OPTIMAL || tstatus === LOCALLY_SOLVED || tstatus === ALMOST_OPTIMAL || - tstatus === ALMOST_LOCALLY_SOLVED) + (tstatus === OPTIMAL || tstatus === LOCALLY_SOLVED || tstatus === ALMOST_OPTIMAL || + tstatus === ALMOST_LOCALLY_SOLVED) end import InfiniteOpt: JuMP, GeneralVariableRef diff --git a/ext/MTKPyomoDynamicOptExt.jl b/ext/MTKPyomoDynamicOptExt.jl index cb0aafc432..16209bfc79 100644 --- a/ext/MTKPyomoDynamicOptExt.jl +++ b/ext/MTKPyomoDynamicOptExt.jl @@ -108,7 +108,8 @@ function MTK.add_constraint!(pmodel::PyomoDynamicOptModel, cons; n_idxs = 1) else cons.lhs - cons.rhs ≤ 0 end - expr = Symbolics.substitute(Symbolics.unwrap(expr), SPECIAL_FUNCTIONS_DICT, fold = false) + expr = Symbolics.substitute( + Symbolics.unwrap(expr), SPECIAL_FUNCTIONS_DICT, fold = false) cons_sym = Symbol("cons", hash(cons)) if occursin(Symbolics.unwrap(t_sym), expr) @@ -141,17 +142,17 @@ end function MTK.lowered_integral(m::PyomoDynamicOptModel, arg, lo, hi) @unpack model, model_sym, t_sym, dummy_sym = m total = 0 - dt = Pyomo.pyconvert(Float64, (model.t.at(-1) - model.t.at(1))/(model.steps - 1)) + dt = Pyomo.pyconvert(Float64, (model.t.at(-1) - model.t.at(1)) / (model.steps - 1)) f = Symbolics.build_function(arg, model_sym, t_sym, expression = false) for (i, t) in enumerate(model.t) if Bool(lo < t) && Bool(t < hi) - t_p = model.t.at(i-1) + t_p = model.t.at(i - 1) Δt = min(t - lo, t - t_p) - total += 0.5*Δt*(f(model, t) + f(model, t_p)) + total += 0.5 * Δt * (f(model, t) + f(model, t_p)) elseif Bool(t >= hi) && Bool(t - dt < hi) - t_p = model.t.at(i-1) + t_p = model.t.at(i - 1) Δt = hi - t + dt - total += 0.5*Δt*(f(model, t) + f(model, t_p)) + total += 0.5 * Δt * (f(model, t) + f(model, t_p)) end end PyomoVar(model.tₛ * total) diff --git a/src/bipartite_graph.jl b/src/bipartite_graph.jl index b6665646c9..cf43241d0d 100644 --- a/src/bipartite_graph.jl +++ b/src/bipartite_graph.jl @@ -29,7 +29,7 @@ Base.iterate(u::Unassigned, state) = nothing Base.show(io::IO, ::Unassigned) = printstyled(io, "u"; color = :light_black) -struct Matching{U, V <: AbstractVector} <: AbstractVector{Union{U, Int}} #=> :Unassigned =# +struct Matching{U, V <: AbstractVector} <: AbstractVector{Union{U, Int}} match::V inv_match::Union{Nothing, V} end @@ -608,6 +608,7 @@ function Graphs.incidence_matrix(g::BipartiteGraph, val = true) I = Int[] J = Int[] for i in 𝑠vertices(g), n in 𝑠neighbors(g, i) + push!(I, i) push!(J, n) end diff --git a/src/deprecations.jl b/src/deprecations.jl index aaafa0d283..bda4ff1dac 100644 --- a/src/deprecations.jl +++ b/src/deprecations.jl @@ -62,6 +62,7 @@ for T in [:ODEProblem, :DDEProblem, :SDEProblem, :SDDEProblem, :DAEProblem, end for pType in [SciMLBase.NullParameters, Nothing], uType in [Any, Nothing] + @eval function SciMLBase.$T(sys::System, u0::$uType, tspan, p::$pType; kw...) ctor = string($T) pT = string($(QuoteNode(pType))) @@ -142,6 +143,7 @@ for T in [:NonlinearProblem, :NonlinearLeastSquaresProblem, end end for pType in [SciMLBase.NullParameters, Nothing], uType in [Any, Nothing] + @eval function SciMLBase.$T(sys::System, u0::$uType, p::$pType; kw...) ctor = string($T) pT = string($(QuoteNode(pType))) diff --git a/src/inputoutput.jl b/src/inputoutput.jl index ac022cbe2f..6a17f76d80 100644 --- a/src/inputoutput.jl +++ b/src/inputoutput.jl @@ -120,7 +120,7 @@ function same_or_inner_namespace(u, var) nu == nv || # namespaces are the same startswith(nv, nu) || # or nv starts with nu, i.e., nv is an inner namespace to nu occursin(NAMESPACE_SEPARATOR, string(getname(var))) && - !occursin(NAMESPACE_SEPARATOR, string(getname(u))) # or u is top level but var is internal + !occursin(NAMESPACE_SEPARATOR, string(getname(u))) # or u is top level but var is internal end function inner_namespace(u, var) @@ -129,7 +129,7 @@ function inner_namespace(u, var) nu == nv && return false startswith(nv, nu) || # or nv starts with nu, i.e., nv is an inner namespace to nu occursin(NAMESPACE_SEPARATOR, string(getname(var))) && - !occursin(NAMESPACE_SEPARATOR, string(getname(u))) # or u is top level but var is internal + !occursin(NAMESPACE_SEPARATOR, string(getname(u))) # or u is top level but var is internal end """ diff --git a/src/linearization.jl b/src/linearization.jl index 88fe1d39d3..72ace55e89 100644 --- a/src/linearization.jl +++ b/src/linearization.jl @@ -297,7 +297,8 @@ function (linfun::LinearizationFunction)(u, p, t) error("Number of unknown variables ($(linfun.num_states)) does not match the number of input unknowns ($(length(u)))") integ_cache = (linfun.caches,) integ = MockIntegrator{true}(u, p, t, fun, integ_cache, nothing) - u, p, success = SciMLBase.get_initial_values( + u, p, + success = SciMLBase.get_initial_values( linfun.prob, integ, fun, linfun.initializealg, Val(true); linfun.initialize_kwargs...) if !success @@ -744,7 +745,8 @@ function linearize(sys, inputs, outputs; op = Dict(), t = 0.0, allow_input_derivatives = false, zero_dummy_der = false, kwargs...) - lin_fun, ssys = linearization_function(sys, + lin_fun, + ssys = linearization_function(sys, inputs, outputs; zero_dummy_der, diff --git a/src/modelingtoolkitize/sdeproblem.jl b/src/modelingtoolkitize/sdeproblem.jl index ff8c238559..0f63a35fc1 100644 --- a/src/modelingtoolkitize/sdeproblem.jl +++ b/src/modelingtoolkitize/sdeproblem.jl @@ -26,7 +26,8 @@ function modelingtoolkitize( odefn = ODEFunction{SciMLBase.isinplace(prob)}( prob.f.f; mass_matrix = prob.f.mass_matrix, sys = prob.f.sys) odeprob = ODEProblem(odefn, prob.u0, prob.tspan, prob.p) - sys, vars, params = modelingtoolkitize( + sys, vars, + params = modelingtoolkitize( odeprob; u_names, p_names, return_symbolic_u0_p = true, name = gensym(:MTKizedSDE), kwargs...) t = get_iv(sys) diff --git a/src/problems/bvproblem.jl b/src/problems/bvproblem.jl index 6344dfaa0b..1c193bca51 100644 --- a/src/problems/bvproblem.jl +++ b/src/problems/bvproblem.jl @@ -14,7 +14,8 @@ # for initialization. _op = has_alg_eqs(sys) ? op : merge(Dict(op), Dict(guesses)) - fode, u0, p = process_SciMLProblem( + fode, u0, + p = process_SciMLProblem( ODEFunction{iip, spec}, sys, _op; guesses, t = tspan !== nothing ? tspan[1] : tspan, check_compatibility = false, cse, checkbounds, time_dependent_init = false, expression, kwargs...) diff --git a/src/problems/daeproblem.jl b/src/problems/daeproblem.jl index 6134923d4d..7ac5b0b721 100644 --- a/src/problems/daeproblem.jl +++ b/src/problems/daeproblem.jl @@ -67,7 +67,9 @@ end check_complete(sys, DAEProblem) check_compatibility && check_compatible_system(DAEProblem, sys) - f, du0, u0, p = process_SciMLProblem(DAEFunction{iip, spec}, sys, op; + f, du0, + u0, + p = process_SciMLProblem(DAEFunction{iip, spec}, sys, op; t = tspan !== nothing ? tspan[1] : tspan, check_length, eval_expression, eval_module, check_compatibility, implicit_dae = true, expression, kwargs...) diff --git a/src/problems/ddeproblem.jl b/src/problems/ddeproblem.jl index 45af3edddb..c86b8a6c2a 100644 --- a/src/problems/ddeproblem.jl +++ b/src/problems/ddeproblem.jl @@ -46,7 +46,8 @@ end check_complete(sys, DDEProblem) check_compatibility && check_compatible_system(DDEProblem, sys) - f, u0, p = process_SciMLProblem(DDEFunction{iip, spec}, sys, op; + f, u0, + p = process_SciMLProblem(DDEFunction{iip, spec}, sys, op; t = tspan !== nothing ? tspan[1] : tspan, check_length, cse, checkbounds, eval_expression, eval_module, check_compatibility, symbolic_u0 = true, expression, u0_constructor, kwargs...) diff --git a/src/problems/discreteproblem.jl b/src/problems/discreteproblem.jl index c4ab069ebe..820819d282 100644 --- a/src/problems/discreteproblem.jl +++ b/src/problems/discreteproblem.jl @@ -45,7 +45,8 @@ end dvs = unknowns(sys) op = to_varmap(op, dvs) add_toterms!(op; replace = true) - f, u0, p = process_SciMLProblem(DiscreteFunction{iip, spec}, sys, op; + f, u0, + p = process_SciMLProblem(DiscreteFunction{iip, spec}, sys, op; t = tspan !== nothing ? tspan[1] : tspan, check_compatibility, expression, kwargs...) diff --git a/src/problems/implicitdiscreteproblem.jl b/src/problems/implicitdiscreteproblem.jl index 9e500d7354..9c85153c4f 100644 --- a/src/problems/implicitdiscreteproblem.jl +++ b/src/problems/implicitdiscreteproblem.jl @@ -50,7 +50,8 @@ end dvs = unknowns(sys) op = to_varmap(op, dvs) add_toterms!(op; replace = true) - f, u0, p = process_SciMLProblem( + f, u0, + p = process_SciMLProblem( ImplicitDiscreteFunction{iip, spec}, sys, op; t = tspan !== nothing ? tspan[1] : tspan, check_compatibility, expression, kwargs...) diff --git a/src/problems/intervalnonlinearproblem.jl b/src/problems/intervalnonlinearproblem.jl index 2dc929fe25..d7221632d3 100644 --- a/src/problems/intervalnonlinearproblem.jl +++ b/src/problems/intervalnonlinearproblem.jl @@ -33,7 +33,8 @@ function SciMLBase.IntervalNonlinearProblem( u0map = unknowns(sys) .=> uspan[1] op = anydict([unknowns(sys)[1] => uspan[1]]) merge!(op, to_varmap(parammap, parameters(sys))) - f, u0, p = process_SciMLProblem(IntervalNonlinearFunction, sys, op; + f, u0, + p = process_SciMLProblem(IntervalNonlinearFunction, sys, op; check_compatibility, expression, kwargs...) kwargs = process_kwargs(sys; kwargs...) diff --git a/src/problems/jumpproblem.jl b/src/problems/jumpproblem.jl index 28de46c35f..96df4d4b96 100644 --- a/src/problems/jumpproblem.jl +++ b/src/problems/jumpproblem.jl @@ -22,7 +22,8 @@ build_initializeprob = false, checkbounds, cse, check_length = false, kwargs...) else - _, u0, p = process_SciMLProblem(EmptySciMLFunction{iip}, sys, op; + _, u0, + p = process_SciMLProblem(EmptySciMLFunction{iip}, sys, op; t = tspan === nothing ? nothing : tspan[1], check_length = false, build_initializeprob = false, kwargs...) observedfun = ObservedFunctionCache(sys; eval_expression, eval_module, @@ -32,7 +33,8 @@ prob = ODEProblem{true}(df, u0, tspan, p; kwargs...) end else - _f, u0, p = process_SciMLProblem(EmptySciMLFunction{iip}, sys, op; + _f, u0, + p = process_SciMLProblem(EmptySciMLFunction{iip}, sys, op; t = tspan === nothing ? nothing : tspan[1], check_length = false, build_initializeprob = false, cse, kwargs...) f = DiffEqBase.DISCRETE_INPLACE_DEFAULT diff --git a/src/problems/nonlinearproblem.jl b/src/problems/nonlinearproblem.jl index 1ed7cda750..093e8db762 100644 --- a/src/problems/nonlinearproblem.jl +++ b/src/problems/nonlinearproblem.jl @@ -63,7 +63,8 @@ end end check_compatibility && check_compatible_system(NonlinearProblem, sys) - f, u0, p = process_SciMLProblem(NonlinearFunction{iip, spec}, sys, op; + f, u0, + p = process_SciMLProblem(NonlinearFunction{iip, spec}, sys, op; check_length, check_compatibility, expression, kwargs...) kwargs = process_kwargs(sys; kwargs...) @@ -79,7 +80,8 @@ end check_complete(sys, NonlinearLeastSquaresProblem) check_compatibility && check_compatible_system(NonlinearLeastSquaresProblem, sys) - f, u0, p = process_SciMLProblem(NonlinearFunction{iip}, sys, op; + f, u0, + p = process_SciMLProblem(NonlinearFunction{iip}, sys, op; check_length, expression, kwargs...) kwargs = process_kwargs(sys; kwargs...) diff --git a/src/problems/odeproblem.jl b/src/problems/odeproblem.jl index 6726322907..0cdfcf003c 100644 --- a/src/problems/odeproblem.jl +++ b/src/problems/odeproblem.jl @@ -70,7 +70,8 @@ end check_complete(sys, ODEProblem) check_compatibility && check_compatible_system(ODEProblem, sys) - f, u0, p = process_SciMLProblem(ODEFunction{iip, spec}, sys, op; + f, u0, + p = process_SciMLProblem(ODEFunction{iip, spec}, sys, op; t = tspan !== nothing ? tspan[1] : tspan, check_length, eval_expression, eval_module, expression, check_compatibility, kwargs...) @@ -88,7 +89,8 @@ end check_complete(sys, SteadyStateProblem) check_compatibility && check_compatible_system(SteadyStateProblem, sys) - f, u0, p = process_SciMLProblem(ODEFunction{iip}, sys, op; + f, u0, + p = process_SciMLProblem(ODEFunction{iip}, sys, op; steady_state = true, check_length, check_compatibility, expression, time_dependent_init = false, kwargs...) diff --git a/src/problems/optimizationproblem.jl b/src/problems/optimizationproblem.jl index 36b34a8867..978f973d83 100644 --- a/src/problems/optimizationproblem.jl +++ b/src/problems/optimizationproblem.jl @@ -23,7 +23,8 @@ function SciMLBase.OptimizationFunction{iip}(sys::System; _grad = nothing end if hess - _hess, hess_prototype = generate_cost_hessian( + _hess, + hess_prototype = generate_cost_hessian( sys; expression, wrap_gfw = Val{true}, eval_expression, eval_module, checkbounds, cse, sparse, simplify, return_sparsity = true, kwargs...) @@ -40,7 +41,8 @@ function SciMLBase.OptimizationFunction{iip}(sys::System; cons = generate_cons(sys; expression, wrap_gfw = Val{true}, eval_expression, eval_module, checkbounds, cse, kwargs...) if cons_j - _cons_j, cons_jac_prototype = generate_constraint_jacobian( + _cons_j, + cons_jac_prototype = generate_constraint_jacobian( sys; expression, wrap_gfw = Val{true}, eval_expression, eval_module, checkbounds, cse, simplify, sparse = cons_sparse, return_sparsity = true, kwargs...) @@ -48,7 +50,8 @@ function SciMLBase.OptimizationFunction{iip}(sys::System; _cons_j = cons_jac_prototype = nothing end if cons_h - _cons_h, cons_hess_prototype = generate_constraint_hessian( + _cons_h, + cons_hess_prototype = generate_constraint_hessian( sys; expression, wrap_gfw = Val{true}, eval_expression, eval_module, checkbounds, cse, simplify, sparse = cons_sparse, return_sparsity = true, kwargs...) @@ -92,7 +95,8 @@ function SciMLBase.OptimizationProblem{iip}( check_complete(sys, OptimizationProblem) check_compatibility && check_compatible_system(OptimizationProblem, sys) - f, u0, p = process_SciMLProblem(OptimizationFunction{iip}, sys, op; + f, u0, + p = process_SciMLProblem(OptimizationFunction{iip}, sys, op; check_compatibility, tofloat = false, check_length = false, expression, kwargs...) dvs = unknowns(sys) diff --git a/src/problems/sccnonlinearproblem.jl b/src/problems/sccnonlinearproblem.jl index d2afb7b315..2a44e3de4e 100644 --- a/src/problems/sccnonlinearproblem.jl +++ b/src/problems/sccnonlinearproblem.jl @@ -111,7 +111,8 @@ function SciMLBase.SCCNonlinearProblem{iip}(sys::System, op; eval_expression = f eqs = equations(sys) obs = observed(sys) - _, u0, p = process_SciMLProblem( + _, u0, + p = process_SciMLProblem( EmptySciMLFunction{iip}, sys, op; eval_expression, eval_module, kwargs...) explicitfuns = [] diff --git a/src/problems/sddeproblem.jl b/src/problems/sddeproblem.jl index e1cc00b2a7..3f852577d0 100644 --- a/src/problems/sddeproblem.jl +++ b/src/problems/sddeproblem.jl @@ -48,7 +48,8 @@ end check_complete(sys, SDDEProblem) check_compatibility && check_compatible_system(SDDEProblem, sys) - f, u0, p = process_SciMLProblem(SDDEFunction{iip, spec}, sys, op; + f, u0, + p = process_SciMLProblem(SDDEFunction{iip, spec}, sys, op; t = tspan !== nothing ? tspan[1] : tspan, check_length, cse, checkbounds, eval_expression, eval_module, check_compatibility, sparse, symbolic_u0 = true, expression, u0_constructor, kwargs...) diff --git a/src/problems/sdeproblem.jl b/src/problems/sdeproblem.jl index 1bc47118ff..86f8dd8e84 100644 --- a/src/problems/sdeproblem.jl +++ b/src/problems/sdeproblem.jl @@ -72,7 +72,8 @@ end check_complete(sys, SDEProblem) check_compatibility && check_compatible_system(SDEProblem, sys) - f, u0, p = process_SciMLProblem(SDEFunction{iip, spec}, sys, op; + f, u0, + p = process_SciMLProblem(SDEFunction{iip, spec}, sys, op; t = tspan !== nothing ? tspan[1] : tspan, check_length, eval_expression, eval_module, check_compatibility, sparse, expression, kwargs...) diff --git a/src/structural_transformation/bareiss.jl b/src/structural_transformation/bareiss.jl index 7957427e5d..602f656c27 100644 --- a/src/structural_transformation/bareiss.jl +++ b/src/structural_transformation/bareiss.jl @@ -131,6 +131,7 @@ end function bareiss_update!(zero!, M::StridedMatrix, k, swapto, pivot, prev_pivot) @inbounds for i in (k + 1):size(M, 2), j in (k + 1):size(M, 1) + M[j, i] = exactdiv(M[j, i] * pivot - M[j, k] * M[k, i], prev_pivot) end zero!(M, (k + 1):size(M, 1), k) @@ -269,6 +270,7 @@ function reduce_echelon!(A::AbstractMatrix{T}, rank, d, end @label out @inbounds for i in (rank + 1):m, j in 1:n + A[i, j] = zero(T) end isreduced && return A diff --git a/src/structural_transformation/codegen.jl b/src/structural_transformation/codegen.jl index 9afe7ec5e7..37a5b380ac 100644 --- a/src/structural_transformation/codegen.jl +++ b/src/structural_transformation/codegen.jl @@ -45,6 +45,7 @@ function torn_system_with_nlsolve_jacobian_sparsity(state, var_eq_matching, var_ var_rename = ones(Int64, ndsts(graph)) nlsolve_vars = Int[] for i in nlsolve_scc_idxs, c in var_sccs[i] + append!(nlsolve_vars, c) for v in c var_rename[v] = 0 diff --git a/src/structural_transformation/partial_state_selection.jl b/src/structural_transformation/partial_state_selection.jl index 89447a4e16..ab8e7f0f3d 100644 --- a/src/structural_transformation/partial_state_selection.jl +++ b/src/structural_transformation/partial_state_selection.jl @@ -213,7 +213,8 @@ function tearing_with_dummy_derivatives(structure, dummy_derivatives) can_eliminate[v] = true end end - var_eq_matching, full_var_eq_matching, var_sccs = tear_graph_modia(structure, + var_eq_matching, full_var_eq_matching, + var_sccs = tear_graph_modia(structure, Base.Fix1(isdiffed, (structure, dummy_derivatives)), Union{Unassigned, SelectedState}; varfilter = Base.Fix1(getindex, can_eliminate)) diff --git a/src/structural_transformation/symbolics_tearing.jl b/src/structural_transformation/symbolics_tearing.jl index af0973f1dc..d9a49a4e15 100644 --- a/src/structural_transformation/symbolics_tearing.jl +++ b/src/structural_transformation/symbolics_tearing.jl @@ -338,6 +338,7 @@ function generate_derivative_variables!( # We need the inverse mapping of `var_sccs` to update it efficiently later. v_to_scc = Vector{NTuple{2, Int}}(undef, ndsts(graph)) for (i, scc) in enumerate(var_sccs), (j, v) in enumerate(scc) + v_to_scc[v] = (i, j) end # Pairs of `(x_t, dx)` added below @@ -475,7 +476,7 @@ function find_duplicate_dd(dv, solvable_graph, diff_to_var, linear_eqs, mm) if length(nzs) == 2 && (abs(nzs[1]) == 1 && nzs[1] == -nzs[2]) && (v_t = rvs[1] == dv ? rvs[2] : rvs[1]; - diff_to_var[v_t] === nothing) + diff_to_var[v_t] === nothing) @assert dv in rvs return eq, v_t end @@ -1064,7 +1065,11 @@ function tearing_reassemble(state::TearingState, var_eq_matching::Matching, var_sccs = generate_derivative_variables!( state, neweqs, var_eq_matching, full_var_eq_matching, var_sccs; mm, iv, D) - neweqs, solved_eqs, eq_ordering, var_ordering, nelim_eq, nelim_var = generate_system_equations!( + neweqs, solved_eqs, + eq_ordering, + var_ordering, + nelim_eq, + nelim_var = generate_system_equations!( state, neweqs, var_eq_matching, full_var_eq_matching, var_sccs, extra_eqs_vars; simplify, iv, D) @@ -1123,6 +1128,7 @@ function add_additional_history!( # We need the inverse mapping of `var_sccs` to update it efficiently later. v_to_scc = Vector{NTuple{2, Int}}(undef, ndsts(graph)) for (i, scc) in enumerate(var_sccs), (j, v) in enumerate(scc) + v_to_scc[v] = (i, j) end @@ -1346,7 +1352,8 @@ function dummy_derivative(sys, state = TearingState(sys); simplify = false, p end end - var_eq_matching, full_var_eq_matching, var_sccs, can_eliminate, summary = dummy_derivative_graph!( + var_eq_matching, full_var_eq_matching, var_sccs, + can_eliminate, summary = dummy_derivative_graph!( state, jac; state_priority, kwargs...) tearing_reassemble(state, var_eq_matching, full_var_eq_matching, var_sccs; diff --git a/src/structural_transformation/tearing.jl b/src/structural_transformation/tearing.jl index 31ebb8370a..67933ffe0e 100644 --- a/src/structural_transformation/tearing.jl +++ b/src/structural_transformation/tearing.jl @@ -74,6 +74,7 @@ function free_equations(graph, vars_scc, var_eq_matching, varfilter::F) where {F ne = nsrcs(graph) seen_eqs = falses(ne) for vars in vars_scc, var in vars + varfilter(var) || continue ieq = var_eq_matching[var] if ieq isa Int diff --git a/src/structural_transformation/utils.jl b/src/structural_transformation/utils.jl index 02a9763226..9cac634227 100644 --- a/src/structural_transformation/utils.jl +++ b/src/structural_transformation/utils.jl @@ -170,6 +170,7 @@ function sorted_incidence_matrix(ts::TransformationState, val = true; only_algeq varidx = 0 eqidx = 0 for vs in var_scc, v in vs + eq = var_eq_matching[v] if eq !== unassigned eqsmap[eq] = (eqidx += 1) @@ -255,9 +256,9 @@ function find_eq_solvables!(state::TearingState, ieq, to_rm = Int[], coeffs = no # if any of the variables in `a` are present in fullvars (taking into account arrays) if any( v -> any(isequal(v), fullvars) || - symbolic_type(v) == ArraySymbolic() && - Symbolics.shape(v) != Symbolics.Unknown() && - any(x -> any(isequal(x), fullvars), collect(v)), + symbolic_type(v) == ArraySymbolic() && + Symbolics.shape(v) != Symbolics.Unknown() && + any(x -> any(isequal(x), fullvars), collect(v)), vars(a)) continue end diff --git a/src/systems/abstractsystem.jl b/src/systems/abstractsystem.jl index 941739f127..18485e4a9f 100644 --- a/src/systems/abstractsystem.jl +++ b/src/systems/abstractsystem.jl @@ -731,46 +731,46 @@ function unflatten_parameters!(buffer, params, all_ps) end for prop in [:eqs - :tag - :noise_eqs - :iv - :unknowns - :ps - :tspan - :brownians - :jumps - :name - :description - :var_to_name - :defaults - :guesses - :observed - :systems - :constraints - :bcs - :domain - :ivs - :dvs - :connector_type - :preface - :initializesystem - :initialization_eqs - :schedule - :tearing_state - :metadata - :gui_metadata - :is_initializesystem - :is_discrete - :parameter_dependencies - :assertions - :ignored_connections - :parent - :is_dde - :tstops - :index_cache - :isscheduled - :costs - :consolidate] + :tag + :noise_eqs + :iv + :unknowns + :ps + :tspan + :brownians + :jumps + :name + :description + :var_to_name + :defaults + :guesses + :observed + :systems + :constraints + :bcs + :domain + :ivs + :dvs + :connector_type + :preface + :initializesystem + :initialization_eqs + :schedule + :tearing_state + :metadata + :gui_metadata + :is_initializesystem + :is_discrete + :parameter_dependencies + :assertions + :ignored_connections + :parent + :is_dde + :tstops + :index_cache + :isscheduled + :costs + :consolidate] fname_get = Symbol(:get_, prop) fname_has = Symbol(:has_, prop) @eval begin @@ -1626,7 +1626,8 @@ function full_equations(sys::AbstractSystem; simplify = false) subs = get_substitutions(sys) neweqs = map(equations(sys)) do eq if iscall(eq.lhs) && operation(eq.lhs) isa Union{Shift, Differential} - return substitute_and_simplify(eq.lhs, subs, simplify) ~ substitute_and_simplify( + return substitute_and_simplify(eq.lhs, subs, simplify) ~ + substitute_and_simplify( eq.rhs, subs, simplify) else @@ -2364,7 +2365,8 @@ function component_post_processing(expr, isconnector) if $isconnector $Setfield.@set!(res.connector_type=$connector_type(res)) end - $Setfield.@set!(res.gui_metadata=$GUIMetadata($GlobalRef(@__MODULE__, name))) + $Setfield.@set!(res.gui_metadata=$GUIMetadata($GlobalRef( + @__MODULE__, name))) else res end @@ -2812,8 +2814,8 @@ function process_parameter_equations(sys::AbstractSystem) if all(varsbuf) do sym is_parameter(sys, sym) || symbolic_type(sym) == ArraySymbolic() && - is_sized_array_symbolic(sym) && - all(Base.Fix1(is_parameter, sys), collect(sym)) + is_sized_array_symbolic(sym) && + all(Base.Fix1(is_parameter, sys), collect(sym)) end if !isparameter(eq.lhs) throw(ArgumentError(""" diff --git a/src/systems/alias_elimination.jl b/src/systems/alias_elimination.jl index fb4fedc920..4d769cf440 100644 --- a/src/systems/alias_elimination.jl +++ b/src/systems/alias_elimination.jl @@ -215,6 +215,7 @@ function find_linear_variables(graph, linear_equations, var_to_diff, irreducible eqs = get(var_to_lineq, v, nothing) eqs === nothing && continue for eq in eqs, v′ in 𝑠neighbors(graph, eq) + if linear_variables[v′] linear_variables[v′] = false push!(stack, v′) @@ -224,6 +225,7 @@ function find_linear_variables(graph, linear_equations, var_to_diff, irreducible end end for eq in linear_equations, v in 𝑠neighbors(graph, eq) + linear_variables[v] = true vlineqs = get!(() -> BitSet(), var_to_lineq, v) push!(vlineqs, eq) diff --git a/src/systems/analysis_points.jl b/src/systems/analysis_points.jl index 3b65dd6669..e5be7c6c6b 100644 --- a/src/systems/analysis_points.jl +++ b/src/systems/analysis_points.jl @@ -238,8 +238,10 @@ end Return all the namespaces in `name`. Namespaces should be separated by `.` or `$NAMESPACE_SEPARATOR`. """ -namespace_hierarchy(name::Symbol) = map( - Symbol, split(string(name), ('.', NAMESPACE_SEPARATOR))) +function namespace_hierarchy(name::Symbol) + map( + Symbol, split(string(name), ('.', NAMESPACE_SEPARATOR))) +end """ $(TYPEDSIGNATURES) @@ -311,15 +313,19 @@ by `fn`. `root`, and returns a 2-tuple consisting of the modified version of `target` and a tuple of the extra variables added. """ -modify_nested_subsystem(fn, root::AbstractSystem, target::AbstractSystem) = modify_nested_subsystem( - fn, root, nameof(target)) +function modify_nested_subsystem(fn, root::AbstractSystem, target::AbstractSystem) + modify_nested_subsystem( + fn, root, nameof(target)) +end """ $(TYPEDSIGNATURES) Apply the modification to the system containing the namespaced analysis point `target`. """ -modify_nested_subsystem(fn, root::AbstractSystem, target::AnalysisPoint) = modify_nested_subsystem( - fn, root, @view namespace_hierarchy(nameof(target))[1:(end - 1)]) +function modify_nested_subsystem(fn, root::AbstractSystem, target::AnalysisPoint) + modify_nested_subsystem( + fn, root, @view namespace_hierarchy(nameof(target))[1:(end - 1)]) +end """ $(TYPEDSIGNATURES) @@ -327,8 +333,10 @@ Apply the modification to the nested subsystem of `root` whose namespaced name m the provided name `target`. The namespace separator in `target` should be `.` or `$NAMESPACE_SEPARATOR`. The `target` may include `nameof(root)` as the first namespace. """ -modify_nested_subsystem(fn, root::AbstractSystem, target::Symbol) = modify_nested_subsystem( - fn, root, namespace_hierarchy(target)) +function modify_nested_subsystem(fn, root::AbstractSystem, target::Symbol) + modify_nested_subsystem( + fn, root, namespace_hierarchy(target)) +end """ $(TYPEDSIGNATURES) @@ -388,8 +396,10 @@ end Given a system `sys` and analysis point `ap`, return the index in `get_eqs(sys)` containing an equation which has as it's RHS an analysis point with name `nameof(ap)`. """ -analysis_point_index(sys::AbstractSystem, ap::AnalysisPoint) = analysis_point_index( - sys, nameof(ap)) +function analysis_point_index(sys::AbstractSystem, ap::AnalysisPoint) + analysis_point_index( + sys, nameof(ap)) +end """ $(TYPEDSIGNATURES) @@ -629,7 +639,8 @@ function apply_transformation(tf::PerturbOutput, sys::AbstractSystem) tf.with_output || return ap_sys, (new_var,) # add output variable, equation, default - out_var, out_def = get_analysis_variable( + out_var, + out_def = get_analysis_variable( ap_ivar, nameof(ap), get_iv(sys); perturb = false) push!(ap_sys_eqs, out_var ~ ap_ivar + wrap(new_var)) push!(unks, out_var) @@ -679,7 +690,8 @@ function apply_transformation(tf::AddVariable, sys::AbstractSystem) # add equations involving new variable ap_ivar = ap_var(ap.input) - new_var, new_def = get_analysis_variable( + new_var, + new_def = get_analysis_variable( ap_ivar, tf.name, get_iv(sys); perturb = false) # add variable unks = copy(get_unknowns(ap_sys)) @@ -728,7 +740,8 @@ end function apply_transformation(cst::ComplementarySensitivityTransform, sys::AbstractSystem) sys, (u,) = apply_transformation(GetInput(cst.ap), sys) - sys, (du,) = apply_transformation( + sys, + (du,) = apply_transformation( AddVariable( cst.ap, Symbol(namespace_hierarchy(nameof(cst.ap))[end], :_comp_sens_du)), sys) @@ -906,7 +919,8 @@ for f in [:get_sensitivity, :get_comp_sensitivity, :get_looptransfer] utility_fun = Symbol(f, :_function) @eval function $f( sys, ap, args...; loop_openings = [], system_modifier = identity, kwargs...) - lin_fun, ssys = $(utility_fun)( + lin_fun, + ssys = $(utility_fun)( sys, ap, args...; loop_openings, system_modifier, kwargs...) ModelingToolkit.linearize(ssys, lin_fun), ssys end @@ -974,7 +988,8 @@ end function linearization_function(sys::AbstractSystem, inputs::Union{Symbol, Vector{Symbol}, AnalysisPoint, Vector{AnalysisPoint}}, outputs; loop_openings = [], system_modifier = identity, kwargs...) - sys, input_vars, output_vars = linearization_ap_transform( + sys, input_vars, + output_vars = linearization_ap_transform( sys, inputs, outputs, loop_openings) return linearization_function(system_modifier(sys), input_vars, output_vars; kwargs...) end diff --git a/src/systems/callbacks.jl b/src/systems/callbacks.jl index a12dad9d72..e443f4a0c7 100644 --- a/src/systems/callbacks.jl +++ b/src/systems/callbacks.jl @@ -526,7 +526,7 @@ function Base.:(==)(e1::AbstractCallback, e2::AbstractCallback) (is_discrete(e1) === is_discrete(e2)) || return false (isequal(e1.conditions, e2.conditions) && isequal(e1.affect, e2.affect) && isequal(e1.initialize, e2.initialize) && isequal(e1.finalize, e2.finalize)) && - isequal(e1.reinitializealg, e2.reinitializealg) || + isequal(e1.reinitializealg, e2.reinitializealg) || return false is_discrete(e1) || (isequal(e1.affect_neg, e2.affect_neg) && isequal(e1.rootfind, e2.rootfind)) @@ -817,8 +817,8 @@ function compile_equational_affect( obseqs, Dict([p => unPre(p) for p in parameters(affsys)])) rhss = map(x -> x.rhs, update_eqs) lhss = map(x -> x.lhs, update_eqs) - is_p = [lhs ∈ Set(ps_to_update) for lhs in lhss] - is_u = [lhs ∈ Set(dvs_to_update) for lhs in lhss] + is_p = [lhs in Set(ps_to_update) for lhs in lhss] + is_u = [lhs in Set(dvs_to_update) for lhs in lhss] dvs = unknowns(sys) ps = parameters(sys) t = get_iv(sys) @@ -835,11 +835,13 @@ function compile_equational_affect( _ps = reorder_parameters(sys, ps) integ = gensym(:MTKIntegrator) - u_up, u_up! = build_function_wrapper(sys, (@view rhss[is_u]), dvs, _ps..., t; + u_up, + u_up! = build_function_wrapper(sys, (@view rhss[is_u]), dvs, _ps..., t; wrap_code = add_integrator_header(sys, integ, :u), expression = Val{false}, outputidxs = u_idxs, wrap_mtkparameters, cse = false, eval_expression, eval_module) - p_up, p_up! = build_function_wrapper(sys, (@view rhss[is_p]), dvs, _ps..., t; + p_up, + p_up! = build_function_wrapper(sys, (@view rhss[is_p]), dvs, _ps..., t; wrap_code = add_integrator_header(sys, integ, :p), expression = Val{false}, outputidxs = p_idxs, wrap_mtkparameters, cse = false, eval_expression, eval_module) @@ -854,9 +856,8 @@ function compile_equational_affect( end end else - return let dvs_to_update = dvs_to_update, - affsys = affsys, ps_to_update = ps_to_update, aff = aff, sys = sys, - reset_jumps = reset_jumps + return let dvs_to_update = dvs_to_update, affsys = affsys, + ps_to_update = ps_to_update, aff = aff, sys = sys, reset_jumps = reset_jumps dvs_to_access = unknowns(affsys) ps_to_access = [unPre(p) for p in parameters(affsys)] diff --git a/src/systems/clock_inference.jl b/src/systems/clock_inference.jl index 42fe28f7c7..97b6be27ab 100644 --- a/src/systems/clock_inference.jl +++ b/src/systems/clock_inference.jl @@ -202,7 +202,8 @@ function split_system(ci::ClockInference{S}) where {S} if continuous_id != 0 tss[continuous_id], tss[end] = tss[end], tss[continuous_id] inputs[continuous_id], inputs[end] = inputs[end], inputs[continuous_id] - id_to_clock[continuous_id], id_to_clock[end] = id_to_clock[end], + id_to_clock[continuous_id], + id_to_clock[end] = id_to_clock[end], id_to_clock[continuous_id] continuous_id = lastindex(tss) end diff --git a/src/systems/codegen.jl b/src/systems/codegen.jl index 17ff652c41..e7f1b1d75f 100644 --- a/src/systems/codegen.jl +++ b/src/systems/codegen.jl @@ -753,7 +753,8 @@ function generate_constraint_jacobian( simplify = false, sparse = false, kwargs...) dvs = unknowns(sys) ps = reorder_parameters(sys) - jac, sparsity = calculate_constraint_jacobian( + jac, + sparsity = calculate_constraint_jacobian( sys; simplify, sparse, return_sparsity = true) res = build_function_wrapper(sys, jac, dvs, ps...; expression = Val{true}, kwargs...) fn = maybe_compile_function( @@ -807,7 +808,8 @@ function generate_constraint_hessian( simplify = false, sparse = false, kwargs...) dvs = unknowns(sys) ps = reorder_parameters(sys) - hess, sparsity = calculate_constraint_hessian( + hess, + sparsity = calculate_constraint_hessian( sys; simplify, sparse, return_sparsity = true) res = build_function_wrapper(sys, hess, dvs, ps...; expression = Val{true}, kwargs...) fn = maybe_compile_function( diff --git a/src/systems/connectors.jl b/src/systems/connectors.jl index cf2f24a7d1..46a2191d96 100644 --- a/src/systems/connectors.jl +++ b/src/systems/connectors.jl @@ -728,6 +728,7 @@ function Base.merge(csets::AbstractVector{<:ConnectionSet}, allouter = false) union_find = IntDisjointSets(0) prev_id = Ref(-1) for cset in csets, (j, s) in enumerate(cset.set) + v = allouter ? withtrueouter(s) : s id = let ele2idx = ele2idx, idx2ele = idx2ele get!(ele2idx, v) do @@ -810,7 +811,8 @@ function domain_defaults(sys, domain_csets) for p in parameters(m.sys.namespace) d_p = get(ns_s_def, p, nothing) if d_p !== nothing - def[parameters(m.sys.namespace, p)] = parameters(s.sys.namespace, + def[parameters(m.sys.namespace, + p)] = parameters(s.sys.namespace, parameters(s.sys.sys, d_p)) end diff --git a/src/systems/diffeqs/basic_transformations.jl b/src/systems/diffeqs/basic_transformations.jl index 430260f60a..0dddfcd309 100644 --- a/src/systems/diffeqs/basic_transformations.jl +++ b/src/systems/diffeqs/basic_transformations.jl @@ -439,7 +439,8 @@ All accumulation variables have a default of zero. function add_accumulations(sys::System, vars::Vector{<:Pair}) eqs = get_eqs(sys) avars = map(first, vars) - if (ints = intersect(avars, unknowns(sys)); !isempty(ints)) + ints = intersect(avars, unknowns(sys)) + if !isempty(ints) error("$ints already exist in the system!") end D = Differential(get_iv(sys)) diff --git a/src/systems/if_lifting.jl b/src/systems/if_lifting.jl index c2eba035ff..aba9dbdcbf 100644 --- a/src/systems/if_lifting.jl +++ b/src/systems/if_lifting.jl @@ -375,34 +375,69 @@ end const CONDITION_SIMPLIFIER = Rewriters.Fixpoint(Rewriters.Postwalk(Rewriters.Chain([ # simple boolean laws (@rule (!!(~x)) => (~x)) - (@rule ((~x) & true) => (~x)) - (@rule ((~x) & false) => false) - (@rule ((~x) | true) => true) - (@rule ((~x) | false) => (~x)) - (@rule ((~x) & !(~x)) => false) - (@rule ((~x) | !(~x)) => true) + (@rule ((~x) & + true) => (~x)) + (@rule ((~x) & + false) => false) + (@rule ((~x) | + true) => true) + (@rule ((~x) | + false) => (~x)) + (@rule ((~x) & + !(~x)) => false) + (@rule ((~x) | + !(~x)) => true) # reversed order of the above, because it matters and `@acrule` refuses to do its job - (@rule (true & (~x)) => (~x)) - (@rule (false & (~x)) => false) - (@rule (true | (~x)) => true) - (@rule (false | (~x)) => (~x)) - (@rule (!(~x) & (~x)) => false) - (@rule (!(~x) | (~x)) => true) + (@rule (true & + (~x)) => (~x)) + (@rule (false & + (~x)) => false) + (@rule (true | + (~x)) => true) + (@rule (false | + (~x)) => (~x)) + (@rule (!(~x) & + (~x)) => false) + (@rule (!(~x) | + (~x)) => true) # idempotent - (@rule ((~x) & (~x)) => (~x)) - (@rule ((~x) | (~x)) => (~x)) + (@rule ((~x) & + (~x)) => (~x)) + (@rule ((~x) | + (~x)) => (~x)) # ifelse with determined branches - (@rule ifelse((~x), true, false) => (~x)) - (@rule ifelse((~x), false, true) => !(~x)) + (@rule ifelse( + (~x), + true, + false) => (~x)) + (@rule ifelse( + (~x), + false, + true) => !(~x)) # ifelse with identical branches - (@rule ifelse((~x), (~y), (~y)) => (~y)) - (@rule ifelse((~x), (~y), !(~y)) => ((~x) & - (~y))) - (@rule ifelse((~x), !(~y), (~y)) => ((~x) & - !(~y))) + (@rule ifelse( + (~x), + (~y), + (~y)) => (~y)) + (@rule ifelse( + (~x), + (~y), + !(~y)) => ((~x) & + (~y))) + (@rule ifelse( + (~x), + !(~y), + (~y)) => ((~x) & + !(~y))) # ifelse with determined condition - (@rule ifelse(true, (~x), (~y)) => (~x)) - (@rule ifelse(false, (~x), (~y)) => (~y))]))) + (@rule ifelse( + true, + (~x), + (~y)) => (~x)) + (@rule ifelse( + false, + (~x), + (~y)) => (~y))]))) """ If lifting converts (nested) if statements into a series of continuous events + a logically equivalent if statement + parameters. @@ -499,7 +534,8 @@ function IfLifting(sys::System) for var in new_cond_vars condition = generate_condition(cw, var) - up_affect, down_affect = generate_affects( + up_affect, + down_affect = generate_affects( cw, var, new_cond_vars, new_cond_vars_graph) cb = SymbolicContinuousCallback([condition], up_affect; affect_neg = down_affect, initialize = up_affect, rootfind = SciMLBase.RightRootFind) diff --git a/src/systems/index_cache.jl b/src/systems/index_cache.jl index c66c562e9c..29e6f2db89 100644 --- a/src/systems/index_cache.jl +++ b/src/systems/index_cache.jl @@ -246,9 +246,11 @@ function IndexCache(sys::AbstractSystem) return idxs, buffer_sizes end - const_idxs, const_buffer_sizes = get_buffer_sizes_and_idxs( + const_idxs, + const_buffer_sizes = get_buffer_sizes_and_idxs( ParamIndexMap, constant_buffers) - nonnumeric_idxs, nonnumeric_buffer_sizes = get_buffer_sizes_and_idxs( + nonnumeric_idxs, + nonnumeric_buffer_sizes = get_buffer_sizes_and_idxs( NonnumericMap, nonnumeric_buffers) tunable_idxs = TunableIndexMap() diff --git a/src/systems/model_parsing.jl b/src/systems/model_parsing.jl index 9d293fd40d..33316f1cf0 100644 --- a/src/systems/model_parsing.jl +++ b/src/systems/model_parsing.jl @@ -85,12 +85,17 @@ function _model_macro(mod, fullname::Union{Expr, Symbol}, expr, isconnector) push!(exprs.args, arg) elseif arg.head == :if MLStyle.@match arg begin - Expr(:if, condition, x) => begin + Expr(:if, + condition, + x) => begin parse_conditional_model_statements(comps, dict, eqs, exprs, kwargs, mod, ps, vs, where_types, parse_top_level_branch(condition, x.args)...) end - Expr(:if, condition, x, y) => begin + Expr(:if, + condition, + x, + y) => begin parse_conditional_model_statements(comps, dict, eqs, exprs, kwargs, mod, ps, vs, where_types, parse_top_level_branch(condition, x.args, y)...) @@ -273,7 +278,9 @@ Base.@nospecializeinfer function parse_variable_def!( # Parses: `par6(t)[1:3]::BigFloat` # Recursively called by: `par2(t)::Int` # Recursively called by: `par3(t)::BigFloat = 1.0` - Expr(:(::), a, type) => begin + Expr(:(::), + a, + type) => begin type = getfield(mod, type) parse_variable_def!( dict, mod, a, varclass, kwargs, where_types; def, type, meta) @@ -283,7 +290,9 @@ Base.@nospecializeinfer function parse_variable_def!( # Recursively called by: `j(t) = jval, [description = "j(t)"]` # Recursively called by: `par2(t)::Int` # Recursively called by: `par3(t)::BigFloat = 1.0` - Expr(:call, a, b) => begin + Expr(:call, + a, + b) => begin var = generate_var!(dict, a, b, varclass, mod; type) update_kwargs_and_metadata!(dict, kwargs, a, def, type, varclass, where_types, meta) @@ -307,7 +316,8 @@ Base.@nospecializeinfer function parse_variable_def!( # `e2[1:2], [description = "e2"]` # `h2(t)[1:2], [description = "h2(t)"]` Expr(:tuple, Expr(:(=), Expr(:ref, a, indices...), default_val), meta_val) || - Expr(:tuple, Expr(:(=), Expr(:(::), Expr(:ref, a, indices...), type), default_val), meta_val) || + Expr(:tuple, Expr(:(=), Expr(:(::), Expr(:ref, a, indices...), type), default_val), + meta_val) || Expr(:tuple, Expr(:(::), Expr(:ref, a, indices...), type), meta_val) || Expr(:tuple, Expr(:ref, a, indices...), meta_val) => begin (@isdefined type) || (type = Real) @@ -341,7 +351,9 @@ Base.@nospecializeinfer function parse_variable_def!( # `d2[1:2] = 2` # `l(t)[1:2, 1:3] = 2, [description = "l is more than 1D"]` Expr(:(=), Expr(:(::), Expr(:ref, a, indices...), type), def_n_meta) || - Expr(:(=), Expr(:ref, a, indices...), def_n_meta) => begin + Expr(:(=), + Expr(:ref, a, indices...), + def_n_meta) => begin (@isdefined type) || (type = Real) varname = Meta.isexpr(a, :call) ? a.args[1] : a if Meta.isexpr(def_n_meta, :tuple) @@ -401,7 +413,9 @@ Base.@nospecializeinfer function parse_variable_def!( # `b2(t)[1:2]` # `a2[1:2]` Expr(:(::), Expr(:ref, a, indices...), type) || - Expr(:ref, a, indices...) => begin + Expr(:ref, + a, + indices...) => begin (@isdefined type) || (type = Real) varname = a isa Expr && a.head == :call ? a.args[1] : a if varclass == :parameters @@ -431,10 +445,13 @@ Base.@nospecializeinfer function parse_variable_def!( # Parses: `k = kval, [description = "k"]` # Parses: `par0::Bool = true` # Parses: `par3(t)::BigFloat = 1.0` - Expr(:(=), a, b) => begin + Expr(:(=), + a, + b) => begin Base.remove_linenums!(b) def, meta = parse_default(mod, b) - var, def, _ = parse_variable_def!( + var, def, + _ = parse_variable_def!( dict, mod, a, varclass, kwargs, where_types; def, type, meta) varclass_dict = dict[varclass] isa Vector ? Ref(dict[varclass][1]) : Ref(dict[varclass]) @@ -449,9 +466,12 @@ Base.@nospecializeinfer function parse_variable_def!( # Parses: `e, [description = "e"]` # Parses: `h(t), [description = "h(t)"]` # Parses: `par2(t)::Int` - Expr(:tuple, a, b) => begin + Expr(:tuple, + a, + b) => begin meta = parse_metadata(mod, b) - var, def, _ = parse_variable_def!( + var, def, + _ = parse_variable_def!( dict, mod, a, varclass, kwargs, where_types; type, meta) varclass_dict = dict[varclass] isa Vector ? Ref(dict[varclass][1]) : Ref(dict[varclass]) @@ -690,7 +710,9 @@ function parse_structural_parameters!(exprs, sps, dict, mod, body, kwargs) Base.remove_linenums!(body) for arg in body.args MLStyle.@match arg begin - Expr(:(=), Expr(:(::), a, type), b) => begin + Expr(:(=), + Expr(:(::), a, type), + b) => begin type = getfield(mod, type) b = _type_check!(get_var(mod, b), a, type, :structural_parameters) push!(sps, a) @@ -698,7 +720,9 @@ function parse_structural_parameters!(exprs, sps, dict, mod, body, kwargs) dict[:structural_parameters][a] = dict[:kwargs][a] = Dict( :value => b, :type => type) end - Expr(:(=), a, b) => begin + Expr(:(=), + a, + b) => begin push!(sps, a) push!(kwargs, Expr(:kw, a, b)) dict[:structural_parameters][a] = dict[:kwargs][a] = Dict(:value => b) @@ -742,7 +766,8 @@ function extend_args!(a, b, dict, expr, kwargs, has_param = false) push!(kwargs, Expr(:kw, x, y)) dict[:kwargs][x] = Dict(:value => y) end - Expr(:parameters, x...) => begin + Expr(:parameters, + x...) => begin has_param = true extend_args!(a, arg, dict, expr, kwargs, has_param) end @@ -806,7 +831,9 @@ function parse_extend!(exprs, ext, dict, mod, body, kwargs) push!(exprs, expr) body = deepcopy(body) MLStyle.@match body begin - Expr(:(=), a, b) => begin + Expr(:(=), + a, + b) => begin if Meta.isexpr(b, :(=)) vars = a if !Meta.isexpr(vars, :tuple) @@ -820,7 +847,9 @@ function parse_extend!(exprs, ext, dict, mod, body, kwargs) error("When explicitly destructing in `@extend` please use the syntax: `@extend a, b = oneport = OnePort()`.") end end - Expr(:call, a′, _...) => begin + Expr(:call, + a′, + _...) => begin a = Symbol(Symbol("#mtkmodel"), :__anonymous__, a′) b = body if (model = getproperty(mod, b.args[1])) isa Model @@ -876,7 +905,8 @@ convert_units(::Unitful.FreeUnits, value::Num) = value convert_units(::DynamicQuantities.Quantity, value::Num) = value function parse_variable_arg(dict, mod, arg, varclass, kwargs, where_types) - vv, def, metadata_with_exprs = parse_variable_def!( + vv, def, + metadata_with_exprs = parse_variable_def!( dict, mod, arg, varclass, kwargs, where_types) if !(vv isa Tuple) name = getname(vv) @@ -931,7 +961,8 @@ function handle_conditional_vars!( :constants => Any[Dict{Symbol, Dict{Symbol, Any}}()], :variables => Any[Dict{Symbol, Dict{Symbol, Any}}()]) for _arg in arg.args - name, ex = parse_variable_arg( + name, + ex = parse_variable_arg( conditional_dict, mod, _arg, varclass, kwargs, where_types) push!(conditional_branch.args, ex) push!(conditional_branch.args, :(push!($varclass, $name))) @@ -997,7 +1028,9 @@ function parse_variables!(exprs, vs, dict, mod, body, varclass, kwargs, where_ty for arg in body.args arg isa LineNumberNode && continue MLStyle.@match arg begin - Expr(:if, condition, x) => begin + Expr(:if, + condition, + x) => begin conditional_expr = Expr(:if, condition, Expr(:block)) conditional_dict = handle_conditional_vars!(x, conditional_expr.args[2], @@ -1008,7 +1041,10 @@ function parse_variables!(exprs, vs, dict, mod, body, varclass, kwargs, where_ty push!(expr.args, conditional_expr) push_conditional_dict!(dict, condition, conditional_dict, nothing, varclass) end - Expr(:if, condition, x, y) => begin + Expr(:if, + condition, + x, + y) => begin conditional_expr = Expr(:if, condition, Expr(:block)) conditional_dict = handle_conditional_vars!(x, conditional_expr.args[2], @@ -1016,7 +1052,8 @@ function parse_variables!(exprs, vs, dict, mod, body, varclass, kwargs, where_ty varclass, kwargs, where_types) - conditional_y_expr, conditional_y_tuple = handle_y_vars(y, + conditional_y_expr, + conditional_y_tuple = handle_y_vars(y, conditional_dict, mod, varclass, @@ -1044,7 +1081,8 @@ function handle_y_vars(y, dict, mod, varclass, kwargs, where_types) varclass, kwargs, where_types) - _y_expr, _conditional_dict = handle_y_vars( + _y_expr, + _conditional_dict = handle_y_vars( y.args[end], dict, mod, varclass, kwargs, where_types) push!(conditional_y_expr.args, _y_expr) (:elseif, y.args[1], conditional_dict, _conditional_dict) @@ -1090,13 +1128,15 @@ function parse_equations!(exprs, eqs, dict, body) Base.remove_linenums!(body) for arg in body.args MLStyle.@match arg begin - Expr(:if, condition, x) => begin + Expr(:if, condition, + x) => begin ifexpr = Expr(:if) eq_entry = handle_if_x_equations!(condition, dict, ifexpr, x) push!(exprs, ifexpr) push!(dict[:equations], (:if, condition, eq_entry)) end - Expr(:if, condition, x, y) => begin + Expr(:if, condition, x, + y) => begin ifexpr = Expr(:if) xeq_entry = handle_if_x_equations!(condition, dict, ifexpr, x) yeq_entry = handle_if_y_equations!(ifexpr, y, dict) @@ -1181,19 +1221,25 @@ function parse_consolidate!(body, dict) end end +function _try_validate_uri(body::String) + try + Base.isvalid(URI(body)) + catch e + false + end +end + function parse_icon!(body::String, dict, icon, mod) icon_dir = get(ENV, "MTK_ICONS_DIR", joinpath(DEPOT_PATH[1], "mtk_icons")) + iconpath = abspath(joinpath(icon_dir, body)) + _body = lstrip(body) dict[:icon] = icon[] = if isfile(body) URI("file:///" * abspath(body)) - elseif (iconpath = abspath(joinpath(icon_dir, body)); isfile(iconpath)) + elseif isfile(iconpath) URI("file:///" * abspath(iconpath)) - elseif try - Base.isvalid(URI(body)) - catch e - false - end + elseif _try_validate_uri(body) URI(body) - elseif (_body = lstrip(body); startswith(_body, r"<\?xml| begin + x::Symbol || + Expr(:kw, x) => begin varname, _varname = _rename(a, x) b.args[i] = Expr(:kw, x, _varname) push!(varexpr.args, :((if $varname !== nothing @@ -1240,7 +1287,9 @@ function component_args!(a, b, varexpr, kwargs; index_name = nothing) Expr(:parameters, x...) => begin component_args!(a, arg, varexpr, kwargs) end - Expr(:kw, x, y) => begin + Expr(:kw, + x, + y) => begin varname, _varname = _rename(a, x) b.args[i] = Expr(:kw, x, _varname) if isnothing(index_name) @@ -1269,7 +1318,9 @@ function _parse_components!(body, kwargs) arg = body.args[end] MLStyle.@match arg begin - Expr(:(=), a, Expr(:comprehension, Expr(:generator, b, Expr(:(=), c, d)))) => begin + Expr(:(=), + a, + Expr(:comprehension, Expr(:generator, b, Expr(:(=), c, d)))) => begin array_varexpr = Expr(:block) push!(comp_names, :($a...)) @@ -1280,10 +1331,14 @@ function _parse_components!(body, kwargs) expr = _named_idxs(a, d, :($c -> $b); extra_args = array_varexpr) end - Expr(:(=), a, Expr(:comprehension, Expr(:generator, b, Expr(:filter, e, Expr(:(=), c, d))))) => begin + Expr(:(=), + a, + Expr(:comprehension, Expr(:generator, b, Expr(:filter, e, Expr(:(=), c, d))))) => begin error("List comprehensions with conditional statements aren't supported.") end - Expr(:(=), a, Expr(:comprehension, Expr(:generator, b, Expr(:(=), c, d), e...))) => begin + Expr(:(=), + a, + Expr(:comprehension, Expr(:generator, b, Expr(:(=), c, d), e...))) => begin # Note that `e` is of the form `Tuple{Expr(:(=), c, d)}` error("More than one index isn't supported while building component array") end @@ -1291,7 +1346,9 @@ function _parse_components!(body, kwargs) # TODO: Do we need this? error("Multiple `@components` block detected within a single block") end - Expr(:(=), a, Expr(:for, Expr(:(=), c, d), b)) => begin + Expr(:(=), + a, + Expr(:for, Expr(:(=), c, d), b)) => begin Base.remove_linenums!(b) array_varexpr = Expr(:block) push!(array_varexpr.args, b.args[1:(end - 1)]...) @@ -1367,15 +1424,20 @@ function parse_components!(exprs, cs, dict, compbody, kwargs) Base.remove_linenums!(compbody) for arg in compbody.args MLStyle.@match arg begin - Expr(:if, condition, x) => begin + Expr(:if, condition, + x) => begin handle_conditional_components(condition, dict, exprs, kwargs, x) end - Expr(:if, condition, x, y) => begin + Expr(:if, + condition, + x, + y) => begin handle_conditional_components(condition, dict, exprs, kwargs, x, y) end # Either the arg is top level component declaration or an invalid cause - both are handled by `_parse_components` _ => begin - comp_names, comps, expr_vec, varexpr = _parse_components!(:(begin + comp_names, comps, expr_vec, + varexpr = _parse_components!(:(begin $arg end), kwargs) @@ -1403,7 +1465,9 @@ define_blocks(branch) = [Expr(branch), Expr(branch), Expr(branch), Expr(branch)] Base.@nospecializeinfer function parse_top_level_branch( condition, x, y = nothing, branch::Symbol = :if) @nospecialize - blocks::Vector{Union{Expr, Nothing}} = component_blk, equations_blk, parameter_blk, variable_blk = define_blocks(branch) + blocks::Vector{Union{Expr, + Nothing}} = component_blk, equations_blk, parameter_blk, + variable_blk = define_blocks(branch) for arg in x if arg.args[1] == Symbol("@components") diff --git a/src/systems/nonlinear/homotopy_continuation.jl b/src/systems/nonlinear/homotopy_continuation.jl index 77e39e4839..96c00411ad 100644 --- a/src/systems/nonlinear/homotopy_continuation.jl +++ b/src/systems/nonlinear/homotopy_continuation.jl @@ -543,7 +543,8 @@ function HomotopyContinuationProblem{iip, spec}( if !iscomplete(sys) error("A completed `System` is required. Call `complete` or `mtkcompile` on the system before creating a `HomotopyContinuationProblem`") end - f, u0, p = process_SciMLProblem( + f, u0, + p = process_SciMLProblem( HomotopyNonlinearFunction{iip, spec}, sys, op; kwargs...) kwargs = filter_kwargs(kwargs) diff --git a/src/systems/nonlinear/initializesystem.jl b/src/systems/nonlinear/initializesystem.jl index 5dfd69b72f..81c1d0160c 100644 --- a/src/systems/nonlinear/initializesystem.jl +++ b/src/systems/nonlinear/initializesystem.jl @@ -544,7 +544,8 @@ function SciMLBase.remake_initialization_data( # is called is because `Initial` parameters are calculated from the corresponding # state values. history_fn = is_time_dependent(sys) && !is_markovian(sys) ? Returns(newu0) : nothing - new_initu0, new_initp = reconstruct_fn( + new_initu0, + new_initp = reconstruct_fn( ProblemState(; u = newu0, p = newp, t = t0, h = history_fn), oldinitprob) if oldinitprob.f.resid_prototype === nothing newf = oldinitprob.f @@ -621,7 +622,8 @@ function SciMLBase.remake_initialization_data( u0map = anydict() pmap = anydict() - missing_unknowns, missing_pars = build_operating_point!(sys, op, + missing_unknowns, + missing_pars = build_operating_point!(sys, op, u0map, pmap, defs, dvs, ps) floatT = float_type_from_varmap(op) u0_constructor = p_constructor = identity diff --git a/src/systems/optimal_control_interface.jl b/src/systems/optimal_control_interface.jl index ea19afea86..ba4234727b 100644 --- a/src/systems/optimal_control_interface.jl +++ b/src/systems/optimal_control_interface.jl @@ -93,7 +93,7 @@ function PyomoCollocation end function warn_overdetermined(sys, op) cstrs = constraints(sys) - init_conds = filter(x -> value(x) ∈ Set(unknowns(sys)), [k for (k,v) in op]) + init_conds = filter(x -> value(x) ∈ Set(unknowns(sys)), [k for (k, v) in op]) if !isempty(cstrs) (length(cstrs) + length(init_conds) > length(unknowns(sys))) && @warn "The control problem is overdetermined. The total number of conditions (# constraints + # fixed initial values given by op) exceeds the total number of states. The solvers will default to doing a nonlinear least-squares optimization." @@ -135,7 +135,8 @@ is_explicit(tableau) = tableau isa DiffEqBase.ExplicitRKTableau initialization_data = nothing, cse = true, kwargs...) where {iip, specialize} - f, _, _ = generate_control_function( + f, _, + _ = generate_control_function( sys, inputs, disturbance_inputs; eval_module, cse, kwargs...) f = f[1] @@ -227,11 +228,11 @@ end ########################## ### MODEL CONSTRUCTION ### ########################## -function process_DynamicOptProblem(prob_type::Type{<:AbstractDynamicOptProblem}, model_type, sys::System, op, tspan; - dt = nothing, - steps = nothing, - guesses = Dict(), kwargs...) - +function process_DynamicOptProblem( + prob_type::Type{<:AbstractDynamicOptProblem}, model_type, sys::System, op, tspan; + dt = nothing, + steps = nothing, + guesses = Dict(), kwargs...) warn_overdetermined(sys, op) ctrls = unbound_inputs(sys) states = unknowns(sys) @@ -239,10 +240,11 @@ function process_DynamicOptProblem(prob_type::Type{<:AbstractDynamicOptProblem}, stidxmap = Dict([v => i for (i, v) in enumerate(states)]) op = Dict([default_toterm(value(k)) => v for (k, v) in op]) u0_idxs = has_alg_eqs(sys) ? collect(1:length(states)) : - [stidxmap[default_toterm(k)] for (k, v) in op if haskey(stidxmap, k)] + [stidxmap[default_toterm(k)] for (k, v) in op if haskey(stidxmap, k)] _op = has_alg_eqs(sys) ? op : merge(Dict(op), Dict(guesses)) - f, u0, p = process_SciMLProblem(ODEInputFunction, sys, _op; + f, u0, + p = process_SciMLProblem(ODEInputFunction, sys, _op; t = tspan !== nothing ? tspan[1] : tspan, kwargs...) model_tspan, steps, is_free_t = process_tspan(tspan, dt, steps) warn_overdetermined(sys, op) @@ -302,7 +304,7 @@ function set_variable_bounds!(m, sys, pmap, tf) end end -is_free_final(model) = model.is_free_final +is_free_final(model) = model.is_free_final function add_cost_function!(model, sys, tspan, pmap) jcosts = cost(sys) @@ -335,14 +337,15 @@ function substitute_integral(model, expr, tspan) Symbolics.substitute(expr, intmap) end -function process_integral_bounds(model, integral_span, tspan) +function process_integral_bounds(model, integral_span, tspan) if is_free_final(model) && isequal(integral_span, tspan) integral_span = (0, 1) elseif is_free_final(model) error("Free final time problems cannot handle partial timespans.") else (lo, hi) = integral_span - (lo < tspan[1] || hi > tspan[2]) && error("Integral bounds are beyond the timespan.") + (lo < tspan[1] || hi > tspan[2]) && + error("Integral bounds are beyond the timespan.") integral_span end end @@ -353,12 +356,14 @@ function substitute_model_vars(model, sys, exprs, tspan) c_ops = [operation(unwrap(ct)) for ct in unbound_inputs(sys)] t = get_iv(sys) - exprs = map(c -> Symbolics.fast_substitute(c, whole_t_map(model, t, x_ops, c_ops)), exprs) + exprs = map( + c -> Symbolics.fast_substitute(c, whole_t_map(model, t, x_ops, c_ops)), exprs) (ti, tf) = tspan if symbolic_type(tf) === ScalarSymbolic() _tf = model.tₛ + ti - exprs = map(c -> Symbolics.fast_substitute(c, free_t_map(model, tf, x_ops, c_ops)), exprs) + exprs = map( + c -> Symbolics.fast_substitute(c, free_t_map(model, tf, x_ops, c_ops)), exprs) exprs = map(c -> Symbolics.fast_substitute(c, Dict(tf => _tf)), exprs) end exprs = map(c -> Symbolics.fast_substitute(c, fixed_t_map(model, x_ops, c_ops)), exprs) @@ -392,7 +397,9 @@ function fixed_t_map end function add_user_constraints!(model, sys, tspan, pmap) jconstraints = get_constraints(sys) (isnothing(jconstraints) || isempty(jconstraints)) && return nothing - cons_dvs, cons_ps = process_constraint_system(jconstraints, Set(unknowns(sys)), parameters(sys), get_iv(sys); validate = false) + cons_dvs, + cons_ps = process_constraint_system( + jconstraints, Set(unknowns(sys)), parameters(sys), get_iv(sys); validate = false) is_free_final(model) && check_constraint_vars(cons_dvs) @@ -421,12 +428,13 @@ function add_equational_constraints!(model, sys, pmap, tspan) end function set_objective! end -objective_value(sol::DynamicOptSolution) = objective_value(sol.model) +objective_value(sol::DynamicOptSolution) = objective_value(sol.model) function substitute_differentials(model, sys, eqs) t = get_iv(sys) D = Differential(t) - diffsubmap = Dict([D(lowered_var(model, :U, i, t)) => lowered_derivative(model, i) for i in 1:length(unknowns(sys))]) + diffsubmap = Dict([D(lowered_var(model, :U, i, t)) => lowered_derivative(model, i) + for i in 1:length(unknowns(sys))]) eqs = map(c -> Symbolics.substitute(c, diffsubmap), eqs) end @@ -466,9 +474,10 @@ function successful_solve end - kwargs are used for other options. For example, the `plugin_options` and `solver_options` will propagated to the Opti object in CasADi. """ -function DiffEqBase.solve(prob::AbstractDynamicOptProblem, solver::AbstractCollocation; verbose = false, kwargs...) +function DiffEqBase.solve(prob::AbstractDynamicOptProblem, + solver::AbstractCollocation; verbose = false, kwargs...) solved_model = prepare_and_optimize!(prob, solver; verbose, kwargs...) - + ts = get_t_values(solved_model) Us = get_U_values(solved_model) Vs = get_V_values(solved_model) diff --git a/src/systems/parameter_buffer.jl b/src/systems/parameter_buffer.jl index bd77c8c519..b2f356002b 100644 --- a/src/systems/parameter_buffer.jl +++ b/src/systems/parameter_buffer.jl @@ -49,7 +49,8 @@ function MTKParameters( u0map = anydict() pmap = anydict() - missing_unknowns, missing_pars = build_operating_point!(sys, op, + missing_unknowns, + missing_pars = build_operating_point!(sys, op, u0map, pmap, defs, dvs, ps) if t0 !== nothing @@ -313,9 +314,9 @@ function SciMLStructures.replace!(::SciMLStructures.Initials, p::MTKParameters, end for (Portion, field, recurse) in [(SciMLStructures.Discrete, :discrete, 1) - (SciMLStructures.Constants, :constant, 1) - (Nonnumeric, :nonnumeric, 1) - (SciMLStructures.Caches, :caches, 1)] + (SciMLStructures.Constants, :constant, 1) + (Nonnumeric, :nonnumeric, 1) + (SciMLStructures.Caches, :caches, 1)] @eval function SciMLStructures.canonicalize(::$Portion, p::MTKParameters) as_vector = buffer_to_arraypartition(p.$field) repack = let p = p @@ -617,14 +618,20 @@ end initials = $similar(oldbuf.initials, $initialsT) copyto!(initials, oldbuf.initials) discretes = $(Expr(:tuple, - (:($similar(oldbuf.discrete[$i], $(discretesT[i]))) for i in 1:length(discretesT))...)) - $((:($copyto!(discretes[$i], oldbuf.discrete[$i])) for i in 1:length(discretesT))...) + (:($similar(oldbuf.discrete[$i], $(discretesT[i]))) + for i in 1:length(discretesT))...)) + $((:($copyto!(discretes[$i], oldbuf.discrete[$i])) + for i in 1:length(discretesT))...) constants = $(Expr(:tuple, - (:($similar(oldbuf.constant[$i], $(constantsT[i]))) for i in 1:length(constantsT))...)) - $((:($copyto!(constants[$i], oldbuf.constant[$i])) for i in 1:length(constantsT))...) + (:($similar(oldbuf.constant[$i], $(constantsT[i]))) + for i in 1:length(constantsT))...)) + $((:($copyto!(constants[$i], oldbuf.constant[$i])) + for i in 1:length(constantsT))...) nonnumerics = $(Expr(:tuple, - (:($similar(oldbuf.nonnumeric[$i], $(nonnumericT[i]))) for i in 1:length(nonnumericT))...)) - $((:($copyto!(nonnumerics[$i], oldbuf.nonnumeric[$i])) for i in 1:length(nonnumericT))...) + (:($similar(oldbuf.nonnumeric[$i], $(nonnumericT[i]))) + for i in 1:length(nonnumericT))...)) + $((:($copyto!(nonnumerics[$i], oldbuf.nonnumeric[$i])) + for i in 1:length(nonnumericT))...) caches = copy.(oldbuf.caches) newbuf = MTKParameters( tunables, initials, discretes, constants, nonnumerics, caches) @@ -643,13 +650,16 @@ end push!(expr.args, :(initials = $similar_type($I, $initialsT)(initials))) push!(expr.args, :(discretes = $(Expr(:tuple, - (:($similar_type($(fieldtype(D, i)), $(discretesT[i]))(discretes[$i])) for i in 1:length(discretesT))...)))) + (:($similar_type($(fieldtype(D, i)), $(discretesT[i]))(discretes[$i])) + for i in 1:length(discretesT))...)))) push!(expr.args, :(constants = $(Expr(:tuple, - (:($similar_type($(fieldtype(C, i)), $(constantsT[i]))(constants[$i])) for i in 1:length(constantsT))...)))) + (:($similar_type($(fieldtype(C, i)), $(constantsT[i]))(constants[$i])) + for i in 1:length(constantsT))...)))) push!(expr.args, :(nonnumerics = $(Expr(:tuple, - (:($similar_type($(fieldtype(C, i)), $(nonnumericT[i]))(nonnumerics[$i])) for i in 1:length(nonnumericT))...)))) + (:($similar_type($(fieldtype(C, i)), $(nonnumericT[i]))(nonnumerics[$i])) + for i in 1:length(nonnumericT))...)))) push!(expr.args, :(newbuf = MTKParameters( tunables, initials, discretes, constants, nonnumerics, caches))) diff --git a/src/systems/problem_utils.jl b/src/systems/problem_utils.jl index 2016b1efd8..35f3b1ea16 100644 --- a/src/systems/problem_utils.jl +++ b/src/systems/problem_utils.jl @@ -1148,8 +1148,8 @@ function maybe_build_initialization_problem( return (; initialization_data = SciMLBase.OverrideInitData( - initializeprob, update_initializeprob!, initializeprobmap, - initializeprobpmap; metadata = meta, is_update_oop = Val(true))) + initializeprob, update_initializeprob!, initializeprobmap, + initializeprobpmap; metadata = meta, is_update_oop = Val(true))) end """ @@ -1251,7 +1251,8 @@ function process_SciMLProblem( u0map = anydict() pmap = anydict() - missing_unknowns, missing_pars = build_operating_point!(sys, op, + missing_unknowns, + missing_pars = build_operating_point!(sys, op, u0map, pmap, defs, dvs, ps) floatT = Bool @@ -1367,7 +1368,7 @@ function process_SciMLProblem( kwargs = merge(kwargs, (; resid_prototype = u0_constructor(calculate_resid_prototype( - length(eqs), u0, p)))) + length(eqs), u0, p)))) end f = constructor(sys; u0 = u0, p = p, @@ -1466,7 +1467,8 @@ function SymbolicTstops( end end rps = reorder_parameters(sys) - tstops, _ = build_function_wrapper(sys, tstops, + tstops, + _ = build_function_wrapper(sys, tstops, rps..., t0, t1; diff --git a/src/systems/system.jl b/src/systems/system.jl index cb330b382f..cda7d60b65 100644 --- a/src/systems/system.jl +++ b/src/systems/system.jl @@ -384,7 +384,8 @@ function System(eqs::Vector{Equation}, iv, dvs, ps, brownians = []; if length(unique_sysnames) != length(sysnames) throw(NonUniqueSubsystemsError(sysnames, unique_sysnames)) end - continuous_events, discrete_events = create_symbolic_events( + continuous_events, + discrete_events = create_symbolic_events( continuous_events, discrete_events, eqs, iv) if iv === nothing && (!isempty(continuous_events) || !isempty(discrete_events)) diff --git a/src/systems/systems.jl b/src/systems/systems.jl index 4a6c1ccbb4..e7d21fe30b 100644 --- a/src/systems/systems.jl +++ b/src/systems/systems.jl @@ -109,6 +109,7 @@ function __mtkcompile(sys::AbstractSystem; simplify = false, dvar2eq[fullvars[dv]] = only(deqs) end for (j, bj) in enumerate(brown_vars), i in 𝑑neighbors(graph, bj) + push!(Is, i) push!(Js, j) eq = new_eqs[i] diff --git a/src/systems/systemstructure.jl b/src/systems/systemstructure.jl index 7e00f6ef95..a5a5920e56 100644 --- a/src/systems/systemstructure.jl +++ b/src/systems/systemstructure.jl @@ -19,7 +19,8 @@ using SparseArrays function quick_cancel_expr(expr) Rewriters.Postwalk(quick_cancel, - similarterm = (x, f, args; kws...) -> maketerm(typeof(x), f, args, + similarterm = (x, f, args; + kws...) -> maketerm(typeof(x), f, args, SymbolicUtils.metadata(x), kws...))(expr) end @@ -264,8 +265,8 @@ end function symbolic_contains(var, set) var in set || symbolic_type(var) == ArraySymbolic() && - Symbolics.shape(var) != Symbolics.Unknown() && - all(x -> x in set, Symbolics.scalarize(var)) + Symbolics.shape(var) != Symbolics.Unknown() && + all(x -> x in set, Symbolics.scalarize(var)) end function TearingState(sys; quick_cancel = false, check = true, sort_eqs = true) @@ -505,6 +506,7 @@ function TearingState(sys; quick_cancel = false, check = true, sort_eqs = true) # build incidence graph graph = BipartiteGraph(neqs, nvars, Val(false)) for (ie, vars) in enumerate(symbolic_incidence), v in vars + jv = var2idx[v] add_edge!(graph, ie, jv) end diff --git a/src/variables.jl b/src/variables.jl index 104616f6e4..8d93a00982 100644 --- a/src/variables.jl +++ b/src/variables.jl @@ -33,7 +33,8 @@ ModelingToolkit.dump_variable_metadata(p) """ function dump_variable_metadata(var) uvar = unwrap(var) - variable_source, name = Symbolics.getmetadata( + variable_source, + name = Symbolics.getmetadata( uvar, VariableSource, (:unknown, :unknown)) type = symtype(uvar) if type <: AbstractArray diff --git a/test/basic_transformations.jl b/test/basic_transformations.jl index 6880f77479..bafb5cf9e2 100644 --- a/test/basic_transformations.jl +++ b/test/basic_transformations.jl @@ -127,7 +127,7 @@ end p = [v => 10.0] prob = ODEProblem(Mx, [u0; p], (0.0, 20.0)) # 1 = dy/dx = (dy/dt)/(dx/dt) means equal initial horizontal and vertical velocities sol = solve(prob, Tsit5(); reltol = 1e-5) - @test all(isapprox.(sol[Mx.y], sol[Mx.x - g * (Mx.t)^2 / 2]; atol = 1e-10)) # compare to analytical solution (x(t) = v*t, y(t) = v*t - g*t^2/2) + @test all(isapprox.(sol[Mx.y], sol[Mx.x - g * (Mx.t) ^ 2 / 2]; atol = 1e-10)) # compare to analytical solution (x(t) = v*t, y(t) = v*t - g*t^2/2) end @testset "Change independent variable (free fall with 2nd order horizontal equation)" begin @@ -140,7 +140,7 @@ end u0 = [Mx.y => 0.0, Dx(Mx.y) => 1.0, Mx.t => 0.0, Mx.xˍt => 10.0] prob = ODEProblem(Mx, u0, (0.0, 20.0)) # 1 = dy/dx = (dy/dt)/(dx/dt) means equal initial horizontal and vertical velocities sol = solve(prob, Tsit5(); reltol = 1e-5) - @test all(isapprox.(sol[Mx.y], sol[Mx.x - g * (Mx.t)^2 / 2]; atol = 1e-10)) # compare to analytical solution (x(t) = v*t, y(t) = v*t - g*t^2/2) + @test all(isapprox.(sol[Mx.y], sol[Mx.x - g * (Mx.t) ^ 2 / 2]; atol = 1e-10)) # compare to analytical solution (x(t) = v*t, y(t) = v*t - g*t^2/2) end @testset "Change independent variable (crazy 3rd order nonlinear system)" begin @@ -158,8 +158,9 @@ end # Compare to pen-and-paper result @variables x xˍt(x) xˍt(x) y(x) t(x) Dx = Differential(x) - areequivalent(eq1, eq2) = isequal(expand(eq1.lhs - eq2.lhs), 0) && - isequal(expand(eq1.rhs - eq2.rhs), 0) + areequivalent(eq1, + eq2) = isequal(expand(eq1.lhs - eq2.lhs), 0) && + isequal(expand(eq1.rhs - eq2.rhs), 0) eq1lhs = xˍt^3 * (Dx^3)(y) + xˍt^2 * Dx(y) * (Dx^2)(xˍt) + xˍt * Dx(y) * (Dx(xˍt))^2 + 3 * xˍt^2 * (Dx^2)(y) * Dx(xˍt) @@ -194,8 +195,9 @@ end @test Set(equations(M2)) == Set([ t ~ √(x), xˍt ~ 2t, - xˍt * Dx(y) ~ 1fc(t) + 2fc(x) + 3fc(y) + - 1callme(f, t) + 2callme(f, x) + 3callme(f, y) + xˍt * Dx(y) ~ + 1fc(t) + 2fc(x) + 3fc(y) + + 1callme(f, t) + 2callme(f, x) + 3callme(f, y) ]) _f = LinearInterpolation([1.0, 1.0], [-100.0, +100.0]) # constant value 1 @@ -230,7 +232,7 @@ end prob = ODEProblem(Mx, u0, (0.0, 20.0)) # 1 = dy/dx = (dy/dt)/(dx/dt) means equal initial horizontal and vertical velocities sol = solve(prob, Tsit5(); reltol = 1e-5) # compare to analytical solution (x(t) = v*t, y(t) = v*t - g*t^2/2) - @test all(isapprox.(sol[Mx.y], sol[Mx.x - g * (Mx.t_units)^2 / 2]; atol = 1e-10)) + @test all(isapprox.(sol[Mx.y], sol[Mx.x - g * (Mx.t_units) ^ 2 / 2]; atol = 1e-10)) end @testset "Change independent variable, no equations" begin diff --git a/test/code_generation.jl b/test/code_generation.jl index 15de194fd8..ac5a3dec77 100644 --- a/test/code_generation.jl +++ b/test/code_generation.jl @@ -17,7 +17,8 @@ using ModelingToolkit: t_nounits as t, D_nounits as D sys, x + y[1] + p1 + p2[1] + p3 * t, [x], [p1, p2, p3]; expression = Val(false)) @test fn1(u0, p, 0.0) == 5.0 - fn3_oop, fn3_iip = generate_custom_function( + fn3_oop, + fn3_iip = generate_custom_function( sys, [x + y[2], y[3] + p2[2], p1 + p3, 3t]; expression = Val(false)) buffer = zeros(4) @@ -41,7 +42,8 @@ using ModelingToolkit: t_nounits as t, D_nounits as D sys, x + y[1] + p1 + p2[1] + p3, [x], [p1, p2, p3]; expression = Val(false)) @test fn1(u0, p) == 6.0 - fn3_oop, fn3_iip = generate_custom_function( + fn3_oop, + fn3_iip = generate_custom_function( sys, [x + y[2], y[3] + p2[2], p1 + p3]; expression = Val(false)) buffer = zeros(3) diff --git a/test/dde.jl b/test/dde.jl index fb92bf61eb..3fd5aa94d6 100644 --- a/test/dde.jl +++ b/test/dde.jl @@ -35,8 +35,9 @@ sol2 = solve(prob2, alg, reltol = 1e-7, abstol = 1e-10) @variables x₀(t) x₁(t) x₂(..) tau = 1 eqs = [D(x₀) ~ (v0 / (1 + beta0 * (x₂(t - tau)^2))) * (p0 - q0) * x₀ - d0 * x₀ - D(x₁) ~ (v0 / (1 + beta0 * (x₂(t - tau)^2))) * (1 - p0 + q0) * x₀ + - (v1 / (1 + beta1 * (x₂(t - tau)^2))) * (p1 - q1) * x₁ - d1 * x₁ + D(x₁) ~ + (v0 / (1 + beta0 * (x₂(t - tau)^2))) * (1 - p0 + q0) * x₀ + + (v1 / (1 + beta1 * (x₂(t - tau)^2))) * (p1 - q1) * x₁ - d1 * x₁ D(x₂(t)) ~ (v1 / (1 + beta1 * (x₂(t - tau)^2))) * (1 - p1 + q1) * x₁ - d2 * x₂(t)] @mtkcompile sys = System(eqs, t) @test ModelingToolkit.is_dde(sys) @@ -174,7 +175,7 @@ end @testset "Issue#3165 DDEs with non-tunables" begin @variables x(..) = 1.0 - @parameters w=1.0 [tunable = false] τ=0.5 + @parameters w=1.0 [tunable=false] τ=0.5 eqs = [D(x(t)) ~ -w * x(t - τ)] @named sys = System(eqs, t) diff --git a/test/dep_graphs.jl b/test/dep_graphs.jl index 04e6bf159a..1fa166f1b7 100644 --- a/test/dep_graphs.jl +++ b/test/dep_graphs.jl @@ -70,21 +70,21 @@ import ModelingToolkit: value @testset "Case $i" for (i, test_case) in enumerate([test_case_1, test_case_2]) (; # filter out vrjs in making graphs - eqs, # eq to vars they depend on - eq_sdeps, - eq_sidepsf, - eq_sidepsb, # eq to params they depend on - eq_pdeps, - eq_pidepsf, - eq_pidepsb, # var to eqs that modify them - s_eqdepsf, - s_eqdepsb, - var_eq_ne, # eq to eqs that depend on them - eq_eqdeps, - eq_eq_ne, # var to vars that depend on them - var_vardeps, - var_var_ne -) = test_case + eqs, # eq to vars they depend on + eq_sdeps, + eq_sidepsf, + eq_sidepsb, # eq to params they depend on + eq_pdeps, + eq_pidepsf, + eq_pidepsb, # var to eqs that modify them + s_eqdepsf, + s_eqdepsb, + var_eq_ne, # eq to eqs that depend on them + eq_eqdeps, + eq_eq_ne, # var to vars that depend on them + var_vardeps, + var_var_ne + ) = test_case deps = equation_dependencies(js; eqs) @test length(deps) == length(eq_sdeps) @test all([issetequal(a, b) for (a, b) in zip(eq_sdeps, deps)]) diff --git a/test/direct.jl b/test/direct.jl index 70e1babe3f..ce2a3f5785 100644 --- a/test/direct.jl +++ b/test/direct.jl @@ -118,7 +118,8 @@ Jiip(J2, [1.0, 2.0, 3.0], [1.0, 2.0, 3.0], 1.0) s∂ = sparse(∂) @test nnz(s∂) == 8 -Joop, Jiip = eval.(ModelingToolkit.build_function(s∂, [x, y, z], [σ, ρ, β], t, +Joop, +Jiip = eval.(ModelingToolkit.build_function(s∂, [x, y, z], [σ, ρ, β], t, linenumbers = true)) J = Joop([1.0, 2.0, 3.0], [1.0, 2.0, 3.0], 1.0) @test length(nonzeros(s∂)) == 8 @@ -147,7 +148,8 @@ function test_worldage() eqs = [σ * (y - x), x * (ρ - z) - y, x * y - β * z] - f, f_iip = ModelingToolkit.build_function(eqs, [x, y, z], [σ, ρ, β]; + f, + f_iip = ModelingToolkit.build_function(eqs, [x, y, z], [σ, ρ, β]; expression = Val{false}) out = [1.0, 2, 3] o1 = f([1.0, 2, 3], [1.0, 2, 3]) diff --git a/test/downstream/analysis_points.jl b/test/downstream/analysis_points.jl index 21c33b2068..a9ccf0ac2d 100644 --- a/test/downstream/analysis_points.jl +++ b/test/downstream/analysis_points.jl @@ -290,7 +290,8 @@ end connect(F.output, sys_inner.add.input1)] sys_outer = System(eqs, t, systems = [F, sys_inner, r], name = :outer) - matrices, _ = get_sensitivity( + matrices, + _ = get_sensitivity( sys_outer, [sys_outer.inner.plant_input, sys_outer.inner.plant_output]) Ps = tf(1, [1, 1]) |> ss @@ -305,7 +306,8 @@ end @test tf(G[1, 2]) ≈ tf(-CS.feedback(Cs, Ps)) @test tf(G[2, 1]) ≈ tf(CS.feedback(Ps, Cs)) - matrices, _ = get_comp_sensitivity( + matrices, + _ = get_comp_sensitivity( sys_outer, [sys_outer.inner.plant_input, sys_outer.inner.plant_output]) G = CS.ss(matrices...) |> sminreal @@ -327,7 +329,8 @@ end @test tf(L[1, 1]) ≈ -tf(Ps * Cs) # Calling looptransfer like below is not the intended way, but we can work out what it should return if we did so it remains a valid test - matrices, _ = get_looptransfer( + matrices, + _ = get_looptransfer( sys_outer, [sys_outer.inner.plant_input, sys_outer.inner.plant_output]) L = CS.ss(matrices...) |> sminreal @test tf(L[1, 1]) ≈ tf(0) @@ -335,7 +338,8 @@ end @test sminreal(L[1, 2]) ≈ ss(-1) @test tf(L[2, 1]) ≈ tf(Ps) - matrices, _ = linearize( + matrices, + _ = linearize( sys_outer, [sys_outer.inner.plant_input], [nameof(sys_inner.plant_output)]) G = CS.ss(matrices...) |> sminreal @test tf(G) ≈ tf(CS.feedback(Ps, Cs)) diff --git a/test/downstream/linearization_dd.jl b/test/downstream/linearization_dd.jl index 3e94a02469..44659f03ff 100644 --- a/test/downstream/linearization_dd.jl +++ b/test/downstream/linearization_dd.jl @@ -43,9 +43,11 @@ ps = poles(G) @test minimum(abs, ps) < 1e-6 @test minimum(abs, complex(0, 1.3777260367206716) .- ps) < 1e-10 -lsys, syss = linearize(model, lin_inputs, lin_outputs, allow_symbolic = true, op = op, +lsys, +syss = linearize(model, lin_inputs, lin_outputs, allow_symbolic = true, op = op, allow_input_derivatives = true, zero_dummy_der = true) -lsyss, sysss = ModelingToolkit.linearize_symbolic(model, lin_inputs, lin_outputs; +lsyss, +sysss = ModelingToolkit.linearize_symbolic(model, lin_inputs, lin_outputs; allow_input_derivatives = true) dummyder = setdiff(unknowns(sysss), unknowns(model)) diff --git a/test/downstream/test_disturbance_model.jl b/test/downstream/test_disturbance_model.jl index c7b9d833d0..642ff85f99 100644 --- a/test/downstream/test_disturbance_model.jl +++ b/test/downstream/test_disturbance_model.jl @@ -149,10 +149,14 @@ sol = solve(prob, Tsit5()) ## Generate function for an augmented Unscented Kalman Filter ===================== # temp = open_loop(model_with_disturbance, :d) outputs = [P.inertia1.phi, P.inertia2.phi, P.inertia1.w, P.inertia2.w] -f, x_sym, p_sym, io_sys = ModelingToolkit.generate_control_function( +f, x_sym, +p_sym, +io_sys = ModelingToolkit.generate_control_function( model_with_disturbance, [:u], [:d1, :d2, :dy], split = false) -f, x_sym, p_sym, io_sys = ModelingToolkit.generate_control_function( +f, x_sym, +p_sym, +io_sys = ModelingToolkit.generate_control_function( model_with_disturbance, [:u], [:d1, :d2, :dy], disturbance_argument = true, split = false) diff --git a/test/dq_units.jl b/test/dq_units.jl index 615de1d642..ea1103db57 100644 --- a/test/dq_units.jl +++ b/test/dq_units.jl @@ -39,18 +39,18 @@ System(eqs, t, name = :sys, checks = false) # connection validation @connector function Pin(; name) - sts = @variables(v(t)=1.0, [unit = u"V"], - i(t)=1.0, [unit = u"A", connect = Flow]) + sts = @variables(v(t)=1.0, [unit=u"V"], + i(t)=1.0, [unit=u"A", connect=Flow]) System(Equation[], t, sts, []; name = name) end @connector function OtherPin(; name) - sts = @variables(v(t)=1.0, [unit = u"mV"], - i(t)=1.0, [unit = u"mA", connect = Flow]) + sts = @variables(v(t)=1.0, [unit=u"mV"], + i(t)=1.0, [unit=u"mA", connect=Flow]) System(Equation[], t, sts, []; name = name) end @connector function LongPin(; name) - sts = @variables(v(t)=1.0, [unit = u"V"], - i(t)=1.0, [unit = u"A", connect = Flow], + sts = @variables(v(t)=1.0, [unit=u"V"], + i(t)=1.0, [unit=u"A", connect=Flow], x(t)=1.0) System(Equation[], t, sts, []; name = name) end diff --git a/test/extensions/dynamic_optimization.jl b/test/extensions/dynamic_optimization.jl index 3884bc1aa9..f9c6a4277f 100644 --- a/test/extensions/dynamic_optimization.jl +++ b/test/extensions/dynamic_optimization.jl @@ -47,7 +47,7 @@ const M = ModelingToolkit @test ≈(csol2.sol.u, osol2.u, rtol = 0.001) pprob = PyomoDynamicOptProblem(sys, [u0map; parammap], tspan, dt = 0.01) psol = solve(pprob, PyomoCollocation("ipopt", BackwardEuler())) - @test all([≈(psol.sol(t), osol2(t), rtol = 1e-2) for t in 0.:0.01:1.]) + @test all([≈(psol.sol(t), osol2(t), rtol = 1e-2) for t in 0.0:0.01:1.0]) # With a constraint u0map = Pair[] @@ -55,7 +55,8 @@ const M = ModelingToolkit constr = [x(0.6) ~ 3.5, x(0.3) ~ 7.0] @mtkcompile lksys = System(eqs, t; constraints = constr) - jprob = JuMPDynamicOptProblem(lksys, [u0map; parammap], tspan; guesses = guess, dt = 0.01) + jprob = JuMPDynamicOptProblem( + lksys, [u0map; parammap], tspan; guesses = guess, dt = 0.01) jsol = solve(jprob, JuMPCollocation(Ipopt.Optimizer, constructTsitouras5())) @test jsol.sol(0.6; idxs = x(t)) ≈ 3.5 @test jsol.sol(0.3; idxs = x(t)) ≈ 7.0 @@ -74,7 +75,8 @@ const M = ModelingToolkit iprob = InfiniteOptDynamicOptProblem( lksys, [u0map; parammap], tspan; guesses = guess, dt = 0.01) - isol = solve(iprob, InfiniteOptCollocation(Ipopt.Optimizer, InfiniteOpt.OrthogonalCollocation(3))) # 48.564 ms, 9.58 MiB + isol = solve(iprob, + InfiniteOptCollocation(Ipopt.Optimizer, InfiniteOpt.OrthogonalCollocation(3))) # 48.564 ms, 9.58 MiB sol = isol.sol @test sol(0.6; idxs = x(t)) ≈ 3.5 @test sol(0.3; idxs = x(t)) ≈ 7.0 @@ -84,14 +86,17 @@ const M = ModelingToolkit @mtkcompile lksys = System(eqs, t; constraints = constr) iprob = InfiniteOptDynamicOptProblem( lksys, [u0map; parammap], tspan; guesses = guess, dt = 0.01) - isol = solve(iprob, InfiniteOptCollocation(Ipopt.Optimizer, InfiniteOpt.OrthogonalCollocation(3))) + isol = solve(iprob, + InfiniteOptCollocation(Ipopt.Optimizer, InfiniteOpt.OrthogonalCollocation(3))) @test all(u -> u > [1, 1], isol.sol.u) - jprob = JuMPDynamicOptProblem(lksys, [u0map; parammap], tspan; guesses = guess, dt = 0.01) + jprob = JuMPDynamicOptProblem( + lksys, [u0map; parammap], tspan; guesses = guess, dt = 0.01) jsol = solve(jprob, JuMPCollocation(Ipopt.Optimizer, constructRadauIA3())) @test all(u -> u > [1, 1], jsol.sol.u) - pprob = PyomoDynamicOptProblem(lksys, [u0map; parammap], tspan; guesses = guess, dt = 0.01) + pprob = PyomoDynamicOptProblem( + lksys, [u0map; parammap], tspan; guesses = guess, dt = 0.01) psol = solve(pprob, PyomoCollocation("ipopt", MidpointEuler())) @test all(u -> u > [1, 1], psol.sol.u) @@ -110,9 +115,9 @@ function is_bangbang(input_sol, lbounds, ubounds, rtol = 1e-4) end function ctrl_to_spline(inputsol, splineType) - us = reduce(vcat, inputsol.u) - ts = reduce(vcat, inputsol.t) - splineType(us, ts) + us = reduce(vcat, inputsol.u) + ts = reduce(vcat, inputsol.t) + splineType(us, ts) end @testset "Linear systems" begin @@ -131,7 +136,8 @@ end tspan = (0.0, 1.0) parammap = [u(t) => 0.0] jprob = JuMPDynamicOptProblem(block, [u0map; parammap], tspan; dt = 0.01) - jsol = solve(jprob, JuMPCollocation(Ipopt.Optimizer, constructVerner8()), verbose = true) + jsol = solve( + jprob, JuMPCollocation(Ipopt.Optimizer, constructVerner8()), verbose = true) # Linear systems have bang-bang controls @test is_bangbang(jsol.input_sol, [-1.0], [1.0]) # Test reached final position. @@ -156,7 +162,7 @@ end isol = solve(iprob, InfiniteOptCollocation(Ipopt.Optimizer)) @test is_bangbang(isol.input_sol, [-1.0], [1.0]) @test ≈(isol.sol[x(t)][end], 0.25, rtol = 1e-5) - + pprob = PyomoDynamicOptProblem(block, [u0map; parammap], tspan; dt = 0.01) psol = solve(pprob, PyomoCollocation("ipopt", BackwardEuler())) @test is_bangbang(psol.input_sol, [-1.0], [1.0]) @@ -165,7 +171,7 @@ end spline = ctrl_to_spline(isol.input_sol, ConstantInterpolation) oprob = ODEProblem(block_ode, [u0map; u_interp => spline], tspan) @test ≈(isol.sol.u, osol.u, rtol = 0.05) - @test all([≈(psol.sol(t), osol(t), rtol = 0.05) for t in 0.:0.01:1.]) + @test all([≈(psol.sol(t), osol(t), rtol = 0.05) for t in 0.0:0.01:1.0]) ################### ### Bee example ### @@ -210,7 +216,7 @@ end @test ≈(osol.u, csol.sol.u, rtol = 0.01) osol2 = solve(oprob, ImplicitEuler(); dt = 0.01, adaptive = false) @test ≈(osol2.u, isol.sol.u, rtol = 0.01) - @test all([≈(psol.sol(t), osol2(t), rtol = 0.01) for t in 0.:0.01:4.]) + @test all([≈(psol.sol(t), osol2(t), rtol = 0.01) for t in 0.0:0.01:4.0]) end @testset "Rocket launch" begin @@ -240,7 +246,8 @@ end jsol = solve(jprob, JuMPCollocation(Ipopt.Optimizer, constructRadauIIA5())) @test jsol.sol[h(t)][end] > 1.012 - cprob = CasADiDynamicOptProblem(rocket, [u0map; pmap], (ts, te); dt = 0.001, cse = false) + cprob = CasADiDynamicOptProblem( + rocket, [u0map; pmap], (ts, te); dt = 0.001, cse = false) csol = solve(cprob, CasADiCollocation("ipopt")) @test csol.sol[h(t)][end] > 1.012 @@ -378,7 +385,8 @@ end @test csol.sol.u[end] ≈ [π, 0, 0, 0] iprob = InfiniteOptDynamicOptProblem(cartpole, [u0map; pmap], tspan; dt = 0.04) - isol = solve(iprob, InfiniteOptCollocation(Ipopt.Optimizer, InfiniteOpt.OrthogonalCollocation(2))) + isol = solve(iprob, + InfiniteOptCollocation(Ipopt.Optimizer, InfiniteOpt.OrthogonalCollocation(2))) @test isol.sol.u[end] ≈ [π, 0, 0, 0] pprob = PyomoDynamicOptProblem(cartpole, [u0map; pmap], tspan; dt = 0.04) diff --git a/test/extensions/homotopy_continuation.jl b/test/extensions/homotopy_continuation.jl index a4af57c8fe..20ff262132 100644 --- a/test/extensions/homotopy_continuation.jl +++ b/test/extensions/homotopy_continuation.jl @@ -192,7 +192,7 @@ end disallowed_y = [7, 5, 4] @test all(!isapprox(sol[x]; atol = 1e-8), disallowed_x) @test all(!isapprox(sol[y]; atol = 1e-8), disallowed_y) - @test abs(sol[x^2 - 4x + y]) >= 1e-8 + @test abs(sol[x ^ 2 - 4x + y]) >= 1e-8 p = parameter_values(prob) for val in disallowed_x diff --git a/test/extensions/test_infiniteopt.jl b/test/extensions/test_infiniteopt.jl index e1e4143bb7..fab44d0cf9 100644 --- a/test/extensions/test_infiniteopt.jl +++ b/test/extensions/test_infiniteopt.jl @@ -53,13 +53,13 @@ m = InfiniteModel(optimizer_with_attributes(Ipopt.Optimizer, guess_xs = [t -> pi, t -> 0.1][permutation] guess_us = [t -> 0.1] InfiniteOpt.@variables(m, - begin - # state variables - (lb[i] <= x[i = 1:nx] <= ub[i], Infinite(τ), start = guess_xs[i]) # state variables - -10 <= u[i = 1:nu] <= 10, Infinite(τ), (start = guess_us[i]) # control variables - 0 <= tf <= 10, (start = 5) # Final time - 0.2 <= L <= 0.6, (start = 0.4) # Length parameter - end) +begin + # state variables + (lb[i] <= x[i = 1:nx] <= ub[i], Infinite(τ), start = guess_xs[i]) # state variables + -10 <= u[i = 1:nu] <= 10, Infinite(τ), (start = guess_us[i]) # control variables + 0 <= tf <= 10, (start = 5) # Final time + 0.2 <= L <= 0.6, (start = 0.4) # Length parameter +end) # Trace the dynamics x0 = ModelingToolkit.get_u0(io_sys, [model.θ => 0, model.ω => 0]) diff --git a/test/fmi/fmi.jl b/test/fmi/fmi.jl index 85e1830650..de5b8a1dac 100644 --- a/test/fmi/fmi.jl +++ b/test/fmi/fmi.jl @@ -177,7 +177,7 @@ end end function build_sspace_model(sspace) - @variables u(t)=1.0 x(t)=1.0 y(t) [guess = 1.0] + @variables u(t)=1.0 x(t)=1.0 y(t) [guess=1.0] @mtkcompile sys = System( [sspace.u ~ u, D(u) ~ t, D(x) ~ sspace.x + sspace.y, y^2 ~ sspace.y + sspace.x], t; systems = [sspace] diff --git a/test/initial_values.jl b/test/initial_values.jl index d4d768c0fe..5412d4d58a 100644 --- a/test/initial_values.jl +++ b/test/initial_values.jl @@ -7,7 +7,7 @@ using SymbolicIndexingInterface @variables x(t)[1:3]=[1.0, 2.0, 3.0] y(t) z(t)[1:2] -@mtkcompile sys=System([D(x) ~ t * x], t) simplify=false +@mtkcompile sys=System([D(x)~t*x], t) simplify=false reorderer = getsym(sys, x) @test reorderer(get_u0(sys, [])) == [1.0, 2.0, 3.0] @test reorderer(get_u0(sys, [x => [2.0, 3.0, 4.0]])) == [2.0, 3.0, 4.0] @@ -15,10 +15,10 @@ reorderer = getsym(sys, x) @test get_u0(sys, [2.0, 3.0, 4.0]) == [2.0, 3.0, 4.0] @mtkcompile sys=System([ - D(x) ~ 3x, - D(y) ~ t, - D(z[1]) ~ z[2] + t, - D(z[2]) ~ y + z[1] + D(x)~3x, + D(y)~t, + D(z[1])~z[2]+t, + D(z[2])~y+z[1] ], t) simplify=false @test_throws ModelingToolkit.MissingVariablesError get_u0(sys, []) @@ -119,7 +119,7 @@ end @testset "split=false systems with all parameter defaults" begin @variables x(t) = 1.0 @parameters p=1.0 q=2.0 r=3.0 - @mtkcompile sys=System(D(x) ~ p * x + q * t + r, t) split=false + @mtkcompile sys=System(D(x)~p*x+q*t+r, t) split=false prob = @test_nowarn ODEProblem(sys, [], (0.0, 1.0)) @test prob.p isa Vector{Float64} end @@ -254,7 +254,7 @@ end @testset "Array initials and scalar parameters with `split = false`" begin @variables x(t)[1:2] @parameters p - @mtkcompile sys=System([D(x[1]) ~ x[1], D(x[2]) ~ x[2] + p], t) split=false + @mtkcompile sys=System([D(x[1])~x[1], D(x[2])~x[2]+p], t) split=false ps = Set(parameters(sys; initial_parameters = true)) @test length(ps) == 5 for i in 1:2 diff --git a/test/initializationsystem.jl b/test/initializationsystem.jl index 99e73e5b17..0f7f762911 100644 --- a/test/initializationsystem.jl +++ b/test/initializationsystem.jl @@ -600,7 +600,8 @@ end # the correct solver. # `rhss` allows adding terms to the end of equations (only 2 equations allowed) to influence # the system type (brownian vars to turn it into an SDE). - @testset "$Problem with $(SciMLBase.parameterless_type(alg)) and $ctor ctor" for ((Problem, alg, rhss), (ctor, expectedT)) in Iterators.product( + @testset "$Problem with $(SciMLBase.parameterless_type(alg)) and $ctor ctor" for ( + (Problem, alg, rhss), (ctor, expectedT)) in Iterators.product( [ (ODEProblem, Tsit5(), zeros(2)), (SDEProblem, ImplicitEM(), [a, b]), @@ -778,7 +779,7 @@ end @testset "No initialization for variables" begin @variables x=1.0 y=0.0 z=0.0 - @parameters σ=10.0 ρ=26.0 β=8 / 3 + @parameters σ=10.0 ρ=26.0 β=8/3 eqs = [0 ~ σ * (y - x), 0 ~ x * (ρ - z) - y, @@ -823,7 +824,7 @@ end @test SciMLBase.successful_retcode(sol) end - @parameters p=2.0 q=missing [guess = 1.0] c=1.0 + @parameters p=2.0 q=missing [guess=1.0] c=1.0 @variables x=1.0 z=3.0 # eqs = [0 ~ p * (y - x), @@ -884,7 +885,8 @@ end @brownian a b x = _x(t) - @testset "$Problem with $(SciMLBase.parameterless_type(typeof(alg)))" for (System, Problem, alg, rhss) in [ + @testset "$Problem with $(SciMLBase.parameterless_type(typeof(alg)))" for ( + System, Problem, alg, rhss) in [ (ModelingToolkit.System, ODEProblem, Tsit5(), zeros(2)), (ModelingToolkit.System, SDEProblem, ImplicitEM(), [a, b]), (ModelingToolkit.System, DDEProblem, MethodOfSteps(Tsit5()), [_x(t - 0.1), 0.0]), @@ -913,7 +915,8 @@ end @brownian a x = _x(t) - @testset "$Problem with $(SciMLBase.parameterless_type(typeof(alg)))" for (System, Problem, alg, rhss) in [ + @testset "$Problem with $(SciMLBase.parameterless_type(typeof(alg)))" for ( + System, Problem, alg, rhss) in [ (ModelingToolkit.System, ODEProblem, Tsit5(), 0), (ModelingToolkit.System, SDEProblem, ImplicitEM(), a), (ModelingToolkit.System, DDEProblem, MethodOfSteps(Tsit5()), _x(t - 0.1)), @@ -936,7 +939,8 @@ end @brownian a b x = _x(t) - @testset "$Problem with $(SciMLBase.parameterless_type(typeof(alg)))" for (Problem, alg, rhss) in [ + @testset "$Problem with $(SciMLBase.parameterless_type(typeof(alg)))" for ( + Problem, alg, rhss) in [ (ODEProblem, Tsit5(), zeros(2)), (SDEProblem, ImplicitEM(), [a, b]), (DDEProblem, MethodOfSteps(Tsit5()), [_x(t - 0.1), 0.0]), @@ -968,7 +972,8 @@ end @brownian a x = _x(t) - @testset "$Problem with $(SciMLBase.parameterless_type(typeof(alg)))" for (System, Problem, alg, rhss) in [ + @testset "$Problem with $(SciMLBase.parameterless_type(typeof(alg)))" for ( + System, Problem, alg, rhss) in [ (ModelingToolkit.System, ODEProblem, Tsit5(), 0), (ModelingToolkit.System, SDEProblem, ImplicitEM(), a), (ModelingToolkit.System, DDEProblem, MethodOfSteps(Tsit5()), _x(t - 0.1)), @@ -1025,7 +1030,8 @@ end @brownian a x = _x(t) - @testset "$Problem with $(SciMLBase.parameterless_type(typeof(alg)))" for (System, Problem, alg, rhss) in [ + @testset "$Problem with $(SciMLBase.parameterless_type(typeof(alg)))" for ( + System, Problem, alg, rhss) in [ (ModelingToolkit.System, ODEProblem, Tsit5(), 0), (ModelingToolkit.System, SDEProblem, ImplicitEM(), a), (ModelingToolkit.System, DDEProblem, MethodOfSteps(Tsit5()), _x(t - 0.1)), @@ -1191,7 +1197,7 @@ end @testset "DAEProblem initialization" begin @variables x(t) [guess = 1.0] y(t) [guess = 1.0] - @parameters p=missing [guess = 1.0] q=missing [guess = 1.0] + @parameters p=missing [guess=1.0] q=missing [guess=1.0] @mtkcompile sys = System( [D(x) ~ p * y + q, x^3 + y^3 ~ 5], t; initialization_eqs = [p^2 + q^3 ~ 3]) @@ -1608,7 +1614,8 @@ end tspan = (0.0, 100.0) getter = getsym(sys, Initial.(unknowns(sys))) prob = ODEProblem(sys, [u0; p], tspan; guesses = [w2 => 3.0]) - new_u0, new_p, _ = SciMLBase.get_initial_values( + new_u0, new_p, + _ = SciMLBase.get_initial_values( prob, prob, prob.f, SciMLBase.OverrideInit(), Val(true); nlsolve_alg = NewtonRaphson(), abstol = 1e-6, reltol = 1e-3) @test getter(prob) != getter(new_p) diff --git a/test/input_output_handling.jl b/test/input_output_handling.jl index ed1fb5fc69..8b966ca814 100644 --- a/test/input_output_handling.jl +++ b/test/input_output_handling.jl @@ -159,13 +159,14 @@ end for split in [true, false] simplify = true - @variables x(t)=0 u(t)=0 [input = true] + @variables x(t)=0 u(t)=0 [input=true] eqs = [ D(x) ~ -x + u ] @named sys = System(eqs, t) - f, dvs, ps, io_sys = ModelingToolkit.generate_control_function( + f, dvs, + ps, io_sys = ModelingToolkit.generate_control_function( sys, [u]; simplify, split) @test isequal(dvs[], x) @@ -177,13 +178,15 @@ end @test f[1](x, u, p, 1) ≈ -x + u # With disturbance inputs - @variables x(t)=0 u(t)=0 [input = true] d(t)=0 + @variables x(t)=0 u(t)=0 [input=true] d(t)=0 eqs = [ D(x) ~ -x + u + d^2 ] @named sys = System(eqs, t) - f, dvs, ps, io_sys = ModelingToolkit.generate_control_function( + f, dvs, + ps, + io_sys = ModelingToolkit.generate_control_function( sys, [u], [d]; simplify, split) @test isequal(dvs[], x) @@ -195,13 +198,15 @@ end @test f[1](x, u, p, 1) ≈ -x + u ## With added d argument - @variables x(t)=0 u(t)=0 [input = true] d(t)=0 + @variables x(t)=0 u(t)=0 [input=true] d(t)=0 eqs = [ D(x) ~ -x + u + d^2 ] @named sys = System(eqs, t) - f, dvs, ps, io_sys = ModelingToolkit.generate_control_function( + f, dvs, + ps, + io_sys = ModelingToolkit.generate_control_function( sys, [u], [d]; simplify, split, disturbance_argument = true) @@ -409,7 +414,8 @@ f, augmented_sys, dvs, p = ModelingToolkit.add_input_disturbance(model, ins) augmented_sys = complete(augmented_sys) -matrices, ssys = linearize(augmented_sys, +matrices, +ssys = linearize(augmented_sys, [ augmented_sys.u, augmented_sys.input.u[2], @@ -429,7 +435,7 @@ matrices = ModelingToolkit.reorder_unknowns( # @test sminreal(G[1, 3]) ≈ sminreal(P[1,1])*dist @testset "Observed functions with inputs" begin - @variables x(t)=0 u(t)=0 [input = true] + @variables x(t)=0 u(t)=0 [input=true] eqs = [ D(x) ~ -x + u ] @@ -454,7 +460,7 @@ end end @testset "With callable symbolic" begin - @variables x(t)=0 u(t)=0 [input = true] + @variables x(t)=0 u(t)=0 [input=true] @parameters p(::Real) = (x -> 2x) eqs = [D(x) ~ -x + p(u)] @named sys = System(eqs, t) diff --git a/test/inputoutput.jl b/test/inputoutput.jl index 28f112425f..5a503be9e4 100644 --- a/test/inputoutput.jl +++ b/test/inputoutput.jl @@ -29,11 +29,11 @@ sys = connected collapsed_eqs = [ D(lorenz1.x) ~ (lorenz1.σ * (lorenz1.y - lorenz1.x) + - (lorenz2.x + lorenz2.y - lorenz2.z)), + (lorenz2.x + lorenz2.y - lorenz2.z)), D(lorenz1.y) ~ lorenz1.x * (lorenz1.ρ - lorenz1.z) - lorenz1.y, D(lorenz1.z) ~ lorenz1.x * lorenz1.y - (lorenz1.β * lorenz1.z), D(lorenz2.x) ~ (lorenz2.σ * (lorenz2.y - lorenz2.x) + - (lorenz1.x + lorenz1.y - lorenz1.z)), + (lorenz1.x + lorenz1.y - lorenz1.z)), D(lorenz2.y) ~ lorenz2.x * (lorenz2.ρ - lorenz2.z) - lorenz2.y, D(lorenz2.z) ~ lorenz2.x * lorenz2.y - (lorenz2.β * lorenz2.z)] diff --git a/test/jacobiansparsity.jl b/test/jacobiansparsity.jl index 7173ab9e1a..4c8c47feab 100644 --- a/test/jacobiansparsity.jl +++ b/test/jacobiansparsity.jl @@ -12,13 +12,17 @@ function brusselator_2d_loop(du, u, p, t) x, y = xyd_brusselator[I[1]], xyd_brusselator[I[2]] ip1, im1, jp1, jm1 = lim(i + 1, N), lim(i - 1, N), lim(j + 1, N), lim(j - 1, N) - du[i, j, 1] = alpha * (u[im1, j, 1] + u[ip1, j, 1] + u[i, jp1, 1] + u[i, jm1, 1] - - 4u[i, j, 1]) + - B + u[i, j, 1]^2 * u[i, j, 2] - (A + 1) * u[i, j, 1] + - brusselator_f(x, y, t) - du[i, j, 2] = alpha * (u[im1, j, 2] + u[ip1, j, 2] + u[i, jp1, 2] + u[i, jm1, 2] - - 4u[i, j, 2]) + - A * u[i, j, 1] - u[i, j, 1]^2 * u[i, j, 2] + du[i, + j, + 1] = alpha * (u[im1, j, 1] + u[ip1, j, 1] + u[i, jp1, 1] + u[i, jm1, 1] - + 4u[i, j, 1]) + + B + u[i, j, 1]^2 * u[i, j, 2] - (A + 1) * u[i, j, 1] + + brusselator_f(x, y, t) + du[i, + j, + 2] = alpha * (u[im1, j, 2] + u[ip1, j, 2] + u[i, jp1, 2] + u[i, jm1, 2] - + 4u[i, j, 2]) + + A * u[i, j, 1] - u[i, j, 1]^2 * u[i, j, 2] end end diff --git a/test/linearize.jl b/test/linearize.jl index ddb90f3eac..ce8da1708f 100644 --- a/test/linearize.jl +++ b/test/linearize.jl @@ -4,7 +4,7 @@ using CommonSolve: solve # r is an input, and y is an output. @independent_variables t @variables x(t)=0 y(t)=0 u(t)=0 r(t)=0 -@variables x(t)=0 y(t)=0 u(t)=0 r(t)=0 [input = true] +@variables x(t)=0 y(t)=0 u(t)=0 r(t)=0 [input=true] @parameters kp = 1 D = Differential(t) @@ -115,7 +115,8 @@ Nd = 10 @named pid = LimPID(; k, Ti, Td, Nd) @unpack reference, measurement, ctr_output = pid -lsys0, ssys = linearize(pid, [reference.u, measurement.u], [ctr_output.u]; +lsys0, +ssys = linearize(pid, [reference.u, measurement.u], [ctr_output.u]; op = Dict(reference.u => 0.0, measurement.u => 0.0)) @unpack int, der = pid desired_order = [int.x, der.x] @@ -126,7 +127,8 @@ lsys = ModelingToolkit.reorder_unknowns(lsys0, unknowns(ssys), desired_order) @test lsys.C == [400 -4000] @test lsys.D == [4400 -4400] -lsyss0, ssys2 = ModelingToolkit.linearize_symbolic(pid, [reference.u, measurement.u], +lsyss0, +ssys2 = ModelingToolkit.linearize_symbolic(pid, [reference.u, measurement.u], [ctr_output.u]) lsyss = ModelingToolkit.reorder_unknowns(lsyss0, unknowns(ssys2), desired_order) diff --git a/test/modelingtoolkitize.jl b/test/modelingtoolkitize.jl index 3f2db92344..30cde7c3b7 100644 --- a/test/modelingtoolkitize.jl +++ b/test/modelingtoolkitize.jl @@ -14,15 +14,20 @@ function brusselator_2d_loop(du, u, p, t) @inbounds for I in CartesianIndices((N, N)) i, j = Tuple(I) x, y = xyd_brusselator[I[1]], xyd_brusselator[I[2]] - ip1, im1, jp1, jm1 = limit(i + 1, N), limit(i - 1, N), limit(j + 1, N), + ip1, im1, jp1, + jm1 = limit(i + 1, N), limit(i - 1, N), limit(j + 1, N), limit(j - 1, N) - du[i, j, 1] = alpha * (u[im1, j, 1] + u[ip1, j, 1] + u[i, jp1, 1] + u[i, jm1, 1] - - 4u[i, j, 1]) + - B + u[i, j, 1]^2 * u[i, j, 2] - (A + 1) * u[i, j, 1] + - brusselator_f(x, y, t) - du[i, j, 2] = alpha * (u[im1, j, 2] + u[ip1, j, 2] + u[i, jp1, 2] + u[i, jm1, 2] - - 4u[i, j, 2]) + - A * u[i, j, 1] - u[i, j, 1]^2 * u[i, j, 2] + du[i, + j, + 1] = alpha * (u[im1, j, 1] + u[ip1, j, 1] + u[i, jp1, 1] + u[i, jm1, 1] - + 4u[i, j, 1]) + + B + u[i, j, 1]^2 * u[i, j, 2] - (A + 1) * u[i, j, 1] + + brusselator_f(x, y, t) + du[i, + j, + 2] = alpha * (u[im1, j, 2] + u[ip1, j, 2] + u[i, jp1, 2] + u[i, jm1, 2] - + 4u[i, j, 2]) + + A * u[i, j, 1] - u[i, j, 1]^2 * u[i, j, 2] end end @@ -115,7 +120,9 @@ function SIRD_ac!(du, u, p, t) # initialize this parameter (death probability stratified by age, taken from literature) - δ₁, δ₂, δ₃, δ₄ = [ + δ₁, δ₂, + δ₃, + δ₄ = [ 0.003 / 100, 0.004 / 100, (0.015 + 0.030 + 0.064 + 0.213 + 0.718) / (5 * 100), @@ -255,7 +262,7 @@ sys = modelingtoolkitize(prob) @test [ModelingToolkit.defaults(sys)[s] for s in unknowns(sys)] == u0 @test [ModelingToolkit.defaults(sys)[s] for s in parameters(sys)] == [10, 20] -@parameters sig=10 rho=28.0 beta=8 / 3 +@parameters sig=10 rho=28.0 beta=8/3 @variables x(t)=100 y(t)=1.0 z(t)=1 eqs = [D(x) ~ sig * (y - x), diff --git a/test/mtkparameters.jl b/test/mtkparameters.jl index 1508d95a74..28ab3759ef 100644 --- a/test/mtkparameters.jl +++ b/test/mtkparameters.jl @@ -56,8 +56,8 @@ setp(sys, a)(ps, 1.0) @test getp(sys, a)(ps) == getp(sys, b)(ps) / 2 == getp(sys, c)(ps) / 3 == 1.0 for (portion, values) in [(Tunable(), [1.0, 5.0, 6.0, 7.0]) - (Discrete(), [3.0]) - (Constants(), vcat([0.1, 0.2, 0.3], ones(9), [4.0]))] + (Discrete(), [3.0]) + (Constants(), vcat([0.1, 0.2, 0.3], ones(9), [4.0]))] buffer, repack, alias = canonicalize(portion, ps) @test alias @test sort(collect(buffer)) == values @@ -173,7 +173,7 @@ ps = [p => 1.0] # Value for `d` is missing # scalar parameters only function level1() - @parameters p1=0.5 [tunable = true] p2=1 [tunable = true] p3=3 [tunable = false] p4=3 [tunable = true] y0=1 + @parameters p1=0.5 [tunable=true] p2=1 [tunable=true] p3=3 [tunable=false] p4=3 [tunable=true] y0=1 @variables x(t)=2 y(t)=y0 D = Differential(t) @@ -188,7 +188,7 @@ end # scalar and vector parameters function level2() - @parameters p1=0.5 [tunable = true] (p23[1:2]=[1, 3.0]) [tunable = true] p4=3 [tunable = false] y0=1 + @parameters p1=0.5 [tunable=true] (p23[1:2]=[1, 3.0]) [tunable=true] p4=3 [tunable=false] y0=1 @variables x(t)=2 y(t)=y0 D = Differential(t) @@ -203,7 +203,7 @@ end # scalar and vector parameters with different scalar types function level3() - @parameters p1=0.5 [tunable = true] (p23[1:2]=[1, 3.0]) [tunable = true] p4::Int=3 [tunable = true] y0::Int=1 + @parameters p1=0.5 [tunable=true] (p23[1:2]=[1, 3.0]) [tunable=true] p4::Int=3 [tunable=true] y0::Int=1 @variables x(t)=2 y(t)=y0 D = Differential(t) diff --git a/test/nonlinearsystem.jl b/test/nonlinearsystem.jl index 73cbbe9639..4c887f5740 100644 --- a/test/nonlinearsystem.jl +++ b/test/nonlinearsystem.jl @@ -408,7 +408,7 @@ end @test prob.f.initialization_data.initializeprobmap === nothing sol = solve(prob) @test SciMLBase.successful_retcode(sol) - @test sol.ps[p^3 + q^3]≈sol.ps[4r] atol=1e-10 + @test sol.ps[p ^ 3 + q ^ 3]≈sol.ps[4r] atol=1e-10 @testset "Differential inside expression also substituted" begin @named sys = System([0 ~ y * D(x) + x^2 - p, 0 ~ x * D(y) + y * p], t) diff --git a/test/odesystem.jl b/test/odesystem.jl index 28bc1d7db0..c5733866c3 100644 --- a/test/odesystem.jl +++ b/test/odesystem.jl @@ -816,7 +816,7 @@ Set(unknowns(new_sys2)) == Set([new_sys2.x1, new_sys2.sys1.x1, new_sys2.sub.s1, new_sys2.sub.s2]) let # Issue https://github.com/SciML/ModelingToolkit.jl/issues/2322 - @parameters a=10 b=a / 10 c=a / 20 + @parameters a=10 b=a/10 c=a/20 Dt = D @@ -1383,7 +1383,7 @@ end @mtkcompile sys = System([D(x) ~ x, y^2 ~ x + sum(p)], t) prob = DAEProblem(sys, [D(x) => x, D(y) => D(x) / 2y], (0.0, 1.0)) sol = solve(prob, DFBDF(), abstol = 1e-8, reltol = 1e-8) - @test sol[x]≈sol[y^2 - sum(p)] atol=1e-5 + @test sol[x]≈sol[y ^ 2 - sum(p)] atol=1e-5 end @testset "Symbolic tstops" begin @@ -1483,7 +1483,8 @@ end obsfn_expr = ModelingToolkit.build_explicit_observed_function( sys, 2x + 1, expression = true) @test obsfn_expr isa Expr - obsfn_expr_oop, obsfn_expr_iip = ModelingToolkit.build_explicit_observed_function( + obsfn_expr_oop, + obsfn_expr_iip = ModelingToolkit.build_explicit_observed_function( sys, [x + 1, x + 2, x + t], return_inplace = true, expression = true) @test obsfn_expr_oop isa Expr @test obsfn_expr_iip isa Expr diff --git a/test/reduction.jl b/test/reduction.jl index 97d360fd67..6fee37d806 100644 --- a/test/reduction.jl +++ b/test/reduction.jl @@ -172,7 +172,7 @@ nlprob.f(residual, reducedsol.u, pp) N = 5 @variables xs[1:N] -A = reshape(1:(N^2), N, N) +A = reshape(1:(N ^ 2), N, N) eqs = xs ~ A * xs @named sys′ = System(eqs, [xs], []) sys = mtkcompile(sys′) diff --git a/test/scc_nonlinear_problem.jl b/test/scc_nonlinear_problem.jl index fad30dbcea..70031ff228 100644 --- a/test/scc_nonlinear_problem.jl +++ b/test/scc_nonlinear_problem.jl @@ -185,7 +185,11 @@ import ModelingToolkitStandardLibrary.Hydraulic.IsothermalCompressible as IC @testset "Caching of subexpressions of different types" begin liquid_pressure(rho, rho_0, bulk) = (rho / rho_0 - 1) * bulk gas_pressure(rho, rho_0, p_gas, rho_gas) = rho * ((0 - p_gas) / (rho_0 - rho_gas)) - full_pressure(rho, rho_0, bulk, p_gas, rho_gas) = ifelse( + full_pressure(rho, + rho_0, + bulk, + p_gas, + rho_gas) = ifelse( rho >= rho_0, liquid_pressure(rho, rho_0, bulk), gas_pressure(rho, rho_0, p_gas, rho_gas)) diff --git a/test/sciml_problem_inputs.jl b/test/sciml_problem_inputs.jl index cd081118fb..8c0ea5edc1 100644 --- a/test/sciml_problem_inputs.jl +++ b/test/sciml_problem_inputs.jl @@ -111,6 +111,7 @@ end # Simulates problems for all input types, checking that identical solutions are found. # test failure. for u0 in u0_alts, p in p_alts + oprob = remake(base_oprob; u0, p) @test base_sol == solve(oprob, Tsit5(); saveat = 1.0) eprob = remake(base_eprob; u0, p) @@ -124,6 +125,7 @@ end base_sol = solve(base_nlprob, NewtonRaphson()) # Solves problems for all input types, checking that identical solutions are found. for u0 in u0_alts, p in p_alts + nlprob = remake(base_nlprob; u0, p) @test base_sol == solve(nlprob, NewtonRaphson()) end @@ -140,6 +142,7 @@ end # Simulates problems for all input types, checking that identical solutions are found. # test failure. for u0 in u0_alts, p in p_alts + ssprob = remake(base_ssprob; u0, p) @test base_sol == solve(ssprob, DynamicSS(Tsit5())) eprob = remake(base_eprob; u0, p) @@ -206,6 +209,7 @@ end (idiscsys, ImplicitDiscreteProblem{true, SciMLBase.FullSpecialize}) ] @testset "$(typeof(u0)) - $(typeof(p))" for u0 in u0s, p in ps + if u0 isa Vector{Float64} && ctor <: ImplicitDiscreteProblem u0 = ones(2) end @@ -224,6 +228,7 @@ end (optsys, OptimizationProblem{true}) ] @testset "$(typeof(u0)) - $(typeof(p))" for u0 in u0s, p in ps + @test_warn ["deprecated"] ctor(sys, u0, p) end end diff --git a/test/stream_connectors.jl b/test/stream_connectors.jl index 6e49cd1fb3..a9d824c548 100644 --- a/test/stream_connectors.jl +++ b/test/stream_connectors.jl @@ -265,7 +265,7 @@ sys_exp = expand_connections(compose(sys, [sp1, sp2])) # array var @connector function VecPin(; name) - sts = @variables v(t)[1:2]=[1.0, 0.0] i(t)[1:2]=1.0 [connect = Flow] + sts = @variables v(t)[1:2]=[1.0, 0.0] i(t)[1:2]=1.0 [connect=Flow] System(Equation[], t, [sts...;], []; name = name) end @@ -286,7 +286,7 @@ sys = expand_connections(compose(simple, [vp1, vp2, vp3])) 0 ~ -vp1.i[2] - vp2.i[2] - vp3.i[2]]) @connector function VectorHeatPort(; name, N = 100, T0 = 0.0, Q0 = 0.0) - @variables (T(t))[1:N]=T0 (Q(t))[1:N]=Q0 [connect = Flow] + @variables (T(t))[1:N]=T0 (Q(t))[1:N]=Q0 [connect=Flow] System(Equation[], t, [T; Q], []; name = name) end diff --git a/test/structural_transformation/index_reduction.jl b/test/structural_transformation/index_reduction.jl index 229b1f5736..2118c8441d 100644 --- a/test/structural_transformation/index_reduction.jl +++ b/test/structural_transformation/index_reduction.jl @@ -68,5 +68,5 @@ let sys, [x => 1, y => 0, D(x) => 0.0, g => 1], (0.0, 10.0), guesses = [λ => 0.0]) sol = solve(prob, Rodas5P()) @test SciMLBase.successful_retcode(sol) - @test sol[x^2 + y^2][end] < 1.1 + @test sol[x ^ 2 + y ^ 2][end] < 1.1 end diff --git a/test/symbolic_events.jl b/test/symbolic_events.jl index 856a97484f..150b58f7f3 100644 --- a/test/symbolic_events.jl +++ b/test/symbolic_events.jl @@ -433,7 +433,8 @@ end compose(System(Equation[], t; name), spring, damper) end - connect_sd(sd, m1, m2) = [ + connect_sd( + sd, m1, m2) = [ sd.spring.x ~ m1.pos - m2.pos, sd.damper.vel ~ m1.vel - m2.vel] sd_force(sd) = -sd.spring.k * sd.spring.x - sd.damper.c * sd.damper.vel @named mass1 = Mass(; m = 1) @@ -1085,8 +1086,10 @@ end [x ~ 0], nothing, initialize = [x ~ 1.5], finalize = f) inited = false finaled = false - a = ModelingToolkit.ImperativeAffect(f = (m, o, ctx, int) -> (inited = true; return (;))) - b = ModelingToolkit.ImperativeAffect(f = (m, o, ctx, int) -> (finaled = true; return (;))) + a = ModelingToolkit.ImperativeAffect(f = ( + m, o, ctx, int) -> (inited = true; return (;))) + b = ModelingToolkit.ImperativeAffect(f = ( + m, o, ctx, int) -> (finaled = true; return (;))) cb2 = ModelingToolkit.SymbolicContinuousCallback( [x ~ 0.1], nothing, initialize = a, finalize = b) @mtkcompile sys = System(D(x) ~ -1, t, [x], []; continuous_events = [cb1, cb2]) diff --git a/test/test_variable_metadata.jl b/test/test_variable_metadata.jl index 482ebf927c..5118a16b8c 100644 --- a/test/test_variable_metadata.jl +++ b/test/test_variable_metadata.jl @@ -117,7 +117,7 @@ d = FakeNormal() ## System interface @independent_variables t Dₜ = Differential(t) -@variables x(t)=0 [bounds = (-10, 10)] u(t)=0 [input = true] y(t)=0 [output = true] +@variables x(t)=0 [bounds=(-10, 10)] u(t)=0 [input=true] y(t)=0 [output=true] @parameters T [bounds = (0, Inf)] @parameters k [tunable = true, bounds = (0, Inf)] @parameters k2 [tunable = false] diff --git a/test/units.jl b/test/units.jl index 77f35877e9..a17dd90575 100644 --- a/test/units.jl +++ b/test/units.jl @@ -66,19 +66,19 @@ System(eqs, t, name = :sys, checks = false) # connection validation @connector function Pin(; name) - sts = @variables(v(t)=1.0, [unit = u"V"], - i(t)=1.0, [unit = u"A", connect = Flow]) + sts = @variables(v(t)=1.0, [unit=u"V"], + i(t)=1.0, [unit=u"A", connect=Flow]) System(Equation[], t, sts, []; name = name) end @connector function OtherPin(; name) - sts = @variables(v(t)=1.0, [unit = u"mV"], - i(t)=1.0, [unit = u"mA", connect = Flow]) + sts = @variables(v(t)=1.0, [unit=u"mV"], + i(t)=1.0, [unit=u"mA", connect=Flow]) System(Equation[], t, sts, []; name = name) end @connector function LongPin(; name) - sts = @variables(v(t)=1.0, [unit = u"V"], - i(t)=1.0, [unit = u"A", connect = Flow], - x(t)=1.0, [unit = NoUnits]) + sts = @variables(v(t)=1.0, [unit=u"V"], + i(t)=1.0, [unit=u"A", connect=Flow], + x(t)=1.0, [unit=NoUnits]) System(Equation[], t, sts, []; name = name) end @named p1 = Pin() diff --git a/test/variable_parsing.jl b/test/variable_parsing.jl index 1ea366d045..60b4e24d64 100644 --- a/test/variable_parsing.jl +++ b/test/variable_parsing.jl @@ -88,7 +88,7 @@ end @test getmetadata(x, VariableUnit) == u @test getmetadata(y, VariableDefaultValue) === 2 -@variables x=[1, 2] [connect = Flow, unit = u] y=2 +@variables x=[1, 2] [connect=Flow, unit=u] y=2 @test getmetadata(x, VariableDefaultValue) == [1, 2] @test getmetadata(x, VariableConnectType) == Flow