179 lines
6.7 KiB
Python
179 lines
6.7 KiB
Python
|
"""Tests for solvers of systems of polynomial equations. """
|
||
|
|
from sympy.core.numbers import (I, Integer, Rational)
|
||
|
|
from sympy.core.singleton import S
|
||
|
|
from sympy.core.symbol import symbols
|
||
|
|
from sympy.functions.elementary.miscellaneous import sqrt
|
||
|
|
from sympy.polys.domains.rationalfield import QQ
|
||
|
|
from sympy.polys.polyerrors import UnsolvableFactorError
|
||
|
|
from sympy.polys.polyoptions import Options
|
||
|
|
from sympy.polys.polytools import Poly
|
||
|
|
from sympy.solvers.solvers import solve
|
||
|
|
from sympy.utilities.iterables import flatten
|
||
|
|
from sympy.abc import x, y, z
|
||
|
|
from sympy.polys import PolynomialError
|
||
|
|
from sympy.solvers.polysys import (solve_poly_system,
|
||
|
|
solve_triangulated,
|
||
|
|
solve_biquadratic, SolveFailed,
|
||
|
|
solve_generic)
|
||
|
|
from sympy.polys.polytools import parallel_poly_from_expr
|
||
|
|
from sympy.testing.pytest import raises
|
||
|
|
|
||
|
|
|
||
|
|
def test_solve_poly_system():
|
||
|
|
assert solve_poly_system([x - 1], x) == [(S.One,)]
|
||
|
|
|
||
|
|
assert solve_poly_system([y - x, y - x - 1], x, y) is None
|
||
|
|
|
||
|
|
assert solve_poly_system([y - x**2, y + x**2], x, y) == [(S.Zero, S.Zero)]
|
||
|
|
|
||
|
|
assert solve_poly_system([2*x - 3, y*Rational(3, 2) - 2*x, z - 5*y], x, y, z) == \
|
||
|
|
[(Rational(3, 2), Integer(2), Integer(10))]
|
||
|
|
|
||
|
|
assert solve_poly_system([x*y - 2*y, 2*y**2 - x**2], x, y) == \
|
||
|
|
[(0, 0), (2, -sqrt(2)), (2, sqrt(2))]
|
||
|
|
|
||
|
|
assert solve_poly_system([y - x**2, y + x**2 + 1], x, y) == \
|
||
|
|
[(-I*sqrt(S.Half), Rational(-1, 2)), (I*sqrt(S.Half), Rational(-1, 2))]
|
||
|
|
|
||
|
|
f_1 = x**2 + y + z - 1
|
||
|
|
f_2 = x + y**2 + z - 1
|
||
|
|
f_3 = x + y + z**2 - 1
|
||
|
|
|
||
|
|
a, b = sqrt(2) - 1, -sqrt(2) - 1
|
||
|
|
|
||
|
|
assert solve_poly_system([f_1, f_2, f_3], x, y, z) == \
|
||
|
|
[(0, 0, 1), (0, 1, 0), (1, 0, 0), (a, a, a), (b, b, b)]
|
||
|
|
|
||
|
|
solution = [(1, -1), (1, 1)]
|
||
|
|
|
||
|
|
assert solve_poly_system([Poly(x**2 - y**2), Poly(x - 1)]) == solution
|
||
|
|
assert solve_poly_system([x**2 - y**2, x - 1], x, y) == solution
|
||
|
|
assert solve_poly_system([x**2 - y**2, x - 1]) == solution
|
||
|
|
|
||
|
|
assert solve_poly_system(
|
||
|
|
[x + x*y - 3, y + x*y - 4], x, y) == [(-3, -2), (1, 2)]
|
||
|
|
|
||
|
|
raises(NotImplementedError, lambda: solve_poly_system([x**3 - y**3], x, y))
|
||
|
|
raises(NotImplementedError, lambda: solve_poly_system(
|
||
|
|
[z, -2*x*y**2 + x + y**2*z, y**2*(-z - 4) + 2]))
|
||
|
|
raises(PolynomialError, lambda: solve_poly_system([1/x], x))
|
||
|
|
|
||
|
|
raises(NotImplementedError, lambda: solve_poly_system(
|
||
|
|
[x-1,], (x, y)))
|
||
|
|
raises(NotImplementedError, lambda: solve_poly_system(
|
||
|
|
[y-1,], (x, y)))
|
||
|
|
|
||
|
|
# solve_poly_system should ideally construct solutions using
|
||
|
|
# CRootOf for the following four tests
|
||
|
|
assert solve_poly_system([x**5 - x + 1], [x], strict=False) == []
|
||
|
|
raises(UnsolvableFactorError, lambda: solve_poly_system(
|
||
|
|
[x**5 - x + 1], [x], strict=True))
|
||
|
|
|
||
|
|
assert solve_poly_system([(x - 1)*(x**5 - x + 1), y**2 - 1], [x, y],
|
||
|
|
strict=False) == [(1, -1), (1, 1)]
|
||
|
|
raises(UnsolvableFactorError,
|
||
|
|
lambda: solve_poly_system([(x - 1)*(x**5 - x + 1), y**2-1],
|
||
|
|
[x, y], strict=True))
|
||
|
|
|
||
|
|
|
||
|
|
def test_solve_generic():
|
||
|
|
NewOption = Options((x, y), {'domain': 'ZZ'})
|
||
|
|
assert solve_generic([x**2 - 2*y**2, y**2 - y + 1], NewOption) == \
|
||
|
|
[(-sqrt(-1 - sqrt(3)*I), Rational(1, 2) - sqrt(3)*I/2),
|
||
|
|
(sqrt(-1 - sqrt(3)*I), Rational(1, 2) - sqrt(3)*I/2),
|
||
|
|
(-sqrt(-1 + sqrt(3)*I), Rational(1, 2) + sqrt(3)*I/2),
|
||
|
|
(sqrt(-1 + sqrt(3)*I), Rational(1, 2) + sqrt(3)*I/2)]
|
||
|
|
|
||
|
|
# solve_generic should ideally construct solutions using
|
||
|
|
# CRootOf for the following two tests
|
||
|
|
assert solve_generic(
|
||
|
|
[2*x - y, (y - 1)*(y**5 - y + 1)], NewOption, strict=False) == \
|
||
|
|
[(Rational(1, 2), 1)]
|
||
|
|
raises(UnsolvableFactorError, lambda: solve_generic(
|
||
|
|
[2*x - y, (y - 1)*(y**5 - y + 1)], NewOption, strict=True))
|
||
|
|
|
||
|
|
|
||
|
|
def test_solve_biquadratic():
|
||
|
|
x0, y0, x1, y1, r = symbols('x0 y0 x1 y1 r')
|
||
|
|
|
||
|
|
f_1 = (x - 1)**2 + (y - 1)**2 - r**2
|
||
|
|
f_2 = (x - 2)**2 + (y - 2)**2 - r**2
|
||
|
|
s = sqrt(2*r**2 - 1)
|
||
|
|
a = (3 - s)/2
|
||
|
|
b = (3 + s)/2
|
||
|
|
assert solve_poly_system([f_1, f_2], x, y) == [(a, b), (b, a)]
|
||
|
|
|
||
|
|
f_1 = (x - 1)**2 + (y - 2)**2 - r**2
|
||
|
|
f_2 = (x - 1)**2 + (y - 1)**2 - r**2
|
||
|
|
|
||
|
|
assert solve_poly_system([f_1, f_2], x, y) == \
|
||
|
|
[(1 - sqrt((2*r - 1)*(2*r + 1))/2, Rational(3, 2)),
|
||
|
|
(1 + sqrt((2*r - 1)*(2*r + 1))/2, Rational(3, 2))]
|
||
|
|
|
||
|
|
query = lambda expr: expr.is_Pow and expr.exp is S.Half
|
||
|
|
|
||
|
|
f_1 = (x - 1 )**2 + (y - 2)**2 - r**2
|
||
|
|
f_2 = (x - x1)**2 + (y - 1)**2 - r**2
|
||
|
|
|
||
|
|
result = solve_poly_system([f_1, f_2], x, y)
|
||
|
|
|
||
|
|
assert len(result) == 2 and all(len(r) == 2 for r in result)
|
||
|
|
assert all(r.count(query) == 1 for r in flatten(result))
|
||
|
|
|
||
|
|
f_1 = (x - x0)**2 + (y - y0)**2 - r**2
|
||
|
|
f_2 = (x - x1)**2 + (y - y1)**2 - r**2
|
||
|
|
|
||
|
|
result = solve_poly_system([f_1, f_2], x, y)
|
||
|
|
|
||
|
|
assert len(result) == 2 and all(len(r) == 2 for r in result)
|
||
|
|
assert all(len(r.find(query)) == 1 for r in flatten(result))
|
||
|
|
|
||
|
|
s1 = (x*y - y, x**2 - x)
|
||
|
|
assert solve(s1) == [{x: 1}, {x: 0, y: 0}]
|
||
|
|
s2 = (x*y - x, y**2 - y)
|
||
|
|
assert solve(s2) == [{y: 1}, {x: 0, y: 0}]
|
||
|
|
gens = (x, y)
|
||
|
|
for seq in (s1, s2):
|
||
|
|
(f, g), opt = parallel_poly_from_expr(seq, *gens)
|
||
|
|
raises(SolveFailed, lambda: solve_biquadratic(f, g, opt))
|
||
|
|
seq = (x**2 + y**2 - 2, y**2 - 1)
|
||
|
|
(f, g), opt = parallel_poly_from_expr(seq, *gens)
|
||
|
|
assert solve_biquadratic(f, g, opt) == [
|
||
|
|
(-1, -1), (-1, 1), (1, -1), (1, 1)]
|
||
|
|
ans = [(0, -1), (0, 1)]
|
||
|
|
seq = (x**2 + y**2 - 1, y**2 - 1)
|
||
|
|
(f, g), opt = parallel_poly_from_expr(seq, *gens)
|
||
|
|
assert solve_biquadratic(f, g, opt) == ans
|
||
|
|
seq = (x**2 + y**2 - 1, x**2 - x + y**2 - 1)
|
||
|
|
(f, g), opt = parallel_poly_from_expr(seq, *gens)
|
||
|
|
assert solve_biquadratic(f, g, opt) == ans
|
||
|
|
|
||
|
|
|
||
|
|
def test_solve_triangulated():
|
||
|
|
f_1 = x**2 + y + z - 1
|
||
|
|
f_2 = x + y**2 + z - 1
|
||
|
|
f_3 = x + y + z**2 - 1
|
||
|
|
|
||
|
|
a, b = sqrt(2) - 1, -sqrt(2) - 1
|
||
|
|
|
||
|
|
assert solve_triangulated([f_1, f_2, f_3], x, y, z) == \
|
||
|
|
[(0, 0, 1), (0, 1, 0), (1, 0, 0)]
|
||
|
|
|
||
|
|
dom = QQ.algebraic_field(sqrt(2))
|
||
|
|
|
||
|
|
assert solve_triangulated([f_1, f_2, f_3], x, y, z, domain=dom) == \
|
||
|
|
[(0, 0, 1), (0, 1, 0), (1, 0, 0), (a, a, a), (b, b, b)]
|
||
|
|
|
||
|
|
|
||
|
|
def test_solve_issue_3686():
|
||
|
|
roots = solve_poly_system([((x - 5)**2/250000 + (y - Rational(5, 10))**2/250000) - 1, x], x, y)
|
||
|
|
assert roots == [(0, S.Half - 15*sqrt(1111)), (0, S.Half + 15*sqrt(1111))]
|
||
|
|
|
||
|
|
roots = solve_poly_system([((x - 5)**2/250000 + (y - 5.0/10)**2/250000) - 1, x], x, y)
|
||
|
|
# TODO: does this really have to be so complicated?!
|
||
|
|
assert len(roots) == 2
|
||
|
|
assert roots[0][0] == 0
|
||
|
|
assert roots[0][1].epsilon_eq(-499.474999374969, 1e12)
|
||
|
|
assert roots[1][0] == 0
|
||
|
|
assert roots[1][1].epsilon_eq(500.474999374969, 1e12)
|