0

I'm trying to solve the following equation for variable x (in dependence on other symbols b, c, d):

from sympy import symbols, sin, cos, solve, solveset, S
  
b, c, d, x = symbols('b c d x', real=True)

expr = sin(x)*sin(c - 2*d - x) - sin(2*b - 2*d - x)*sin(c - x)  # find `x` such that `expr == 0`

There exist additional constraints which apply to the various symbols, namely:

0 < d < b < d+x < c

So I tried to use solve for incorporating these inequality constraints:

result = solve(
    [
        sin(x)*sin(c - 2*d - x) - sin(2*b - 2*d - x)*sin(c - x),
        0 < d,
        d < b,
        b < d+x,
        d+x < c,
    ],
    x,
)

However, I get the following error:

Traceback (most recent call last):
  File "/path/to/venv/lib/python3.9/site-packages/sympy/polys/polyutils.py", line 211, in _parallel_dict_from_expr_if_gens
    monom[indices[base]] = exp
KeyError: sin(c - x)

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "/path/to/venv/lib/python3.9/site-packages/sympy/solvers/inequalities.py", line 818, in _solve_inequality
    p = Poly(expr, s)
  File "/path/to/venv/lib/python3.9/site-packages/sympy/polys/polytools.py", line 181, in __new__
    return cls._from_expr(rep, opt)
  File "/path/to/venv/lib/python3.9/site-packages/sympy/polys/polytools.py", line 310, in _from_expr
    rep, opt = _dict_from_expr(rep, opt)
  File "/path/to/venv/lib/python3.9/site-packages/sympy/polys/polyutils.py", line 368, in _dict_from_expr
    rep, gens = _dict_from_expr_if_gens(expr, opt)
  File "/path/to/venv/lib/python3.9/site-packages/sympy/polys/polyutils.py", line 307, in _dict_from_expr_if_gens
    (poly,), gens = _parallel_dict_from_expr_if_gens((expr,), opt)
  File "/path/to/venv/lib/python3.9/site-packages/sympy/polys/polyutils.py", line 216, in _parallel_dict_from_expr_if_gens
    raise PolynomialError("%s contains an element of "
sympy.polys.polyerrors.PolynomialError: sin(c - x) contains an element of the set of generators.

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "/path/to/venv/lib/python3.9/site-packages/sympy/polys/polyutils.py", line 211, in _parallel_dict_from_expr_if_gens
    monom[indices[base]] = exp
KeyError: sin(c - x)

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "/path/to/venv/lib/python3.9/site-packages/sympy/solvers/inequalities.py", line 246, in reduce_rational_inequalities
    (numer, denom), opt = parallel_poly_from_expr(
  File "/path/to/venv/lib/python3.9/site-packages/sympy/polys/polytools.py", line 4414, in parallel_poly_from_expr
    return _parallel_poly_from_expr(exprs, opt)
  File "/path/to/venv/lib/python3.9/site-packages/sympy/polys/polytools.py", line 4467, in _parallel_poly_from_expr
    reps, opt = _parallel_dict_from_expr(exprs, opt)
  File "/path/to/venv/lib/python3.9/site-packages/sympy/polys/polyutils.py", line 332, in _parallel_dict_from_expr
    reps, gens = _parallel_dict_from_expr_if_gens(exprs, opt)
  File "/path/to/venv/lib/python3.9/site-packages/sympy/polys/polyutils.py", line 216, in _parallel_dict_from_expr_if_gens
    raise PolynomialError("%s contains an element of "
sympy.polys.polyerrors.PolynomialError: sin(c - x) contains an element of the set of generators.

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "/path/to/venv/lib/python3.9/site-packages/sympy/solvers/inequalities.py", line 827, in _solve_inequality
    rv = reduce_rational_inequalities([[ie]], s)
  File "/path/to/venv/lib/python3.9/site-packages/sympy/solvers/inequalities.py", line 249, in reduce_rational_inequalities
    raise PolynomialError(filldedent('''
sympy.polys.polyerrors.PolynomialError: 
only polynomials and rational functions are supported in this context.

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "/tmp/test2.py", line 5, in <module>
    result = solve(
  File "/path/to/venv/lib/python3.9/site-packages/sympy/solvers/solvers.py", line 920, in solve
    return reduce_inequalities(f, symbols=symbols)
  File "/path/to/venv/lib/python3.9/site-packages/sympy/solvers/inequalities.py", line 996, in reduce_inequalities
    rv = _reduce_inequalities(inequalities, symbols)
  File "/path/to/venv/lib/python3.9/site-packages/sympy/solvers/inequalities.py", line 912, in _reduce_inequalities
    other.append(_solve_inequality(Relational(expr, 0, rel), gen))
  File "/path/to/venv/lib/python3.9/site-packages/sympy/solvers/inequalities.py", line 829, in _solve_inequality
    rv = solve_univariate_inequality(ie, s)
  File "/path/to/venv/lib/python3.9/site-packages/sympy/solvers/inequalities.py", line 500, in solve_univariate_inequality
    frange = function_range(e, gen, domain)
  File "/path/to/venv/lib/python3.9/site-packages/sympy/calculus/util.py", line 208, in function_range
    for critical_point in critical_points:
  File "/path/to/venv/lib/python3.9/site-packages/sympy/utilities/iterables.py", line 2881, in roundrobin
    yield nxt()
  File "/path/to/venv/lib/python3.9/site-packages/sympy/sets/sets.py", line 1439, in __iter__
    raise TypeError(
TypeError: The computation had not completed because of the undecidable set membership is found in every candidates.

I'm not sure what that error means. Does it mean that there exists no solution for the given inequality constraints?

(sympy version: 1.10.1)

solveset

I also tried using solveset (however, it doesn't seem to be possible to include the inequality constraints here):

result = solveset(
    sin(x)*sin(c - 2*d - x) - sin(2*b - 2*d - x)*sin(c - x),
    x,
    domain=S.Reals,
)

However, after a long computation, it just returns the original expression:

ConditionSet(x, Eq(-sin(x)*sin(-c + 2*d + x) + sin(c - x)*sin(-2*b + 2*d + x), 0), Reals)

matlab

I also tried using Matlab for solving the equation:

syms b c d x

eq1 = sin(x)*sin(c - 2*d - x) - sin(2*b - 2*d - x)*sin(c - x) == 0;
eq2 = 0 < d;
eq3 = d < b;
eq4 = b < d+x;
eq5 = d+x < c;
result = solve([eq1 eq2 eq3 eq4 eq5], x, "Real",true, "ReturnConditions",true);

However, the result is something entirely different, involving a fourth order polynomial of which I can't make any sense:

>> result.x
 
ans =
 
2*pi*k + atan2(y, 1) - atan2(-y, 1)

>> result.conditions
 
ans =
 
d < b & b + atan2(-y, 1) < d + atan2(y, 1) + 2*pi*k & d + atan2(y, 1) + 2*pi*k < c + atan2(-y, 1) & cos(b)*sin(b) + 2*cos(b)^2*cos(d)*sin(d) ~= cos(d)*sin(d) + 2*cos(b)*cos(d)^2*sin(b) & in(k, 'integer') & ~in(c/pi, 'integer') & in(y, 'real') & 4*y^2*cos(c - 2*d) + 2*y^3*sin(c - 2*d) + sin(2*b - 2*d)*sin(c) + y^2*cos(c - 2*b + 2*d) + 3*y^2*cos(2*b + c - 2*d) + 2*y^3*sin(2*b + c - 2*d) + y^4*sin(2*b - 2*d)*sin(c) == 2*y*sin(2*b + c - 2*d) + 2*y*sin(c - 2*d) & 0 < d
a_guest
  • 34,165
  • 12
  • 64
  • 118
  • An explicit solution here is possible but it would be very complicated resulting from the four roots of the quartic formula and its unlikely that the inequalities could be used to do anything with it. Maybe take a step back and think about what it is you are actually trying to achieve: the explicit solution is unlikely to be the best way. – Oscar Benjamin Apr 27 '22 at 10:45
  • @OscarBenjamin Perhaps, but I'd like obtain a better understanding of the solution. Regarding the Matlab solution, I don't see at all how this fourth order polynomial enters the scene. I mean, I start with an expression involving only trigonometric functions, so why do I end up with a polynomial then? And there's still a discrepancy between the results from Matlab and Sympy. Sympy seems to say that it's impossible to find a solution (though I'm not sure if that is what the error means), but Matlab provides a solution. Only if the fourth order polynomial has no roots, this could be equivalent. – a_guest Apr 27 '22 at 11:11
  • The result from sympy does not imply that there is no solution. The 4th order polynomial comes from eliminating either sin or cos which are quadratic in the original equation. You can do the same using `resultant` in sympy. Either way an explicit form of the solution will be complicated: sufficiently complicated expressions are not really much use for understanding anything. – Oscar Benjamin Apr 27 '22 at 11:32
  • @OscarBenjamin I got the impression that Sympy's error implies that no solution exists and that would have been a valuable piece of information. Though you seem to be saying that the error is about something else. Could you please explain why `solve` raises that error and what it means for the specific example? – a_guest Apr 27 '22 at 12:42
  • The error message is just a bug. What it should say is that the symbolic inequality solver is unable to fully solve the inequalities. – Oscar Benjamin Apr 27 '22 at 12:44

1 Answers1

1

Here is how you can obtain something similar to the Matlab solution using SymPy (note that I redefined the symbols to make the output simpler):

In [19]: a, b, c, d, e, x, y = symbols('a, b, c, d, e, x, y', real=True)

In [20]: expr = sin(x)*sin(a - x) - sin(e - x)*sin(c - x)

In [21]: expr.expand(trig=True)
Out[21]: 
                                        2                                                          
sin(a)⋅sin(x)⋅cos(x) - sin(c)⋅sin(e)⋅cos (x) + sin(c)⋅sin(x)⋅cos(e)⋅cos(x) + sin(e)⋅sin(x)⋅cos(c)⋅c

           2                2                 
os(x) - sin (x)⋅cos(a) - sin (x)⋅cos(c)⋅cos(e)

In [22]: expr2 = resultant(expr.expand(trig=True), cos(x)**2 + sin(x)**2 - 1, cos(x))

In [23]: expr2.collect(sin(x))
Out[23]: 
⎛     2                                                             2       2         2       2    
⎝- sin (a) - 2⋅sin(a)⋅sin(c)⋅cos(e) - 2⋅sin(a)⋅sin(e)⋅cos(c) - 2⋅sin (c)⋅sin (e) - sin (c)⋅cos (e) 

                              2       2   ⎞    2      ⎛   2                                        
+ 2⋅sin(c)⋅sin(e)⋅cos(a) - sin (e)⋅cos (c)⎠⋅sin (x) + ⎝sin (a) + 2⋅sin(a)⋅sin(c)⋅cos(e) + 2⋅sin(a)⋅

                   2       2         2       2                                  2       2         2
sin(e)⋅cos(c) + sin (c)⋅sin (e) + sin (c)⋅cos (e) - 2⋅sin(c)⋅sin(e)⋅cos(a) + sin (e)⋅cos (c) + cos 

                                  2       2   ⎞    4         2       2   
(a) + 2⋅cos(a)⋅cos(c)⋅cos(e) + cos (c)⋅cos (e)⎠⋅sin (x) + sin (c)⋅sin (e)

So we see that what we have is a polynomial of degree 4 in sin(x). This one is a bit nicer to solve than the one that Matlab found because it only has even powers of sin(x). The explicit solutions for sin(x) can be found using roots e.g. the first is:

In [26]: r1, r2, r3, r4 = roots(expr2.subs(sin(x), y), y)

In [27]: r1
Out[27]: 
          _________________________________________________________________________________________
         ╱                         ________________________________________________________________
        ╱                         ╱                                                                
       ╱                         ╱                               2       2      (cos(2⋅a) + cos(2⋅c
      ╱                         ╱   - 8⋅(cos(-a + c + e) + 1)⋅sin (c)⋅sin (e) + ───────────────────
     ╱                        ╲╱                                                                   
-   ╱      - ──────────────────────────────────────────────────────────────────────────────────────
   ╱           ⎛   2                                                           2       2         2 
 ╲╱          2⋅⎝sin (a) + 2⋅sin(a)⋅sin(c)⋅cos(e) + 2⋅sin(a)⋅sin(e)⋅cos(c) + sin (c)⋅sin (e) + sin (

___________________________________________________________________________________________________
__________________________________________________________________________________________         
                                                                                        2          
) + cos(2⋅e) - 3⋅cos(-a + c + e) + cos(a - c + e) + cos(a + c - e) + cos(a + c + e) - 3)           
─────────────────────────────────────────────────────────────────────────────────────────          
                                  4                                                                
───────────────────────────────────────────────────────────────────────────────────────────────────
      2                                  2       2         2                                  2    
c)⋅cos (e) - 2⋅sin(c)⋅sin(e)⋅cos(a) + sin (e)⋅cos (c) + cos (a) + 2⋅cos(a)⋅cos(c)⋅cos(e) + cos (c)⋅

___________________________________________________________________________________________________
                                                                                                   
                                                                                                   
                                                                                                   
                                                                                                   
           -cos(2⋅a) - cos(2⋅c) - cos(2⋅e) + 3⋅cos(-a + c + e) - cos(a - c + e) - cos(a + c - e) - 
──────── + ────────────────────────────────────────────────────────────────────────────────────────
   2   ⎞                                            8⋅(cos(-a + c + e) + 1)                        
cos (e)⎠                                                                                           

___________________
                   
                   
                   
                   
cos(a + c + e) + 3 
──────────────────

Whether or not this is real or satisfies the inequalities will depend on the values of the symbolic parameters so I doubt that the inequalities can be used to choose amongst the 4 expressions for the roots. It is perhaps possible to answer different questions using the polynomial expr2 though.

Oscar Benjamin
  • 12,649
  • 1
  • 12
  • 14