Improve _solve_trig1 (handle rational & symbolic coefficients) #19507

gschintgen · 2020-06-07T13:37:58Z

References to other Issues or PRs

Fixes #18315
Closes #9616
Fixes #16870
Closes #17543 (adds test; bare except was already taken out by earlier commits; now handled by _solve_trig1)
Fixes #11218
Improves: #18427 (exception no longer raised since the equation is now solved by _solve_trig1; the bug in _solve_trig2 is hidden, not fixed. Test added)

Brief description of what is fixed or changed

This PR fixes and improves the main trig (and hyperbolic) solver in solveset.

The initial commit adds some pre-processing to _solve_trig1. The motivation
of this processing is to fix the change of variable exp(I*x) -> y
for equations where some coefficients are non-integers. The actual code
is an enhanced version of the one found in _solve_trig2.

Before:

    _solve_trig1(sin(x/2), x, Complexes)
    _solve_trig1(sin(pi*x), x, Complexes)
    _solve_trig1(sin(sqrt(2)*x), x, Complexes)
    _solve_trig1(sin(a*x), x, Complexes)

all result in ConditionSets due to a failing 'subs' and execution passes
on to _solve_trig2 which already has this pre-processing. (Which is
incomplete, so that uncaught exceptions are raised...)
Also: _solve_trig1(sin(5*x)+cos(10*x), x, Complexes) results in an equation
of degree 20(!) even though a simple change of variable could take this
down to 4.

After:
The equation types above are handled fine. If necessary a ConditionSet is returned. See the added tests for examples.

To sum up:

sometimes huge performance improvements with far nicer solution sets. (e.g. for solveset(sin(10*x)))
NEW: solveset handles symbolic coefficients in trig arguments (e.g. solveset(sin(a**2*x/pi), x)) and returns the consistency conditions as part of the solution set.

Other comments

Release Notes

solvers
- Improved solveset capabilities for solving trigonometric equations, notably rational and symbolic coefficients are now supported.

Fixes sympy#18315 (and others). This commit adds some pre-processing to _solve_trig1. The motivation of this processing is to fix the change of variable exp(I*x) -> y for equations where some coefficients are non-integers. The actual code is an enhanced version of the one found in _solve_trig2. Before: _solve_trig1(sin(x/2), x, Complexes) _solve_trig1(sin(pi*x), x, Complexes) _solve_trig1(sin(sqrt(2)*x), x, Complexes) _solve_trig1(sin(a*x), x, Complexes) all result in ConditionSets due to a failing 'subs' and execution passes on to _solve_trig2 which already has this pre-processing. (Which is incomplete, so that uncaught exceptions are raised...) Also: _solve_trig1(sin(5*x)+cos(10*x), x, Complexes) results in an equation of degree 20(!) even though a simple change of variable could take this down to 4. After: The equation types above are handled fine.

Previously some of the tests couldn't be solved by _solve_trig1 and ended up being solved by _solve_trig2. These are moved around a bit and a comment is added.

sympy-bot · 2020-06-07T13:37:59Z

✅

Hi, I am the SymPy bot (v160). I'm here to help you write a release notes entry. Please read the guide on how to write release notes.

Your release notes are in good order.

Here is what the release notes will look like:

solvers
- Improved solveset capabilities for solving trigonometric equations, notably rational and symbolic coefficients are now supported. (#19507 by @gschintgen)

This will be added to https://github.com/sympy/sympy/wiki/Release-Notes-for-1.7.

Note: This comment will be updated with the latest check if you edit the pull request. You need to reload the page to see it.

Click here to see the pull request description that was parsed.

#### References to other Issues or PRs
Fixes #18315
Closes #9616
Fixes #16870
Closes #17543 (adds test; bare except was already taken out by earlier commits; now handled by _solve_trig1)
Fixes #11218
Improves: #18427 (exception no longer raised since the equation is now solved by `_solve_trig1`; the bug in `_solve_trig2` is hidden, not fixed. Test added)

#### Brief description of what is fixed or changed
This PR fixes and improves the main trig (and hyperbolic) solver in `solveset`.

The initial commit adds some pre-processing to `_solve_trig1`. The motivation
of this processing is to fix the change of variable `exp(I*x) -> y`
 for equations where some coefficients are non-integers. The actual code
is an enhanced version of the one found in `_solve_trig2`.

Before:
```
    _solve_trig1(sin(x/2), x, Complexes)
    _solve_trig1(sin(pi*x), x, Complexes)
    _solve_trig1(sin(sqrt(2)*x), x, Complexes)
    _solve_trig1(sin(a*x), x, Complexes)
```
all result in ConditionSets due to a failing 'subs' and execution passes
on to _solve_trig2 which already has this pre-processing. (Which is
incomplete, so that uncaught exceptions are raised...)
Also: `_solve_trig1(sin(5*x)+cos(10*x), x, Complexes)` results in an equation
of degree 20(!) even though a simple change of variable could take this
down to 4.

After:
The equation types above are handled fine. If necessary a `ConditionSet` is returned. See the added tests for examples.

To sum up:  
- sometimes huge performance improvements with far nicer solution sets. (e.g. for `solveset(sin(10*x))`)  
- NEW: `solveset` handles symbolic coefficients in trig arguments (e.g. `solveset(sin(a**2*x/pi), x)`) and returns the consistency conditions as part of the solution set.
#### Other comments


#### Release Notes

<!-- Write the release notes for this release below. See
https://github.com/sympy/sympy/wiki/Writing-Release-Notes for more information
on how to write release notes. The bot will check your release notes
automatically to see if they are formatted correctly. -->

<!-- BEGIN RELEASE NOTES -->
* solvers
  * Improved `solveset` capabilities for solving trigonometric equations, notably rational and symbolic coefficients are now supported.
<!-- END RELEASE NOTES -->

Update

The release notes on the wiki have been updated.

gschintgen · 2020-06-07T13:40:39Z

(I stopped the tests, they'll fail anyway due to #19496)

gschintgen · 2020-06-07T14:55:16Z

sympy/solvers/solveset.py

+
+    # lcm() and gcd() require more than one argument
+    if len(numerators) > 1:
+        mu = lcm(*denominators)/gcd(*numerators)


Are there any pitfalls lurking here that should be checked for? (e.g. coefficients involving I and their handling by the Polys module.)

The gcd function has a strange 2-arg vs 1-arg iterable syntax. It's best not to use star-unpacking:

In [16]: gcd(S(6), S(2), S(3)) Out[16]: 2 In [17]: gcd([S(6), S(2), S(3)]) Out[17]: 1

Ouch, thanks, wasn't aware of that!

gschintgen · 2020-06-07T14:59:26Z

sympy/solvers/solveset.py

-    if g.has(symbol) or h.has(symbol):
-        return ConditionSet(symbol, Eq(f, 0), domain)
+    if g.has(x) or h.has(x):
+        raise NotImplementedError("change of variable not possible")


returning a ConditionSet would be more in line with the rest of solveset. But that entails further modifications in the upper layers: e.g. add an is_unsolved_conditionset() function to distinguish between ConditionSets as failure condition and ConditionSets as perfectly complete solution sets including validity conditions.

Another possibility is returning a special value. If the intention is always for the exception to be caught in a specific place then I think that a specific local exception class is better as there are many places around the code base that catch or raise NotImplementedError:

$ git grep 'raise NotImplementedError' | wc -l 728 $ git grep 'except NotImplementedError' | wc -l 127

That makes it easy for an exception to end up being unintentionally caught in the wrong place.

Thanks for the advice. Maybe I'll just add is_unsolved_conditionset after all. I guess there are other places in solveset that are not 100% prepared to have ConditionSets returned to them that don't mean "I couldn't solve this". So something like that will be needed anyway and passing along Sets/ConditionSets is the established pattern in solveset.py.

oscarbenjamin · 2020-06-07T19:38:15Z

sympy/solvers/solveset.py

+        try:
+            poly_ar = Poly(ar, symbol)
+        except ValueError:
+            raise ValueError("give up, we can't solve if this is not a polynomial in x")


When does Poly raise ValueError?

If these ValueErrors are intended to be caught above then it would be better to make a special exception class for it like:

class _SolveTrig1Error(Exception): """Raised when solvetrig1 heuristics do not apply""" pass

If the exception class is local to these helper functions then there is no need for it to be part of any hierarchy or public API.

When does Poly raise ValueError?

Hmm, good question. I copied it over from _solve_trig2 and thought it was to harden the code a bit if for some reason a non-polynomial ends up as an argument to a trig function. (This shouldn't happen if called via solveset, but maybe it was intended to make it more self-contained and guard against user error.)
But now I notice that it actually raises PolynomialError for obviously non-polynomial arguments and not ValueError. I'll grep through Poly code and probably change it to PolynomialError.

I introduced _SolveTrig1Error as suggested.

oscarbenjamin · 2020-06-07T20:08:40Z

sympy/solvers/solveset.py

        # avoid spurious intersections with C in solution set
        if domain is S.Complexes:
-            return result
+            return ConditionSet(symbol, cond, result)


I've often wondered looking at these lines what are the spurious intersections with C? It would be better to fix is_subset rather than work around its limitations here (assuming that the fix is straight-forward).

I think that with a blanket .intersect(domain) i.e. .intersect(S.Complexes) there would have to be substantial improvements to is_subset to automatically simplify them away. result could be a rather arbitrary Set expression involving ImageSets, ConditionSets, Unions thereof etc. all with potentially complex algebraic expressions. Trying to enforce mathematical correctness (protect against possible zoo) by slapping on .intersect(Complexes) doesn't seem the right approach to me. It only leads to almost any result having an Intersection(..., Complexes) that I may have to peel off.
Instead, singularities should be avoided/handled from the get go.
Here's an unrelated but instructive example:

In [1]: solveset(exp(x) - y, x, Reals) Out[1]: ℝ ∩ {log(y)}

That's what I consider a "slapped on" .intersect(Reals). Yes, it makes the result mathematically correct, but it gives no useful information at all on acceptable values of y. IMO that call should return

ConditionSet(x, y > 0, FiniteSet(log(y))) Out[2]: {x | x ∊ {log(y)} ∧ (y > 0)}

or, because non-reals raise an exception when tested for > 0 instead of returning False, maybe even:

In [14]: ConditionSet(x, And(Contains(y, Reals), y > 0), FiniteSet(log(y))) Out[14]: {x | x ∊ {log(y)} ∧ ((y ∈ ℝ) ∧ y > 0)} In [15]: _.subs(y, I) Out[15]: ∅

Unfortunately that's quite involved/ugly. I'd prefer if I > 0 did not raise and simply evaluated to False.

I'd prefer if I > 0 did not raise and simply evaluated to False.

You aren't the only one who would prefer that. This is discussed in many places e.g. #18778 (comment)

Unfortunately the guard And(Contains(y, Reals), y>0) doesn't work because there's no real guarantee that subs will sub into the Contains before the Gt.

Trying to enforce mathematical correctness (protect against possible zoo) by slapping on .intersect(Complexes) doesn't seem the right approach to me.

I agree. It shouldn't be hard though in the cases where there is no possible zoo to determine if the result is a subset of C. If the unevaluated intersection is there because it's needed to exclude a genuine zoo then the intersection is needed for correctness.

Using ConditionSet for consistency conditions doesn't seem right to me. I think of ConditionSet is filtering from a set so a condition that doesn't involve the symbols in the expression feels like an abuse somehow.

I would prefer to express the result as (something like):

In [27]: Piecewise(({log(y)}, y>0), (EmptySet, True)) Out[27]: ⎧{log(y)} for y > 0 ⎨ ⎩ ∅ otherwise

For inequalities raising see also https://groups.google.com/forum/#!topic/sympy/W0iKfAGsc2E

I'll look into it. At least this shouldn't be too difficult to fix:

The example that works there works because of this (which I'd like to remove):

sympy/sympy/sets/sets.py

Lines 386 to 391 in eee759c

# Fall back on computing the intersection

# XXX: We shouldn't do this. A query like this should be handled

# without evaluating new Set objects. It should be the other way round

# so that the intersect method uses is_subset for evaluation.

if self.intersect(other) == self:

return True

The example that works there works because of this (which I'd like to remove):

Yes, I've also now tracked it down. It's because there's special handling of S.Reals in

@dispatch(ImageSet, Set) # type: ignore # noqa:F811 def intersection_sets(self, other)

It's even using solveset_real on the imaginary part, so it's more sophisticated than what I expected. For the present case, intersection with C, maybe all that is needed is something like _is_finite_with_finite_vars

sympy/sympy/solvers/solveset.py

Lines 440 to 446 in eee759c

def _is_finite_with_finite_vars(f, domain=S.Complexes):

"""

Return True if the given expression is finite. For symbols that

don't assign a value for `complex` and/or `real`, the domain will

be used to assign a value; symbols that don't assign a value

for `finite` will be made finite. All other assumptions are

left unmodified.

I'll do some experiments later.

Checking imgset.lamda.expr.is_complex should be enough I think. (I did a local proof-of-concept.) But it seems this is not the only deficiency. I think some further checking of .intersect(Reals), .intersect(Complexes), .is_subset(Reals), is_subset(Complexes) is necessary, particularly for ImageSets. For now I'd like to just add an additional comment concerning the "spurious" intersections and leave the is_subset details for another PR.

newly filed issues: #19513, #19514

I think _is_finite_with_finite_vars sounds like the best solution. Checking imgset.lamda.expr.is_complex will depend on the assumptions of bound symbols:

In [238]: x = Symbol('x') In [239]: print(ImageSet(Lambda(x, x**2), Reals).lamda.expr.is_real) None In [240]: x = Symbol('x', real=True) In [241]: print(ImageSet(Lambda(x, x**2), Reals).lamda.expr.is_real) True

It doesn't need to be done in this PR though.

It doesn't need to be done in this PR though.

Ok! I added some comments in the code.

sympy/solvers/tests/test_solveset.py

oscarbenjamin · 2020-06-07T20:21:59Z

sympy/solvers/tests/test_solveset.py

+        ImageSet(Lambda(n, n*pi/2 + pi/16), S.Integers)))
+
+    # This is the only remaining solveset test that actually ends up being solved
+    # by _solve_trig2(). All others are handled by the improved _solve_trig1.
    assert dumeq(solveset_real(2*cos(x)*cos(2*x) - 1, x),


I haven't looked into the general architecture of solveset. Why are there two solvetrig functions?

The comment above makes me wonder if solvetrig2 is actually needed or if it's just compensating for some obvious deficiencies in solvetrig1.

You're right about that. Here's what I wrote in #18315 (comment):

I've also been wondering if this absence of pre-processing is not the actual reason _solve_trig2 is needed at all. AFAIR I stumbled over the PR introducing _solve_trig2 and it was added as an alternative to _solve_trig1 for cases where it could not solve the equation, but maybe those deficiencies were only a symptom of the present issue.

But there seem to be examples where trig2 is still needed.
_solve_trig1 does rewrite to exp, while _solve_trig2 will rewrite to tan and use half-angle formulas. Both do a change of variable to obtain a rational equation. It seems that in some cases the exp-induced polynomials result in CRootOfs while the polynomials resulting from tan-rewrite can be solved with radicals. (Further investigation is necessary to confirm.)

Anyway, on my todo list (not for this PR) is: document both functions, analyze _solve_trig2 necessity, port the improvements to _solve_trig1 here over to _solve_trig2, add support for hyperbolic equations (rewrite to tanh...)

oscarbenjamin · 2020-06-07T20:29:45Z

Looking through the tests here I can certainly see why we need to be able to simplify unions of image sets....

codecov · 2020-06-07T21:38:11Z

Codecov Report

Merging #19507 into master will decrease coverage by 0.002%.
The diff coverage is 81.395%.

@@              Coverage Diff              @@
##            master    #19507       +/-   ##
=============================================
- Coverage   75.690%   75.687%   -0.003%     
=============================================
  Files          653       653               
  Lines       169873    169895       +22     
  Branches     40059     40061        +2     
=============================================
+ Hits        128577    128590       +13     
- Misses       35686     35693        +7     
- Partials      5610      5612        +2

add fix for ConditionSet.dummy_eq()

Issues: sympy#11218, sympy#18427, sympy#17543

oscarbenjamin · 2020-06-08T20:58:35Z

Tests haven't passed yet but this looks good to me

oscarbenjamin · 2020-06-09T22:41:17Z

The coverage report indicates a few lines that don't get hit during the tests. For example what sort of input would lead to PolynomialError? Is it possible to add a test that covers that case?

This looks good to me so I'm happy to merge it if you're done. I often find that adding last tests to try and ensure that everything is covered leads to me discovering possible improvements in my work (but that doesn't apply in all situations).

gschintgen · 2020-06-10T07:14:40Z

The coverage report indicates a few lines that don't get hit during the tests. For example what sort of input would lead to PolynomialError? Is it possible to add a test that covers that case?

Thanks for pointing that out. I think that there's not much to do about those missed lines though: Some conditionals are for hardening against erroneous input and making the code more explicit. (E.g. PolynomialError or degree error). Those should be hard to actually hit (if at all), if execution ends up there via a call to solveset. At the upper levels, there's a call to _is_function_class_equation in order to dispatch to the trig solver. Equations that violate _solve_trig1 assumptions should probably not end up in _solve_trig1 at all, at least with the current implementation of _is_function_class_equation. (One way would be add a test directly targeting _solve_trig1, but I'm not sure that's necessary.)

In other cases, it would have been easy to hit the condition with the previous code, but I'm not so sure now. ("change of variable not possible")

I'd be grateful if it could be merged as is. _solve_trig2 also has a few rough edges, so I intend to revisit this part of the codebase anyway. (And this week will be rather busy.)
Thanks for your thorough review!

oscarbenjamin · 2020-06-10T10:01:33Z

If we think that it shouldn't be possible for PolynomialError to be raised because of checking done beforehand then I think it's best not to catch the exception. The reason is that if a PolynomialError can only indicate a bug somewhere then catching the exception would just hide or obscure that bug. On the other hand if there are reasonable inputs that can lead to the exception then it should be possible to add a test for them.

All the same this is a good improvement so let's merge it so the work can continue.

Thanks!

gschintgen added 7 commits June 6, 2020 12:07

solveset: some tests now handled by _solve_trig1

7dcf014

Previously some of the tests couldn't be solved by _solve_trig1 and ended up being solved by _solve_trig2. These are moved around a bit and a comment is added.

solveset tests: improve coverage of trig/hyperbolic equations

70df553

solveset tests: add test for issue sympy#9616

e0426fc

solveset tests (trig): add test for sympy#16870

9efb5a7

_solve_trig1: handle symbolic coefficients

8671f3c

_solve_trig tests: add tests for symbolic coefficients

7087e34

This was referenced Jun 7, 2020

bug in _solve_trig1 #18315

Closed

solveset(sin(pi*x), x, domain=S.Reals) doesn't give the solution #11218

Closed

solveset complicates and fails to solve sin(x/180*pi) - 1/2 #16870

Closed

solveset(sin(pi*x), x) raises TypeError in _solve_trig2 #18427

Closed

gschintgen commented Jun 7, 2020

View reviewed changes

oscarbenjamin closed this Jun 7, 2020

oscarbenjamin reopened this Jun 7, 2020

oscarbenjamin reviewed Jun 7, 2020

View reviewed changes

sympy/solvers/tests/test_solveset.py Outdated Show resolved Hide resolved

oscarbenjamin reviewed Jun 7, 2020

View reviewed changes

gschintgen added 4 commits June 8, 2020 12:24

_solve_trig1: apply some review suggestions

a99aa64

Merge remote-tracking branch 'upstream/master' into fix-solve-trig1

81954b6

add fix for ConditionSet.dummy_eq()

_solve_trig1: further review suggestions

c50a783

solveset (solve_trig): add tests for issues

3c8b95b

Issues: sympy#11218, sympy#18427, sympy#17543

gschintgen marked this pull request as ready for review June 8, 2020 19:33

gschintgen added the solvers.solveset label Jun 8, 2020

gschintgen mentioned this pull request Jun 8, 2020

solveset(I*cot(8*x - 8*E), x) inconsistent results based on keyboard interrupt #17543

Closed

oscarbenjamin merged commit d9bc559 into sympy:master Jun 10, 2020

gschintgen deleted the fix-solve-trig1 branch June 13, 2020 19:07

gschintgen mentioned this pull request Jan 17, 2021

Solveset returns empty set for cos(2x) = 1/2 #20798

Closed

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Improve _solve_trig1 (handle rational & symbolic coefficients) #19507

Improve _solve_trig1 (handle rational & symbolic coefficients) #19507

gschintgen commented Jun 7, 2020 •

edited

sympy-bot commented Jun 7, 2020 •

edited

gschintgen commented Jun 7, 2020

gschintgen Jun 7, 2020

oscarbenjamin Jun 7, 2020

gschintgen Jun 7, 2020

gschintgen Jun 7, 2020

oscarbenjamin Jun 7, 2020

gschintgen Jun 7, 2020

oscarbenjamin Jun 7, 2020

gschintgen Jun 7, 2020

gschintgen Jun 8, 2020

oscarbenjamin Jun 7, 2020

gschintgen Jun 7, 2020

oscarbenjamin Jun 7, 2020

oscarbenjamin Jun 7, 2020

oscarbenjamin Jun 7, 2020

oscarbenjamin Jun 8, 2020

gschintgen Jun 8, 2020 •

edited

gschintgen Jun 8, 2020

oscarbenjamin Jun 8, 2020

gschintgen Jun 8, 2020

oscarbenjamin Jun 7, 2020

gschintgen Jun 7, 2020

oscarbenjamin commented Jun 7, 2020

codecov bot commented Jun 7, 2020 •

edited

oscarbenjamin commented Jun 8, 2020

oscarbenjamin commented Jun 9, 2020

gschintgen commented Jun 10, 2020

oscarbenjamin commented Jun 10, 2020

	# Fall back on computing the intersection
	# XXX: We shouldn't do this. A query like this should be handled
	# without evaluating new Set objects. It should be the other way round
	# so that the intersect method uses is_subset for evaluation.
	if self.intersect(other) == self:
	return True

	def _is_finite_with_finite_vars(f, domain=S.Complexes):
	"""
	Return True if the given expression is finite. For symbols that
	don't assign a value for `complex` and/or `real`, the domain will
	be used to assign a value; symbols that don't assign a value
	for `finite` will be made finite. All other assumptions are
	left unmodified.

Improve _solve_trig1 (handle rational & symbolic coefficients) #19507

Improve _solve_trig1 (handle rational & symbolic coefficients) #19507

Conversation

gschintgen commented Jun 7, 2020 • edited

References to other Issues or PRs

Brief description of what is fixed or changed

Other comments

Release Notes

sympy-bot commented Jun 7, 2020 • edited

gschintgen commented Jun 7, 2020

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

gschintgen Jun 8, 2020 • edited

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

oscarbenjamin commented Jun 7, 2020

codecov bot commented Jun 7, 2020 • edited

Codecov Report

oscarbenjamin commented Jun 8, 2020

oscarbenjamin commented Jun 9, 2020

gschintgen commented Jun 10, 2020

oscarbenjamin commented Jun 10, 2020

gschintgen commented Jun 7, 2020 •

edited

sympy-bot commented Jun 7, 2020 •

edited

gschintgen Jun 8, 2020 •

edited

codecov bot commented Jun 7, 2020 •

edited