[SciPy-user] how to resolve a cubic equation in scipy?

Johann Cohen-Tanugi cohen@slac.stanford....
Wed May 21 18:02:13 CDT 2008


hi,
just for fun I compared the 3 methods with timeit I got :
[cohen@jarrett ~]$ ipython testSolveCubic.py
using fsolve:
2.44568693982 1.84741111298e-12
time=  0.000438928604126

using roots:
[-1.22284347+2.30637633j -1.22284347-2.30637633j  2.44568694+0.j        ]
time=  0.002277135849

using brentq
2.44568693982
time=  3.31401824951e-05

So brentq is indeed by far the fastest. I am slightly surprised that 
fsolve outperforms roots, but I guess it could be due to the fact that 
the problem is too simple...
Anyway, for completeness here is the script :
from scipy import optimize,roots
import timeit

def myCubicEq(r):
    return 1.2*r**3  + r - 20

print 'using fsolve:'
start=timeit.time.time()
results=optimize.fsolve(myCubicEq, 5) 
stop=timeit.time.time()
print results,myCubicEq(results)
print 'time= ',stop-start

print '\n'
print 'using roots:'
start=timeit.time.time()
results2=roots([1.2, 0, 1, -20])
stop=timeit.time.time()
print results2
print 'time= ',stop-start

print '\n'
print 'using brentq'
start=timeit.time.time()
results3=optimize.brentq(myCubicEq,0,10)
stop=timeit.time.time()
print results3
print 'time= ',stop-start

cheers,
Johann


Anne Archibald wrote:
> 2008/5/16 Pauli Virtanen <pav@iki.fi>:
>   
>> Fri, 16 May 2008 21:03:22 +0800, zhang chi wrote:
>>     
>>>    I want to resolve a cubic equation 1.2r^3  + r - 20 = 0.
>>>       
>>>>> scipy.roots([1.2, 0, 1, -20])
>>>>>           
>> array([-1.22284347+2.30637633j, -1.22284347-2.30637633j,
>>       2.44568694+0.j        ])
>>     
>
> Given that the equation has just one real root, scipy's root-finders
> (e.g. brentq) should be reliable and fast.
>
> Of course, if the OP really only has the one equation to solve, it's
> done. But presumably they're working with a family of related
> cubics...
>
> Anne
> _______________________________________________
> SciPy-user mailing list
> SciPy-user@scipy.org
> http://projects.scipy.org/mailman/listinfo/scipy-user
>   


More information about the SciPy-user mailing list