Merge pull request #2760 from SciML/ChrisRackauckas-patch-6

Update to SymbolicUtils v2

Author: ChrisRackauckas

Project.toml

@@ -107,8 +107,8 @@ SparseArrays = "1"
 SpecialFunctions = "0.7, 0.8, 0.9, 0.10, 1.0, 2"
 StaticArrays = "0.10, 0.11, 0.12, 1.0"
 SymbolicIndexingInterface = "0.3.12"
-SymbolicUtils = "<1.6"
-Symbolics = "5.26"
+SymbolicUtils = "2"
+Symbolics = "5.30.1"
 URIs = "1"
 UnPack = "0.1, 1.0"
 Unitful = "1.1"

src/ModelingToolkit.jl

@@ -9,7 +9,7 @@ using PrecompileTools, Reexport
 end
 
 import SymbolicUtils
-import SymbolicUtils: istree, arguments, operation, similarterm, promote_symtype,
+import SymbolicUtils: iscall, arguments, operation, maketerm, promote_symtype,
                       Symbolic, isadd, ismul, ispow, issym, FnType,
                       @rule, Rewriters, substitute, metadata, BasicSymbolic,
                       Sym, Term
@@ -59,7 +59,7 @@ using Symbolics: _parse_vars, value, @derivatives, get_variables,
                  exprs_occur_in, solve_for, build_expr, unwrap, wrap,
                  VariableSource, getname, variable, Connection, connect,
                  NAMESPACE_SEPARATOR, set_scalar_metadata, setdefaultval,
-                 initial_state, transition, activeState, entry,
+                 initial_state, transition, activeState, entry, hasnode,
                  ticksInState, timeInState, fixpoint_sub, fast_substitute
 import Symbolics: rename, get_variables!, _solve, hessian_sparsity,
                   jacobian_sparsity, isaffine, islinear, _iszero, _isone,

src/clock.jl

@@ -21,7 +21,7 @@ function is_continuous_domain(x)
 end
 
 function get_time_domain(x)
-    if istree(x) && operation(x) isa Operator
+    if iscall(x) && operation(x) isa Operator
         output_timedomain(x)
     else
         getmetadata(x, TimeDomain, nothing)
@@ -130,7 +130,7 @@ is_concrete_time_domain(x) = x isa Union{AbstractClock, Continuous}
     SolverStepClock()
     SolverStepClock(t)
 
-A clock that ticks at each solver step (sometimes referred to as "continuous sample time"). This clock **does generally not have equidistant tick intervals**, instead, the tick interval depends on the adaptive step-size slection of the continuous solver, as well as any continuous event handling. If adaptivity of the solver is turned off and there are no continuous events, the tick interval will be given by the fixed solver time step `dt`. 
+A clock that ticks at each solver step (sometimes referred to as "continuous sample time"). This clock **does generally not have equidistant tick intervals**, instead, the tick interval depends on the adaptive step-size slection of the continuous solver, as well as any continuous event handling. If adaptivity of the solver is turned off and there are no continuous events, the tick interval will be given by the fixed solver time step `dt`.
 
 Due to possibly non-equidistant tick intervals, this clock should typically not be used with discrete-time systems that assume a fixed sample time, such as PID controllers and digital filters.
 """

src/debugging.jl

@@ -28,9 +28,9 @@ end
 
 debug_sub(eq::Equation) = debug_sub(eq.lhs) ~ debug_sub(eq.rhs)
 function debug_sub(ex)
-    istree(ex) || return ex
+    iscall(ex) || return ex
     f = operation(ex)
     args = map(debug_sub, arguments(ex))
     f in LOGGED_FUN ? logged_fun(f, args...) :
-    similarterm(ex, f, args, metadata = metadata(ex))
+    maketerm(typeof(ex), f, args, symtype(t), metadata(ex))
 end

src/discretedomain.jl

@@ -39,7 +39,7 @@ end
 function (D::Shift)(x::Num, allow_zero = false)
     !allow_zero && D.steps == 0 && return x
     vt = value(x)
-    if istree(vt)
+    if iscall(vt)
         op = operation(vt)
         if op isa Sample
             error("Cannot shift a `Sample`. Create a variable to represent the sampled value and shift that instead")

src/inputoutput.jl

@@ -19,8 +19,8 @@ function outputs(sys)
     lhss = [eq.lhs for eq in o]
     unique([filter(isoutput, unknowns(sys))
             filter(isoutput, parameters(sys))
-            filter(x -> istree(x) && isoutput(x), rhss) # observed can return equations with complicated expressions, we are only looking for single Terms
-            filter(x -> istree(x) && isoutput(x), lhss)])
+            filter(x -> iscall(x) && isoutput(x), rhss) # observed can return equations with complicated expressions, we are only looking for single Terms
+            filter(x -> iscall(x) && isoutput(x), lhss)])
 end
 
 """

src/parameters.jl

@@ -15,9 +15,9 @@ function isparameter(x)
         isparameter(p) ||
             (hasmetadata(p, Symbolics.VariableSource) &&
              getmetadata(p, Symbolics.VariableSource)[1] == :parameters)
-    elseif istree(x) && operation(x) isa Symbolic
+    elseif iscall(x) && operation(x) isa Symbolic
         varT === PARAMETER || isparameter(operation(x))
-    elseif istree(x) && operation(x) == (getindex)
+    elseif iscall(x) && operation(x) == (getindex)
         isparameter(arguments(x)[1])
     elseif x isa Symbolic
         varT === PARAMETER

src/structural_transformation/StructuralTransformations.jl

@@ -7,7 +7,7 @@ using Symbolics: unwrap, linear_expansion, fast_substitute
 using SymbolicUtils
 using SymbolicUtils.Code
 using SymbolicUtils.Rewriters
-using SymbolicUtils: similarterm, istree
+using SymbolicUtils: maketerm, iscall
 
 using ModelingToolkit
 using ModelingToolkit: ODESystem, AbstractSystem, var_from_nested_derivative, Differential,

src/structural_transformation/symbolics_tearing.jl

@@ -92,7 +92,7 @@ function full_equations(sys::AbstractSystem; simplify = false)
     @unpack subs = substitutions
     solved = Dict(eq.lhs => eq.rhs for eq in subs)
     neweqs = map(equations(sys)) do eq
-        if istree(eq.lhs) && operation(eq.lhs) isa Union{Shift, Differential}
+        if iscall(eq.lhs) && operation(eq.lhs) isa Union{Shift, Differential}
             return tearing_sub(eq.lhs, solved, simplify) ~ tearing_sub(eq.rhs, solved,
                 simplify)
         else
@@ -568,7 +568,7 @@ function tearing_reassemble(state::TearingState, var_eq_matching,
 
     for eq in obs
         lhs = eq.lhs
-        istree(lhs) || continue
+        iscall(lhs) || continue
         operation(lhs) === getindex || continue
         Symbolics.shape(lhs) !== Symbolics.Unknown() || continue
         arg1 = arguments(lhs)[1]

src/structural_transformation/utils.jl

@@ -445,7 +445,7 @@ function simplify_shifts(var)
         t2 = op2.t
         return simplify_shifts(ModelingToolkit.Shift(t1 === nothing ? t2 : t1, s1 + s2)(vv2))
     else
-        return similarterm(var, operation(var), simplify_shifts.(arguments(var)),
-            Symbolics.symtype(var); metadata = unwrap(var).metadata)
+        return maketerm(typeof(var), operation(var), simplify_shifts.(arguments(var)),
+            Symbolics.symtype(var), unwrap(var).metadata)
     end
 end

src/systems/abstractsystem.jl

@@ -205,7 +205,7 @@ function wrap_array_vars(sys::AbstractSystem, exprs; dvs = unknowns(sys))
     isscalar = !(exprs isa AbstractArray)
     array_vars = Dict{Any, AbstractArray{Int}}()
     for (j, x) in enumerate(dvs)
-        if istree(x) && operation(x) == getindex
+        if iscall(x) && operation(x) == getindex
             arg = arguments(x)[1]
             inds = get!(() -> Int[], array_vars, arg)
             push!(inds, j)
@@ -352,7 +352,7 @@ function SymbolicIndexingInterface.is_variable(sys::AbstractSystem, sym)
     end
     if has_index_cache(sys) && (ic = get_index_cache(sys)) !== nothing
         return is_variable(ic, sym) ||
-               istree(sym) && operation(sym) === getindex &&
+               iscall(sym) && operation(sym) === getindex &&
                is_variable(ic, first(arguments(sym)))
     end
     return any(isequal(sym), variable_symbols(sys)) ||
@@ -378,7 +378,7 @@ function SymbolicIndexingInterface.variable_index(sys::AbstractSystem, sym)
     if has_index_cache(sys) && (ic = get_index_cache(sys)) !== nothing
         return if (idx = variable_index(ic, sym)) !== nothing
             idx
-        elseif istree(sym) && operation(sym) === getindex &&
+        elseif iscall(sym) && operation(sym) === getindex &&
                (idx = variable_index(ic, first(arguments(sym)))) !== nothing
             idx[arguments(sym)[(begin + 1):end]...]
         else
@@ -416,7 +416,7 @@ function SymbolicIndexingInterface.is_parameter(sys::AbstractSystem, sym)
     sym = unwrap(sym)
     if has_index_cache(sys) && (ic = get_index_cache(sys)) !== nothing
         return is_parameter(ic, sym) ||
-               istree(sym) && operation(sym) === getindex &&
+               iscall(sym) && operation(sym) === getindex &&
                is_parameter(ic, first(arguments(sym)))
     end
     if unwrap(sym) isa Int
@@ -441,7 +441,7 @@ function SymbolicIndexingInterface.parameter_index(sys::AbstractSystem, sym)
     if has_index_cache(sys) && (ic = get_index_cache(sys)) !== nothing
         return if (idx = parameter_index(ic, sym)) !== nothing
             idx
-        elseif istree(sym) && operation(sym) === getindex &&
+        elseif iscall(sym) && operation(sym) === getindex &&
                (idx = parameter_index(ic, first(arguments(sym)))) !== nothing
             ParameterIndex(idx.portion, (idx.idx..., arguments(sym)[(begin + 1):end]...))
         else
@@ -745,22 +745,22 @@ function _apply_to_variables(f::F, ex) where {F}
     if isvariable(ex)
         return f(ex)
     end
-    istree(ex) || return ex
-    similarterm(ex, _apply_to_variables(f, operation(ex)),
+    iscall(ex) || return ex
+    maketerm(typeof(ex), _apply_to_variables(f, operation(ex)),
         map(Base.Fix1(_apply_to_variables, f), arguments(ex)),
-        metadata = metadata(ex))
+        symtype(ex), metadata(ex))
 end
 
 abstract type SymScope end
 
 struct LocalScope <: SymScope end
 function LocalScope(sym::Union{Num, Symbolic, Symbolics.Arr{Num}})
     apply_to_variables(sym) do sym
-        if istree(sym) && operation(sym) === getindex
+        if iscall(sym) && operation(sym) === getindex
             args = arguments(sym)
             a1 = setmetadata(args[1], SymScope, LocalScope())
-            similarterm(sym, operation(sym), [a1, args[2:end]...];
-                metadata = metadata(sym))
+            maketerm(typeof(sym), operation(sym), [a1, args[2:end]...],
+                symtype(sym), metadata(sym))
         else
             setmetadata(sym, SymScope, LocalScope())
         end
@@ -772,12 +772,12 @@ struct ParentScope <: SymScope
 end
 function ParentScope(sym::Union{Num, Symbolic, Symbolics.Arr{Num}})
     apply_to_variables(sym) do sym
-        if istree(sym) && operation(sym) === getindex
+        if iscall(sym) && operation(sym) === getindex
             args = arguments(sym)
             a1 = setmetadata(args[1], SymScope,
                 ParentScope(getmetadata(value(args[1]), SymScope, LocalScope())))
-            similarterm(sym, operation(sym), [a1, args[2:end]...];
-                metadata = metadata(sym))
+            maketerm(typeof(sym), operation(sym), [a1, args[2:end]...],
+                symtype(sym), metadata(sym))
         else
             setmetadata(sym, SymScope,
                 ParentScope(getmetadata(value(sym), SymScope, LocalScope())))
@@ -791,12 +791,12 @@ struct DelayParentScope <: SymScope
 end
 function DelayParentScope(sym::Union{Num, Symbolic, Symbolics.Arr{Num}}, N)
     apply_to_variables(sym) do sym
-        if istree(sym) && operation(sym) == getindex
+        if iscall(sym) && operation(sym) == getindex
             args = arguments(sym)
             a1 = setmetadata(args[1], SymScope,
                 DelayParentScope(getmetadata(value(args[1]), SymScope, LocalScope()), N))
-            similarterm(sym, operation(sym), [a1, args[2:end]...];
-                metadata = metadata(sym))
+            maketerm(typeof(sym), operation(sym), [a1, args[2:end]...],
+                symtype(sym), metadata(sym))
         else
             setmetadata(sym, SymScope,
                 DelayParentScope(getmetadata(value(sym), SymScope, LocalScope()), N))
@@ -808,11 +808,11 @@ DelayParentScope(sym::Union{Num, Symbolic, Symbolics.Arr{Num}}) = DelayParentSco
 struct GlobalScope <: SymScope end
 function GlobalScope(sym::Union{Num, Symbolic, Symbolics.Arr{Num}})
     apply_to_variables(sym) do sym
-        if istree(sym) && operation(sym) == getindex
+        if iscall(sym) && operation(sym) == getindex
             args = arguments(sym)
             a1 = setmetadata(args[1], SymScope, GlobalScope())
-            similarterm(sym, operation(sym), [a1, args[2:end]...];
-                metadata = metadata(sym))
+            maketerm(typeof(sym), operation(sym), [a1, args[2:end]...],
+                symtype(sym), metadata(sym))
         else
             setmetadata(sym, SymScope, GlobalScope())
         end
@@ -827,16 +827,16 @@ function renamespace(sys, x)
     x = unwrap(x)
     if x isa Symbolic
         T = typeof(x)
-        if istree(x) && operation(x) isa Operator
-            return similarterm(x, operation(x),
-                Any[renamespace(sys, only(arguments(x)))];
-                metadata = metadata(x))::T
+        if iscall(x) && operation(x) isa Operator
+            return maketerm(typeof(x), operation(x),
+                Any[renamespace(sys, only(arguments(x)))],
+                symtype(x), metadata(x))::T
         end
-        if istree(x) && operation(x) === getindex
+        if iscall(x) && operation(x) === getindex
             args = arguments(x)
-            return similarterm(
-                x, operation(x), vcat(renamespace(sys, args[1]), args[2:end]);
-                metadata = metadata(x))::T
+            return maketerm(
+                typeof(x), operation(x), vcat(renamespace(sys, args[1]), args[2:end]),
+                symtype(x), metadata(x))::T
         end
         let scope = getmetadata(x, SymScope, LocalScope())
             if scope isa LocalScope
@@ -904,7 +904,7 @@ function namespace_expr(
     O = unwrap(O)
     if any(isequal(O), ivs)
         return O
-    elseif istree(O)
+    elseif iscall(O)
         T = typeof(O)
         renamed = let sys = sys, n = n, T = T
             map(a -> namespace_expr(a, sys, n; ivs)::Any, arguments(O))
@@ -913,13 +913,14 @@ function namespace_expr(
             # Use renamespace so the scope is correct, and make sure to use the
             # metadata from the rescoped variable
             rescoped = renamespace(n, O)
-            similarterm(O, operation(rescoped), renamed,
-                metadata = metadata(rescoped))
+            maketerm(typeof(rescoped), operation(rescoped), renamed,
+                symtype(rescoped),
+                metadata(rescoped))
         elseif Symbolics.isarraysymbolic(O)
             # promote_symtype doesn't work for array symbolics
-            similarterm(O, operation(O), renamed, symtype(O), metadata = metadata(O))
+            maketerm(typeof(O), operation(O), renamed, symtype(O), metadata(O))
         else
-            similarterm(O, operation(O), renamed, metadata = metadata(O))
+            maketerm(typeof(O), operation(O), renamed, symtype(O), metadata(O))
         end
     elseif isvariable(O)
         renamespace(n, O)
@@ -1097,7 +1098,7 @@ function time_varying_as_func(x, sys::AbstractTimeDependentSystem)
     # than pass in a value in place of x(t).
     #
     # This is done by just making `x` the argument of the function.
-    if istree(x) &&
+    if iscall(x) &&
        issym(operation(x)) &&
        !(length(arguments(x)) == 1 && isequal(arguments(x)[1], get_iv(sys)))
         return operation(x)
@@ -1125,7 +1126,7 @@ function push_vars!(stmt, name, typ, vars)
     isempty(vars) && return
     vars_expr = Expr(:macrocall, typ, nothing)
     for s in vars
-        if istree(s)
+        if iscall(s)
             f = nameof(operation(s))
             args = arguments(s)
             ex = :($f($(args...)))
@@ -1142,7 +1143,7 @@ function round_trip_expr(t, var2name)
     name = get(var2name, t, nothing)
     name !== nothing && return name
     issym(t) && return nameof(t)
-    istree(t) || return t
+    iscall(t) || return t
     f = round_trip_expr(operation(t), var2name)
     args = map(Base.Fix2(round_trip_expr, var2name), arguments(t))
     return :($f($(args...)))
@@ -2564,8 +2565,10 @@ is_diff_equation(eq2) # false
 """
 function is_diff_equation(eq)
     (eq isa Equation) || (return false)
-    isdefined(eq, :lhs) && occursin(is_derivative, wrap(eq.lhs)) && (return true)
-    isdefined(eq, :rhs) && occursin(is_derivative, wrap(eq.rhs)) && (return true)
+    isdefined(eq, :lhs) && hasnode(is_derivative, wrap(eq.lhs)) &&
+        (return true)
+    isdefined(eq, :rhs) && hasnode(is_derivative, wrap(eq.rhs)) &&
+        (return true)
     return false
 end
 

src/systems/clock_inference.jl

@@ -48,7 +48,7 @@ function substitute_sample_time(eq::Equation, dt)
 end
 
 function substitute_sample_time(ex, dt)
-    istree(ex) || return ex
+    iscall(ex) || return ex
     op = operation(ex)
     args = arguments(ex)
     if op == SampleTime
@@ -63,7 +63,7 @@ function substitute_sample_time(ex, dt)
             end
             new_args[i] = ex_arg
         end
-        similarterm(ex, op, new_args, symtype(ex); metadata = metadata(ex))
+        maketerm(typeof(ex), op, new_args, symtype(ex), metadata(ex))
     end
 end
 
@@ -128,7 +128,7 @@ function resize_or_push!(v, val, idx)
 end
 
 function is_time_domain_conversion(v)
-    istree(v) && (o = operation(v)) isa Operator &&
+    iscall(v) && (o = operation(v)) isa Operator &&
         input_timedomain(o) != output_timedomain(o)
 end
 

src/systems/connectors.jl

@@ -14,7 +14,7 @@ end
 
 function get_connection_type(s)
     s = unwrap(s)
-    if istree(s) && operation(s) === getindex
+    if iscall(s) && operation(s) === getindex
         s = arguments(s)[1]
     end
     getmetadata(s, VariableConnectType, Equality)
@@ -82,7 +82,7 @@ function collect_instream!(set, eq::Equation)
 end
 
 function collect_instream!(set, expr, occurs = false)
-    istree(expr) || return occurs
+    iscall(expr) || return occurs
     op = operation(expr)
     op === instream && (push!(set, expr); occurs = true)
     for a in SymbolicUtils.unsorted_arguments(expr)