fix: fix flatten_equations

Author: AayushSabharwal

src/systems/diffeqs/abstractodesystem.jl

@@ -1513,19 +1513,7 @@ function isisomorphic(sys1::AbstractODESystem, sys2::AbstractODESystem)
 end
 
 function flatten_equations(eqs)
-    mapreduce(vcat, eqs; 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
+    mapreduce(Symbolics.scalarize, vcat, eqs; init = Equation[])
 end
 
 struct InitializationProblem{iip, specialization} end

test/odesystem.jl

@@ -1022,3 +1022,10 @@ prob2 = @test_nowarn ODEProblem(sys, [x => ones(3)], (0.0, 10.0), [p => ones(3,
 sol2 = @test_nowarn solve(prob, Tsit5())
 
 @test sol1 ≈ sol2
+
+# Flatten equations properly
+# Issue#1502
+@variables x(t)[1:3] y(t)
+@parameters p[1:3, 1:3]
+@test_nowarn @mtkbuild sys = ODESystem([D(x) ~ p * x, D(y) ~ x' * p * x], t)
+@test_nowarn ODEProblem(sys, [x => ones(3), y => 2], (0.0, 10.0), [p => ones(3, 3)])