2014-03-30 11 views
21

unter Verwendung einer Gleichung I haben, wie folgt:Lösen einer Gleichung ein Python numerischen Solver in numpy

R - ((1.0 - np.exp(-tau))/(1.0 - np.exp(-a*tau))) = 0.

Ich möchte für tau in dieser Gleichung mit einem numerischen Solver lösen in numpy. Was ist der beste Weg, dies zu tun? Die Werte für R und a in dieser Gleichung variieren für verschiedene Implementierungen dieser Formel, sind jedoch auf bestimmte Werte festgelegt, wenn sie für Tau gelöst werden sollen.

+0

ich mir ziemlich sicher bin, dass diese Gleichung eine analytische Lösung in C – lejlot

+0

habe ich nicht so denken würde, wenn Sie es genau bestimmen können. Ich habe versucht, Newtons Methode mit Python zu verwenden, aber es scheint von der Art und Weise, wie die Gleichung geschrieben wird, abhängig zu sein. – stars83clouds

+0

probiere wolfram's alpha solver, wenn du ein bestimmtes R und ein a verwendest gibt es schöne, saubere Lösungen in C. – lejlot

Antwort

30

In herkömmlicher mathematischer Notation Ihre Gleichung

$$ R = \frac{1 - e^{-\tau}}{1 - e^{-a\cdot\tau}}$$

Die Durchsuchungen Funktion SciPy fsolve für einen Punkt, an dem ein gegebenen Ausdruck Null (eine „Null“ oder „Wurzel“ des Ausdrucks) entspricht. Sie müssen fsolve mit einer ersten Schätzung angeben, die Ihrer gewünschten Lösung "nahe" ist. Eine gute Möglichkeit, eine solche erste Schätzung zu finden, besteht darin, den Ausdruck einfach zu plotten und nach dem Nulldurchgang zu suchen.

#!/usr/bin/python 

import numpy as np 
import matplotlib.pyplot as plt 
from scipy.optimize import fsolve 

# Define the expression whose roots we want to find 

a = 0.5 
R = 1.6 

func = lambda tau : R - ((1.0 - np.exp(-tau))/(1.0 - np.exp(-a*tau))) 

# Plot it 

tau = np.linspace(-0.5, 1.5, 201) 

plt.plot(tau, func(tau)) 
plt.xlabel("tau") 
plt.ylabel("expression value") 
plt.grid() 
plt.show() 

# Use the numerical solver to find the roots 

tau_initial_guess = 0.5 
tau_solution = fsolve(func, tau_initial_guess) 

print "The solution is tau = %f" % tau_solution 
print "at which the value of the expression is %f" % func(tau_solution) 
7

Sie können die Gleichung umschreiben als

eq

  • Für integer a und Nicht-Null R Sie a Lösungen im komplexen Raum erhalten wird;
  • Es gibt analytische Lösungen für a=0,1,...4 (siehe here);

Also im Allgemeinen können Sie eine, mehrere oder keine Lösung haben und einige oder alle von ihnen können komplexe Werte sein. Sie können leicht scipy.root bei dieser Gleichung werfen, aber keine numerische Methode wird garantieren, alle Lösungen zu finden.

im komplexen Raum lösen:

import numpy as np 
from scipy.optimize import root 

def poly(xs, R, a): 
    x = complex(*xs) 
    err = R * x - x + 1 - R 
    return [err.real, err.imag] 

root(poly, x0=[0, 0], args=(1.2, 6)) 
Verwandte Themen