Hilft die Transformation von bei der numerischen Integration?

Ich habe anekdotisch gehört, dass, wenn man versucht, ein Integral der Form numerisch zu machen

\int_{0}^{\infty} f (x) J_{0} (x) d x

$\int_0^\infty f(x) J_0(x)\,\mathrm{d}x$

mit glatt und gutmütig (z. B. nicht selbst stark oszillierend, nicht singulär usw.), dann hilft es Genauigkeit, es als umzuschreiben $f(x)$

\frac{1}{π} \int_{0}^{π} \int_{0}^{\infty} f (x) \cos (x \sin θ) d x d θ

$\frac{1}{\pi}\int_0^\pi \int_0^\infty f(x) \cos(x\sin\theta) \,\mathrm{d}x\,\mathrm{d}\theta$

und führe zuerst das innere Integral numerisch durch. Ich kann keinen Grund sehen, warum ich damit rechnen sollte, aber andererseits ist die Genauigkeit einer numerischen Methode selten offensichtlich.

Natürlich weiß ich, dass der beste Weg, dies tatsächlich zu tun, darin besteht, eine Methode zu verwenden, die für solche oszillatorischen Integrale optimiert ist. Um der Neugier willen beschränke ich mich jedoch auf die Verwendung einer Quadraturregel. Kann jemand bestätigen oder widerlegen, dass diese Transformation die Genauigkeit des Integrals verbessert? Und / oder mich auf eine Quelle verweisen, die das erklärt?

quadrature special-functions David Z
quelle

Integriert über

... Dies ist eine der integralen Definitionen der Bessel-Funktion.

0 \leq θ \leq π

$0\le\theta \le\pi$

David Z

So Ihre Frage: generic Gegeben

-Punkt Quadraturformeln

auf

und

auf

, ist

N

$N$

Q_{\infty}^{N} [\cdot]

$Q_\infty^N[\cdot]$

[0, \infty)

$[0,\infty)$

Q_{π}^{N} [\cdot]

$Q_\pi^N[\cdot]$

[0, π]

$[0,\pi]$

schlechter oder besser als

Q_{\infty}^{N \cdot M} [f J_{0}]

$Q^{N\cdot M }_\infty[f\,J_0]$

Q_{π}^{M} [Q_{\infty}^{N} [f (x) \cos (x \sin θ)]]

$Q^M_\pi[Q_\infty^N[f(x)\,\cos(x\sin\theta)]]$

Stefano M

@StefanoM ja, das stimmt.

David Z

FWIW, eine der effizientesten Methoden zur Bewertung der Bessel-Funktion nullter Ordnung, ist die Trapezregel, die bekanntermaßen sehr genaue Ergebnisse liefert, wenn periodische Integranden über eine Periode integriert werden (sogar besser als die übliche Standard-Gaußsche Quadratur). Also: es könnte helfen, es könnte nicht helfen.

Antworten:

Ich denke nicht, dass es einen Unterschied macht. Sie müssen eine ausreichend hohe Quadratur für das Integral über wählen , damit es der Bessel-Funktion . Ich habe im folgenden Beispiel die Reihenfolge 20 gewählt, aber Sie müssen immer Konvergenz in Bezug auf die exakte Funktion und das Intervall durchführen, über das Sie integrieren. Dann habe ich Konvergenz mit , der Ordnung der Gaußschen Quadratur des Integrals über . Ich habe und Domain . Sie können ändern $\theta$ $J_0$ $n$ $x$ $f(x) = e^{-x} x^2$ $[0, x_\text{max}]$ $x_\text{max}$ unten. Ich habe:

 n      direct         rewritten
 1  0.770878284949  0.770878284949
 2  0.304480978430  0.304480978430
 3  0.356922151260  0.356922151260
 4  0.362576361509  0.362576361509
 5  0.362316789057  0.362316789057
 6  0.362314010897  0.362314010897
 7  0.362314071949  0.362314071949
 8  0.362314072182  0.362314072182
 9  0.362314072179  0.362314072179
10  0.362314072179  0.362314072179

$n=9$

Hier ist der Code:

from scipy.integrate import fixed_quad
from scipy.special import jn
from numpy import exp, pi, sin, cos, array

def gauss(f, a, b, n):
    """Gauss quadrature"""
    return fixed_quad(f, a, b, n=n)[0]

def f(x):
    """Function f(x) to integrate"""
    return exp(-x) * x**2

xmax = 3.

print " n      direct         rewritten"
for n in range(1, 20):
    def inner(theta_array):
        return array([gauss(lambda x: f(x) * cos(x*sin(theta)), 0, xmax, n)
            for theta in theta_array])
    direct = gauss(lambda x: f(x) * jn(0, x), 0, xmax, n)
    rewritten = gauss(inner, 0, pi, 20) / pi
    print "%2d  %.12f  %.12f" % (n, direct, rewritten)

xmax $[0, \infty]$ f(x)rewritten = gauss(inner, 0, pi, 20) / pi

Ondřej Čertík
quelle

Ich vermute, Sie haben Recht, meine eigenen Tests haben ähnliche Ergebnisse gezeigt.

David Z