Scipy's fsolve Alternativen?

Question

Ich versuche numerisch ein System von nicht-linearen Gleichungen zu lösen:

def func(p): 
    x, f = p 
    return (math.exp(-x/O)-f, 
      L - L*((1 - math.exp(-x/O))**W) - x*math.exp(-x/O)) 

und ich zur Zeit scipy.fsolve für sie in der folgenden Art und Weise bin mit:

x, f = fsolve(func, (10, 0.2)) 

Ich bin Sicher ist meine Art der Verwendung von fsolve korrekt: Es funktioniert perfekt für bestimmte Parameterbereiche. Allerdings ist es nicht ganz für einen anderen (zB O = 8, L = 1,67, B = 8) mit dem folgenden Fehler:

RuntimeWarning: The number of calls to function has reached maxfev = 600. 

Ich bin sicher, dass es gut auflösbar ist - es gibt (mindestens) Matlab-Tools, die es tun. Gibt es irgendetwas, was ich falsch mache oder gibt es einen anderen Löser, den ich ausprobieren könnte?

Vielen Dank im Voraus für Hinweise!

Answer 1

Nichtlineare Optimierung und Root-Finding ist leider empfindlich für die Wahl des Startpunktes.

In [19]: def func(p): 
    x, f = p    
    return [np.exp(-x/O) -f, L - L*((1 - np.exp(-x/O))**W) - x*np.exp(-x/O)] 
    ....: 

In [20]: O, L, W = 8, 1.67, 8 

In [21]: res = fsolve(func, [1, 1.2]) 

In [22]: res 
Out[22]: array([ 2.19804447, 0.75975782]) 

In [23]: func(res) 
Out[23]: [-2.2204460492503131e-16, -4.4408920985006262e-15] 

Beachten Sie, dass Ihre Parametrisierung überflüssig scheint: Die Position der Wurzel hängt nicht von L und können Sie alles in Bezug auf exp(-x/O) parametrieren, die wahrscheinlich zu den Solver es viel einfacher machen wird.

BEARBEITEN: Definieren y = exp(-x/O). Dann sagt Ihnen Ihre erste Gleichung, dass die Wurzelfindung tatsächlich eindimensional ist (was bedeutet, dass Sie z. B. brentq verwenden können, was robuster ist). Aber Sie können auch fsolve verwenden:

In [43]: ry = fsolve(lambda y: 1 - (1-y)**W + (O/L)*y*np.log(y), 0.4) 

In [44]: ry 
Out[44]: array([ 0.75975782])