Skip to content

Commit

Permalink
Merge pull request #1259 from n0rbed/master
Browse files Browse the repository at this point in the history
removed solve_multipoly
  • Loading branch information
ChrisRackauckas authored Sep 7, 2024
2 parents fa98d3d + d292e16 commit 97d9cad
Show file tree
Hide file tree
Showing 4 changed files with 18 additions and 87 deletions.
19 changes: 0 additions & 19 deletions ext/SymbolicsNemoExt.jl
Original file line number Diff line number Diff line change
Expand Up @@ -57,32 +57,13 @@ function Symbolics.factor_use_nemo(poly::Num)
return sym_unit, sym_factors
end

# gcd(x^2 - y^2, x^3 - y^3) -> x - y
function Symbolics.gcd_use_nemo(poly1::Num, poly2::Num)
Symbolics.check_polynomial(poly1)
Symbolics.check_polynomial(poly2)
vars1 = Symbolics.get_variables(poly1)
vars2 = Symbolics.get_variables(poly2)
vars = vcat(vars1, vars2)
nemo_ring, nemo_vars = Nemo.polynomial_ring(Nemo.QQ, map(string, vars))
sym_to_nemo = Dict(vars .=> nemo_vars)
nemo_to_sym = Dict(v => k for (k, v) in sym_to_nemo)
nemo_poly1 = Symbolics.substitute(poly1, sym_to_nemo)
nemo_poly2 = Symbolics.substitute(poly2, sym_to_nemo)
nemo_gcd = Nemo.gcd(nemo_poly1, nemo_poly2)
sym_gcd = Symbolics.wrap(nemo_crude_evaluate(nemo_gcd, nemo_to_sym))
return sym_gcd
end


# Helps with precompilation time
PrecompileTools.@setup_workload begin
@variables a b c x y z
expr_with_params = expand((x + b)*(x^2 + 2x + 1)*(x^2 - a))
PrecompileTools.@compile_workload begin
symbolic_solve(expr_with_params, x, dropmultiplicity=false)
symbolic_solve(x^10 - a^10, x, dropmultiplicity=false)
symbolic_solve([x^2 - a^2, x + a], x)
end
end

Expand Down
71 changes: 15 additions & 56 deletions src/solver/main.jl
Original file line number Diff line number Diff line change
Expand Up @@ -165,31 +165,21 @@ function symbolic_solve(expr, x::T; dropmultiplicity = true, warns = true) where
expr = [expr]
expr_univar = false
end
if !expr_univar && x_univar
x = [x]
x_univar = false
end

if x_univar
sols = []
if expr_univar
sols = check_poly_inunivar(expr, x) ?
solve_univar(expr, x, dropmultiplicity=dropmultiplicity) :
ia_solve(expr, x, warns=warns)
isequal(sols, nothing) && return nothing
else
for i in eachindex(expr)
if !check_poly_inunivar(expr[i], x)
warns && @warn("Solve can not solve this input currently")
return nothing
end
end
sols = solve_multipoly(
expr, x, dropmultiplicity=dropmultiplicity, warns=warns)
isequal(sols, nothing) && return nothing
end

sols = check_poly_inunivar(expr, x) ?
solve_univar(expr, x, dropmultiplicity = dropmultiplicity) :
ia_solve(expr, x, warns = warns)
isequal(sols, nothing) && return nothing
sols = map(postprocess_root, sols)
return sols
end

if !expr_univar && !x_univar
if !x_univar
for e in expr
for var in x
if !check_poly_inunivar(e, var)
Expand All @@ -201,11 +191,13 @@ function symbolic_solve(expr, x::T; dropmultiplicity = true, warns = true) where

sols = solve_multivar(expr, x, dropmultiplicity=dropmultiplicity, warns=warns)
isequal(sols, nothing) && return nothing
for sol in sols
sols = convert(Vector{Any}, sols)
for i in eachindex(sols)
for var in x
if haskey(sol, var)
sol[var] = postprocess_root(sol[var])
end
sols[i][var] = postprocess_root(sols[i][var])
end
if length(collect(keys(sols[i]))) == 1
sols[i] = collect(values(sols[i]))[1]
end
end

Expand Down Expand Up @@ -310,39 +302,6 @@ function solve_univar(expression, x; dropmultiplicity=true)
return arr_roots
end

# You can compute the GCD between a system of polynomials by doing the following:
# Get the GCD between the first two polys,
# and get the GCD between this result and the following index,
# say: solve([x^2 - 1, x - 1, (x-1)^20], x)
# the GCD between the first two terms is obviously x-1,
# now we call gcd_use_nemo() on this term, and the following,
# gcd_use_nemo(x - 1, (x-1)^20), which is again x-1.
# now we just need to solve(x-1, x) to get the common root in this
# system of equations.
function solve_multipoly(polys::Vector, x::Num; dropmultiplicity = true, warns = true)
polys = unique(polys)

if length(polys) < 1
warns && @warn("No expressions entered")
return nothing
end
if length(polys) == 1
return solve_univar(polys[1], x, dropmultiplicity = dropmultiplicity)
end

gcd = gcd_use_nemo(polys[1], polys[2])

for i in eachindex(polys)[3:end]
gcd = gcd_use_nemo(gcd, polys[i])
end

if isequal(gcd, 1)
return []
end

return solve_univar(gcd, x, dropmultiplicity = dropmultiplicity)
end

function solve_multivar(eqs::Any, vars::Any; dropmultiplicity = true, warns = true)
throw("Groebner bases engine is required. Execute `using Groebner` to enable this functionality.")
end
4 changes: 0 additions & 4 deletions src/solver/nemo_stuff.jl
Original file line number Diff line number Diff line change
Expand Up @@ -14,7 +14,3 @@ function factor_use_nemo(poly::Any)
throw("Nemo is required. Execute `using Nemo` to enable this functionality.")
end

# gcd(x^2 - y^2, x^3 - y^3) -> x - y
function gcd_use_nemo(poly1::Any, poly2::Any)
throw("Nemo is required. Execute `using Nemo` to enable this functionality.")
end
11 changes: 3 additions & 8 deletions test/solver.jl
Original file line number Diff line number Diff line change
Expand Up @@ -269,7 +269,7 @@ end

@testset "Multivar parametric" begin
@variables x y a
@test isequal(symbolic_solve([x + a, a - 1], x), [])
@test isequal(symbolic_solve([x + a, a - 1], x), [-1])
@test isequal(symbolic_solve([x - a, y + a], [x, y]), [Dict(y => -a, x => a)])
@test isequal(symbolic_solve([x*y - a, x*y + x], [x, y]), [Dict(y => -1, x => -a)])
@test isequal(symbolic_solve([x*y - a, 1 ~ 3], [x, y]), [])
Expand All @@ -283,13 +283,13 @@ end
@test isequal(sol, [Dict(t => -5 / (-1 + u + v), w => 1 - u - v)])

sol = symbolic_solve([x-y, y-z], [x])
@test isequal(sol, [Dict(x=>z)])
@test isequal(sol, [z])

sol = symbolic_solve([x-y, y-z], [x, y])
@test isequal(sol, [Dict(x=>z, y=>z)])

sol = symbolic_solve([x + y - z, y - z], [x])
@test isequal(sol, [Dict(x=>0)])
@test isequal(sol, [0])
end

@testset "Factorisation" begin
Expand All @@ -310,11 +310,6 @@ end
@test isequal(expand(u*prod(factors) - f), 0)
end

@testset "GCD" begin
f1, f2 = x^2 - y^2, x^3 - y^3
@test isequal(x - y, Symbolics.gcd_use_nemo(f1, f2))
end


# Post Process roots #
@testset "Post Process roots" begin
Expand Down

0 comments on commit 97d9cad

Please sign in to comment.