# How to find all roots with solve?

I have the following quantity:

c_p = 1 - (4*sin(th)^2+2*5*sin(th)/pi+(5/(2*pi))^2)


and I am trying to find its roots (the solutions for c_p==0).

When plotting c_p as a function of th between $-\pi$ and $\pi$, I can see the curve crosses the x-axis at four positions. However, solve(c_p==0, th) is only giving two roots: $$\left[theta = -\arcsin\left(\frac{5}{4 \ \pi} + \frac{1}{2}\right),\quad theta = -\arcsin\left(\frac{5}{4 \ \pi} - \frac{1}{2}\right)\right].$$ It appears that solve can only find the roots that are in the domain of the $\arcsin$ function, i.e in the interval $[-\dfrac{\pi}{2},\dfrac{\pi}{2}]$. Is there a way to get the other two roots?

edit retag close merge delete

Sort by » oldest newest most voted

Sympy offers us a solution, but sage is not (yet) ready to parse it totally. But with a little Sympy's doc reading and a little manual unscrewing, one can get it:

sage: [u.args.args._sage_() for u in c_p.solve(theta, algorithm="sympy")]
[2*pi*n + arctan((2*pi - 5)/(pi*sqrt(20/pi - 25/pi^2 + 12))),
pi + 2*pi*n - arctan((2*pi - 5)/(pi*sqrt(20/pi - 25/pi^2 + 12))),
-pi + 2*pi*n + arctan(1/2*(2*pi + 5)/(pi*sqrt(-5/pi - 25/4/pi^2 + 3))),
2*pi + 2*pi*n - arctan(1/2*(2*pi + 5)/(pi*sqrt(-5/pi - 25/4/pi^2 + 3)))]


A bit of fooling around with latex and string splicing allows us to get the $\LaTeX$ form:

$$2 \pi n + \arctan\left(\frac{2 \pi - 5}{\pi \sqrt{\frac{20}{\pi} - \frac{25}{\pi^{2}} + 12}}\right)$$ $$\pi + 2 \pi n - \arctan\left(\frac{2 \pi - 5}{\pi \sqrt{\frac{20}{\pi} - \frac{25}{\pi^{2}} + 12}}\right)$$ $$-\pi + 2 \pi n + \arctan\left(\frac{2 \pi + 5}{2 \pi \sqrt{-\frac{5}{\pi} - \frac{25}{4 \pi^{2}} + 3}}\right)$$ $$2 \pi + 2 \pi n - \arctan\left(\frac{2 \pi + 5}{2 \pi \sqrt{-\frac{5}{\pi} - \frac{25}{4 \pi^{2}} + 3}}\right)$$

Note: Maxima could also get it (via to_poly_solve), but gives pompous pedantic nonsense involving logs of complex expressions, probably due to the domain:complex default...

Note 2: Substituting these values in c_p gets rid of the trig functions easily (via trig_reduce() after declaring n integer...). However, the resulting expressions involve sums ad products of scalar and square roots values ; proving them to be 0 is another possibly nontrivial problem... It can be "checked" numerically, and formally proved for the last two:

sage: L2=[u.args.args._sage_() for u in c_p.solve(theta, algorithm="sympy")]
sage: [c_p.subs(theta==u).trig_reduce().n() for u in L2]
[2.22044604925031e-16,
2.22044604925031e-16,
2.22044604925031e-16,
2.22044604925031e-16]

sage: [bool(c_p.subs(theta==u).trig_reduce()==0) for u in L2]
[False, False, True, True]


EDIT : A few notes for the terminally curious:

fricas and giac give the same results as maxima (i. e. Sage).

I have been unable to obtain no-nonsensical results from to_poly_solve.

One can also screw around Mathematica's answers to get another form of the results (closer to what Sage/Maxima gives) :

sage: var("k", domain="integer")
k
sage: LM=[u.sage(locals={"C":lambda x:k}) for u in mathematica.Solve(c_p
....: ==0, theta)]
sage: LM
[pi + 2*pi*k + arcsin(1/4*(2*pi + 5)/pi),
2*pi*k - arcsin(1/4*(2*pi + 5)/pi),
pi + 2*pi*k - arcsin(1/4*(2*pi - 5)/pi),
2*pi*k + arcsin(1/4*(2*pi - 5)/pi)]


which is:

$$\pi + 2 \pi k + \arcsin\left(\frac{2 \pi + 5}{4 \pi}\right)$$ $$2 \pi k - \arcsin\left(\frac{2 \pi + 5}{4 \pi}\right)$$ $$\pi + 2 \pi k - \arcsin\left(\frac{2 \pi - 5}{4 \pi}\right)$$ $$2 \pi k + \arcsin\left(\frac{2 \pi - 5}{4 \pi}\right)$$

which turns out to be somewhat easier to check/prove:

sage: [c_p.subs(theta==u).trig_reduce().expand() for u in LM]
[0, 0, 0, 0]

more

On what subject ? Getting all Sympy results in Sage u=is a large work, already undertaken notably by someone whose trac login is rws and whose name escapes me presently (getting old isn't fun...). The present occurrence is but a special case...

Do you suggest a "sympy usability" metaticket ?

The Mathematica interface has other problems, most notably the lack of correct Mathematica<-->Sage conversions. E. g., try to use arctan in a Mathematica call...

sage: mathematica.Integrate(arctan(x),x)
Integrate[Arctan[x], x]


This is true in the reverse direction also 'Sage complains of the lack of a function arcTan...)

Since Mathematica isn't free, I doubt that polishing its interface can be considered high-priority. Sympy's can, however...