in scipy.optimize.fsolve mit Gleichungen Verbesserung der Genauigkeit Integration Einbeziehung

Question

Ich versuche, eine Integralgleichung mit dem folgenden Code (irrelevante Teile entfernt) zu lösen:

def _pdf(self, a, b, c, t): 
    pdf = some_pdf(a,b,c,t) 
    return pdf 

def _result(self, a, b, c, flag): 
    return fsolve(lambda t: flag - 1 + quad(lambda tau: self._pdf(a, b, c, tau), 0, t)[0], x0)[0] 

, die eine Wahrscheinlichkeitsdichtefunktion nimmt und findet Ergebnis tau so dass das Integral von PDF von Tau nach unendlich ist gleich flag. Beachten Sie, dass x0 eine (Gleitkomma-) Schätzung des Stamms ist, der an anderer Stelle im Skript definiert ist. Beachte auch, dass flag eine extrem kleine Zahl ist, in der Größenordnung von 1e-9.

In meiner Anwendung findet fsolve nur etwa 50% der Zeit erfolgreich eine Wurzel. Es gibt oft nur x0 zurück, was meine Ergebnisse signifikant beeinflusst. Es gibt keine geschlossene Form für das Integral von pdf, also bin ich gezwungen, mich numerisch zu integrieren und fühle, dass dies eine gewisse Ungenauigkeit mit sich bringt?

EDIT:

Dies wurde seit löst ein anderes Verfahren als das unten beschrieben ist, aber ich mag quadpy Arbeit bekommen und sehen, ob die Ergebnisse überhaupt verbessern. Der spezifische Code, den ich versuche zu bekommen zu arbeiten wie folgt:

import quadpy 
import numpy as np 
from scipy.optimize import * 
from scipy.special import gammaln, kv, gammaincinv, gamma 
from scipy.integrate import quad, simps 

l = 226.02453163 
mu = 0.00212571582056 
nu = 4.86569872444 
flag = 2.5e-09 
estimate = 3 * mu 

def pdf(l, mu, nu, t): 
    return np.exp(np.log(2) + (l + nu - 1 + 1)/2 * np.log(l * nu/mu) + (l + nu - 1 - 1)/2 * np.log(t) + np.log(
     kv(nu - l, 2 * np.sqrt(l * nu/mu * t))) - gammaln(l) - gammaln(nu)) 


def tail_cdf(l, mu, nu, tau): 
    i, error = quadpy.line_segment.adaptive_integrate(
     lambda t: pdf(l, mu, nu, t), [tau, 10000], 1.0e-10 
    ) 
    return i 

result = fsolve(lambda tau: flag - tail_cdf(l, mu, nu, tau[0]), estimate) 

Als ich das laufen erhalte ich eine Behauptung Fehler von assert all(lengths > minimum_interval_length). Ich bin mir nicht sicher, wie ich das beheben soll. jede Hilfe würde sehr geschätzt werden!

Answer 1

Als Beispiel versuchte ich 1/x für die Integration zwischen 1 und alpha, um das Zielintegral 2.0 abzurufen. Diese

import quadpy 
from scipy.optimize import fsolve 

def f(alpha): 
    beta, _ = quadpy.line_segment.adaptive_integrate(
     lambda x: 1.0/x, [1, alpha], 1.0e-10 
     ) 
    return beta 

target = 2.0 
res = fsolve(lambda alpha: target - f(alpha), x0=2.0) 
print(res) 

kehrt richtig 7.3890561.

Das Versagen quadpy Behauptung

assert all(lengths > minimum_interval_length)

Sie Hilfe bekommen, dass die adaptive Integration seine Grenze erreicht: Entweder Ihre Toleranz ein wenig entspannen, oder verringern die minimum_interval_length (siehe here).