Symbolic linear system solving not supported by ModelKit

[luanborelli]

Suppose we want to solve a system $Ax = b$, where $A$ is a polynomial matrix and $b$ is a vector of polynomials. ModelKit seems to be unable to solve such a system. Here's an example:

using HomotopyContinuation

@polyvar x y

A1 = x^2 + y^2
A2 = y + x^2 
A3 = x*y 
A4 = y^2 

b1 = 5*y - 4*x
b2 = -3*x - y 

A = [A1 A2; A3 A4]
b = [b1; b2]

A\b

Output:

ERROR: MethodError: no method matching MultivariatePolynomials.RationalPoly{DynamicPolynomials.Polynomial{true, Int64}, DynamicPolynomials.Polynomial{true, Int64}}(::MultivariatePolynomials.RationalPoly{DynamicPolynomials.Polynomial{true, Int64}, DynamicPolynomials.Polynomial{true, Int64}})

Some of the types have been truncated in the stacktrace for improved reading. To emit complete information
in the stack trace, evaluate `TruncatedStacktraces.VERBOSE[] = true` and re-run the code.

Closest candidates are:
  MultivariatePolynomials.RationalPoly{NT, DT}(::Any, ::Any) where {NT<:(MultivariatePolynomials.AbstractPolynomialLike), DT<:(MultivariatePolynomials.AbstractPolynomialLike)} at D:\Users\b46525\.julia\packages\MultivariatePolynomials\sWAOE\src\rational.jl:5
  MultivariatePolynomials.RationalPoly{NT, DT}(::Bool) where {NT, DT} at D:\Users\b46525\.julia\packages\MultivariatePolynomials\sWAOE\src\rational.jl:10
Stacktrace:
 [1] oneunit(#unused#::Type{MultivariatePolynomials.RationalPoly{DynamicPolynomials.Polynomial{true, Int64}, DynamicPolynomials.Polynomial{true, Int64}}})
   @ Base .\number.jl:358
 [2] lutype(T::Type)
   @ LinearAlgebra D:\Users\b46525\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\LinearAlgebra\src\lu.jl:205
 [3] lu(A::Matrix{DynamicPolynomials.Polynomial{true, Int64}}, pivot::RowMaximum; check::Bool)
   @ LinearAlgebra D:\Users\b46525\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\LinearAlgebra\src\lu.jl:279
 [4] lu(A::Matrix{DynamicPolynomials.Polynomial{true, Int64}}, pivot::RowMaximum) (repeats 2 times)
   @ LinearAlgebra D:\Users\b46525\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\LinearAlgebra\src\lu.jl:278
 [5] \(A::Matrix{DynamicPolynomials.Polynomial{true, Int64}}, B::Vector{DynamicPolynomials.Polynomial{true, Int64}})
   @ LinearAlgebra D:\Users\b46525\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\LinearAlgebra\src\generic.jl:1110
 [6] top-level scope
   @ Untitled-2:16

This is a type of operation that other symbolic systems seem to handle with ease. I have already performed such an operation in Matlab, for example. In Julia, other systems like Symbolics (suggested in #456 for integration --- I endorse the suggestion!) are also capable of solving this kind of problem. Here's an example:

using Symbolics

@variables x y

A1 = x^2 + y^2
A2 = y + x^2 
A3 = x*y 
A4 = y^2 

b1 = 5*y - 4*x
b2 = -3*x - y 

A = [A1 A2; A3 A4]
b = [b1; b2]

A\b

Output:

2-element Vector{Num}:
 (5y + (-(y + x^2)*((-x*y*(5y - 4x)) / (x^2 + y^2) - 3x - y)) / (y^2 + (-x*y*(y + x^2)) / (x^2 + y^2)) - 4x) / (x^2 + y^2)
                                        ((-x*y*(5y - 4x)) / (x^2 + y^2) - 3x - y) / (y^2 + (-x*y*(y + x^2)) / (x^2 + y^2))

Obs.: I have a rather complicated system that depends on this type of operation to be built. I need to solve it using homotopy and I'm trying to rebuild it so that it's ModelKit compatible, but I'm stuck on this operation. Any suggestions on how to get around this problem (at least temporarily) would also be greatly appreciated. I have been trying to parse Symbolics expressions to ModelKit expressions, but without success so far. Thanks.