FAQ
Just in case anybody was nervous about whether Python arithmetic expressions
have been adequately tested, the following works. Cool, huh? (Composed using
Maple 6 via Fortran codegeneration and a text editor)

from math import sqrt
def cmplx (x, y):
return x + 1.0j * y

def cubic_roots (a, b, c, d):
"Roots of a*x**3+b*x**2+c*x+d"
t0 = 1/a*(36*c*b*a-108*d*a**2-8*b**3+12*sqrt(3.0)*sqrt(4*c**3*a-c\
**2*b**2-18*c*b*a*d+27*d**2*a**2+4*d*b**3)*a)**(1.0/3.0)/6-2.0/\
3.0*(3*c*a-b**2)/a/(36*c*b*a-108*d*a**2-8*b**3+12*sqrt(3.0)*sqrt\

(4*c**3*a-c**2*b**2-18*c*b*a*d+27*d**2*a**2+4*d*b**3)*a)**(1.0/3.0)-b/a/3

t1 = -1/a*(36*c*b*a-108*d*a**2-8*b**3+12*sqrt(3.0)*sqrt(4*c**3*a-\
c**2*b**2-18*c*b*a*d+27*d**2*a**2+4*d*b**3)*a)**(1.0/3.0)/12+(3*\
c*a-b**2)/a/(36*c*b*a-108*d*a**2-8*b**3+12*sqrt(3.0)*sqrt(4*c**3*\
a-c**2*b**2-18*c*b*a*d+27*d**2*a**2+4*d*b**3)*a)**(1.0/3.0)/3-b/\
a/3+cmplx(0.0,1.0/2.0)*sqrt(3.0)*(1/a*(36*c*b*a-108*d*a**2-8*b\
**3+12*sqrt(3.0)*sqrt(4*c**3*a-c**2*b**2-18*c*b*a*d+27*d**2*a**2+\
4*d*b**3)*a)**(1.0/3.0)/6+2.0/3.0*(3*c*a-b**2)/a/(36*c*b*a-108\
*d*a**2-8*b**3+12*sqrt(3.0)*sqrt(4*c**3*a-c**2*b**2-18*c*b*a*d+27\
*d**2*a**2+4*d*b**3)*a)**(1.0/3.0))

t2 = -1/a*(36*c*b*a-108*d*a**2-8*b**3+12*sqrt(3.0)*sqrt(4*c**3*a-\
c**2*b**2-18*c*b*a*d+27*d**2*a**2+4*d*b**3)*a)**(1.0/3.0)/12+(3*\
c*a-b**2)/a/(36*c*b*a-108*d*a**2-8*b**3+12*sqrt(3.0)*sqrt(4*c**3*\
a-c**2*b**2-18*c*b*a*d+27*d**2*a**2+4*d*b**3)*a)**(1.0/3.0)/3-b/\
a/3+cmplx(0.0,-1.0/2.0)*sqrt(3.0)*(1/a*(36*c*b*a-108*d*a**2-8*\
b**3+12*sqrt(3.0)*sqrt(4*c**3*a-c**2*b**2-18*c*b*a*d+27*d**2*a**2\
+4*d*b**3)*a)**(1.0/3.0)/6+2.0/3.0*(3*c*a-b**2)/a/(36*c*b*a-108\
*d*a**2-8*b**3+12*sqrt(3.0)*sqrt(4*c**3*a-c**2*b**2-18*c*b*a*d+27\
*d**2*a**2+4*d*b**3)*a)**(1.0/3.0))

return [t0, t1, t2]

if __name__ == "__main__":
print cubic_roots(1.0, -1.0, 1.0, -1.0)

## Related Discussions

 view thread | post posts ‹ prev | 1 of 1 | next ›
Discussion Overview
 group python-list categories python posted Apr 13, '00 at 10:02p active Apr 13, '00 at 10:02p posts 1 users 1 website python.org

### 1 user in discussion

Content

People

Support

Translate

site design / logo © 2022 Grokbase