fix: do not scalarize in system constructors

Author: AayushSabharwal

src/systems/diffeqs/odesystem.jl

@@ -202,23 +202,10 @@ function ODESystem(deqs::AbstractVector{<:Equation}, iv, dvs, ps;
     name === nothing &&
         throw(ArgumentError("The `name` keyword must be provided. Please consider using the `@named` macro"))
     @assert all(control -> any(isequal.(control, ps)), controls) "All controls must also be parameters."
-    deqs = mapreduce(vcat, deqs; init = Equation[]) do eq
-        islhsarr = eq.lhs isa AbstractArray || Symbolics.isarraysymbolic(eq.lhs)
-        isrhsarr = eq.rhs isa AbstractArray || Symbolics.isarraysymbolic(eq.rhs)
-        if islhsarr || isrhsarr
-            islhsarr && isrhsarr ||
-                error("LHS ($(eq.lhs)) and RHS ($(eq.rhs)) must either both be array expressions or both scalar")
-            size(eq.lhs) == size(eq.rhs) ||
-                error("Size of LHS ($(eq.lhs)) and RHS ($(eq.rhs)) must match: got $(size(eq.lhs)) and $(size(eq.rhs))")
-            return collect(eq.lhs) .~ collect(eq.rhs)
-        else
-            eq
-        end
-    end
     iv′ = value(iv)
     ps′ = value.(ps)
     ctrl′ = value.(controls)
-    dvs′ = value.(symbolic_type(dvs) === NotSymbolic() ? dvs : [dvs])
+    dvs′ = value.(dvs)
     dvs′ = filter(x -> !isdelay(x, iv), dvs′)
     if !(isempty(default_u0) && isempty(default_p))
         Base.depwarn(

src/systems/diffeqs/sdesystem.jl

@@ -171,14 +171,6 @@ function SDESystem(deqs::AbstractVector{<:Equation}, neqs::AbstractArray, iv, dv
         gui_metadata = nothing)
     name === nothing &&
         throw(ArgumentError("The `name` keyword must be provided. Please consider using the `@named` macro"))
-    deqs = flatten_equations(deqs)
-    neqs = mapreduce(vcat, neqs) do expr
-        if expr isa AbstractArray || Symbolics.isarraysymbolic(expr)
-            collect(expr)
-        else
-            expr
-        end
-    end
     iv′ = value(iv)
     dvs′ = value.(dvs)
     ps′ = value.(ps)

src/systems/jumps/jumpsystem.jl

@@ -151,6 +151,7 @@ function JumpSystem(eqs, iv, unknowns, ps;
         kwargs...)
     name === nothing &&
         throw(ArgumentError("The `name` keyword must be provided. Please consider using the `@named` macro"))
+    eqs = scalarize.(eqs)
     sysnames = nameof.(systems)
     if length(unique(sysnames)) != length(sysnames)
         throw(ArgumentError("System names must be unique."))

src/systems/systemstructure.jl

@@ -251,7 +251,8 @@ function TearingState(sys; quick_cancel = false, check = true)
     sys = flatten(sys)
     ivs = independent_variables(sys)
     iv = length(ivs) == 1 ? ivs[1] : nothing
-    eqs = copy(equations(sys))
+    # scalarize array equations, without scalarizing arguments to registered functions
+    eqs = flatten_equations(copy(equations(sys)))
     neqs = length(eqs)
     dervaridxs = OrderedSet{Int}()
     var2idx = Dict{Any, Int}()