forked from sympy/sympy
-
Notifications
You must be signed in to change notification settings - Fork 0
Finding roots of polynomials
certik edited this page Feb 8, 2011
·
1 revision
Let's try <math>x^3 + 2x^2 + 8</math>:
In [1]: solve(x**3+2*x**2+8, x)
Out[1]:
⎡ ⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽ 2/3
⎢ 3 ╱ ⎽⎽⎽⎽ 2/3 3 ⎽⎽⎽ ⎛ ⎽⎽⎽⎽⎞ 3
⎢- 18*╲╱ 1044 + 108*╲╱ 93 - 36*3 - 3*╲╱ 3 *�1044 + 108*╲╱ 93 ⎠18*╲╱
⎢──────────────────────────────────────────────────────────────────────, ─────
⎢ ⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽
⎢ 3 ╱ ⎽⎽⎽⎽
⎣ 27*╲╱ 1044 + 108*╲╱ 93
⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽
⎽⎽⎽ ⎽⎽⎽⎽ 3 ⎽⎽⎽ ⎽⎽⎽⎽ 3 ⎽⎽⎽ 2/3 3 ╱ ⎽⎽⎽⎽
3 *╲╱ 93 - 54*ⅈ*╲╱ 3 *╲╱ 31 + 174*╲╱ 3 + 2*3 *╲╱ 1044 + 108*╲╱ 93 + 6
──────────────────────────────────────────────────────────────────────────────
⎛ ⎽⎽
3*�1044 + 108*╲╱ 9
⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽ 2/3
6 ⎽⎽⎽ 3 ╱ ⎽⎽⎽⎽ ⎛ ⎽⎽⎽⎽⎞ 5/6 3
- ⅈ*╲╱ 3 *╲╱ 1044 + 108*╲╱ 93 - 2*�1044 + 108*╲╱ 93 ⎠- 174*ⅈ*3 18*╲╱
- ⅈ*╲╱ 3 *╲╱ 1044 + 108*╲╱ 93 - 2*�1044 + 108*╲╱ 93 ⎠+ 174*ⅈ*3 ⎥
Too messy? Let's latex it:
>>> latex(solve(x**3+2*x**2+8, x))
'$\\begin{bmatrix}\\frac{- 18 \\left(1044 + 108 \\sqrt{93}\\right)^{\\frac{1}{3}} - 36 {3}^{\\frac{2}{3}} - 3 {3}^{\\frac{1}{3}} \\left(1044 + 108 \\sqrt{93}\\right)^{\\frac{2}{3}}}{27 \\left(1044 + 108 \\sqrt{93}\\right)^{\\frac{1}{3}}}, & \\frac{- 174 \\mathbf{\\imath} {3}^{\\frac{5}{6}} - 2 \\left(1044 + 108 \\sqrt{93}\\right)^{\\frac{2}{3}} - 54 \\mathbf{\\imath} {3}^{\\frac{1}{3}} \\sqrt{31} + 2 {3}^{\\frac{2}{3}} \\left(1044 + 108 \\sqrt{93}\\right)^{\\frac{1}{3}} + 18 {3}^{\\frac{1}{3}} \\sqrt{93} + 6 \\mathbf{\\imath} {3}^{\\frac{1}{6}} \\left(1044 + 108 \\sqrt{93}\\right)^{\\frac{1}{3}} + 174 {3}^{\\frac{1}{3}}}{3 \\left(1044 + 108 \\sqrt{93}\\right)^{\\frac{2}{3}}}, & \\frac{- 2 \\left(1044 + 108 \\sqrt{93}\\right)^{\\frac{2}{3}} + 2 {3}^{\\frac{2}{3}} \\left(1044 + 108 \\sqrt{93}\\right)^{\\frac{1}{3}} + 54 \\mathbf{\\imath} {3}^{\\frac{1}{3}} \\sqrt{31} + 18 {3}^{\\frac{1}{3}} \\sqrt{93} - 6 \\mathbf{\\imath} {3}^{\\frac{1}{6}} \\left(1044 + 108 \\sqrt{93}\\right)^{\\frac{1}{3}} + 174 {3}^{\\frac{1}{3}} + 174 \\mathbf{\\imath} {3}^{\\frac{5}{6}}}{3 \\left(1044 + 108 \\sqrt{93}\\right)^{\\frac{2}{3}}}\\end{bmatrix}$'
In [1]: p = x**3+2*x**2+8
In [2]: r = solve(p, x)
In [3]: p.subs(x, r[0])
Out[3]:
⎛ ⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽ 2/3⎞
⎜ 3 ╱ ⎽⎽⎽⎽ 2/3 3 ⎽⎽⎽ ⎛ ⎽⎽⎽⎽⎞ ⎟
2*âŽ�- 18*╲╱ 1044 + 108*╲╱ 93 - 36*3 - 3*╲╱ 3 *âŽ�1044 + 108*╲╱ 93 ⎠âŽ
8 + ──────────────────────────────────────────────────────────────────────────
2/3
⎛ ⎽⎽⎽⎽⎞
729*âŽ�1044 + 108*╲╱ 93 âŽ
2 3
⎛ ⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽ 2/3⎞
⎜ 3 ╱ ⎽⎽⎽⎽ 2/3 3 ⎽⎽⎽ ⎛ ⎽⎽⎽⎽⎞ ⎟
âŽ�- 18*╲╱ 1044 + 108*╲╱ 93 - 36*3 - 3*╲╱ 3 *âŽ�1044 + 108*╲╱ 93 ⎠âŽ
─ + ─────────────────────────────────────────────────────────────────────────
⎛ ⎽⎽⎽⎽⎞
19683*âŽ�1044 + 108*╲╱ 93 âŽ
Hm, is this zero? Let's use a "brute force":
In [4]: p.subs(x, r[0]).evalf()
Out[4]: 2.1316282072803e-14
In [5]: p.subs(x, r[1]).evalf()
Out[5]: -1.77635683940025e-15 - 3.10862446895044e-15*â…ˆ
In [6]: p.subs(x, r[2]).evalf()
Out[6]: -3.5527136788005e-15 + 2.66453525910038e-15*â…ˆ
In [1]: solve(x**3 + x**2 - x + 1, x)[1]
Out[1]:
⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽
3 ⎽⎽⎽ ⎽⎽⎽⎽ 3 ⎽⎽⎽ ⎽⎽⎽⎽ 3 ⎽⎽⎽ 2/3 3 ╱ ⎽⎽⎽⎽
9*╲╱ 3 *╲╱ 33 + 27*ⅈ*╲╱ 3 *╲╱ 11 + 57*╲╱ 3 + 4*3 *╲╱ 171 + 27*╲╱ 33 -
──────────────────────────────────────────────────────────────────────────────
⎛ ⎽⎽
6*�171 + 27*╲╱ 3
⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽ 2/3
6 ⎽⎽⎽ 3 ╱ ⎽⎽⎽⎽ ⎛ ⎽⎽⎽⎽⎞ 5/6
12*ⅈ*╲╱ 3 *╲╱ 171 + 27*╲╱ 33 - 2*�171 + 27*╲╱ 33 ⎠+ 57*ⅈ*3
────────────────────────────────────────────────────────────────────
2/3
⎽⎽⎞
3 âŽ
In [1]: p = x**3+x**2+x-1
In [2]: a = solve(p, x)
In [3]: p.subs(x, a[0])
Out[3]:
⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽
3 ⎽⎽⎽ ⎽⎽⎽⎽ 3 ⎽⎽⎽ ⎽⎽⎽⎽ 3 ⎽⎽⎽ 2/3 3 ╱ ⎽⎽
9*╲╱ 3 *╲╱ 33 - 27*ⅈ*╲╱ 3 *╲╱ 11 - 51*╲╱ 3 - 2*3 *╲╱ -153 + 27*╲╱ 3
-1 + ─────────────────────────────────────────────────────────────────────────
⎛
6*�-153 + 2
⎽⎽⎽ ⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽ 2/3
⎽⎽ 6 ⎽⎽⎽ 3 ╱ ⎽⎽⎽⎽ ⎛ ⎽⎽⎽⎽⎞ 5/6
3 - 6*ⅈ*╲╱ 3 *╲╱ -153 + 27*╲╱ 33 - 2*�-153 + 27*╲╱ 33 ⎠+ 51*ⅈ*3
─────────────────────────────────────────────────────────────────────────── +
2/3
⎽⎽⎽⎽⎞
7*╲╱ 33 âŽ
⎛ ⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽
⎜ 3 ⎽⎽⎽ ⎽⎽⎽⎽ 3 ⎽⎽⎽ ⎽⎽⎽⎽ 3 ⎽⎽⎽ 2/3 3 ╱ ⎽⎽⎽⎽
�9*╲╱ 3 *╲╱ 33 - 27*ⅈ*╲╱ 3 *╲╱ 11 - 51*╲╱ 3 - 2*3 *╲╱ -153 + 27*╲╱ 33
──────────────────────────────────────────────────────────────────────────────
⎛
36*�-153 + 27*
2
⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽ 2/3 ⎞ ⎛
6 ⎽⎽⎽ 3 ╱ ⎽⎽⎽⎽ ⎛ ⎽⎽⎽⎽⎞ 5/6⎟ ⎜
- 6*ⅈ*╲╱ 3 *╲╱ -153 + 27*╲╱ 33 - 2*�-153 + 27*╲╱ 33 ⎠+ 51*ⅈ*3 ⎠�9
───────────────────────────────────────────────────────────────────────── + ──
4/3
⎽⎽⎽⎽⎞
╲╱ 33 âŽ
⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽
3 ⎽⎽⎽ ⎽⎽⎽⎽ 3 ⎽⎽⎽ ⎽⎽⎽⎽ 3 ⎽⎽⎽ 2/3 3 ╱ ⎽⎽⎽⎽
- ╲╱ 3 *╲╱ 33 - 27*ⅈ*╲╱ 3 *╲╱ 11 - 51*╲╱ 3 - 2*3 *╲╱ -153 + 27*╲╱ 33 -