[SciPy-User] Simple root solving issues

[josef.pktd at gmail.com]

On Sun, Oct 24, 2010 at 5:38 PM, Gökhan Sever <gokhansever at gmail.com> wrote:
> On Sun, Oct 24, 2010 at 4:29 PM, Warren Weckesser
> <warren.weckesser at enthought.com> wrote:
>>
>> Your function is basically c1 - c2/log(gsd), so yeah, it will have a problem
>> if the initial guess is 1.0.  It works fine if the initial guess is, say,
>> 1.01.
>
> Yes, catch this better now. Interestingly it suffers when x0 < 1.00
> too. Probably log(gsd<1) explodes the function by making the right
> term positive and restricting the function to find a zero.
>
> --
> Gökhan
> _______________________________________________
> SciPy-User mailing list
> SciPy-User at scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-user
>

with a bit of random search
--------
#!/usr/bin/env python

import numpy as np
from scipy.optimize import leastsq, fsolve

def fitfunc(gsd):

    #return dH_dlogDP1 -
    (h1.sum()/((2*np.pi)**0.5*np.log(gsd)))*np.exp(-(np.log(Dp1)-np.log(Dp1))**2./(2.*np.log(gsd)**2.))
    return dH_dlogDP1 - (h1.sum()/((2*np.pi)**0.5*np.log(gsd)))

dH_dlogDP1 = np.array([ 869.11014589])
Dp1 = np.array([ 0.02994996])
h1 = np.array([ 1906.7283])

res = leastsq(fitfunc, x0=1.5)

np.random.seed(5)
lb = -20
ub = 20
for x0 in lb + (ub-lb)*np.random.random(50):
    res2 = fsolve(fitfunc, x0=x0)
    if np.abs(res2 - x0) > 1e-6:
        print x0, res2
----

prints

>>>
0.736719514918 [ 0.09143545]
3.19351240758 [ 2.39943611]
3.997167865 [ 2.39943611]
3.10651432519 [ 2.39943611]
0.618904476216 [ 0.15916137]
1.8582599189 [ 2.39943611]