Całkowanie numeryczne całki wielowymiarowej ze znanymi granicami

Mam (2-wymiarową) niewłaściwą całkę

I = \int_{A} \frac{W (x, y)}{F (x, y)} d x d y

$I=\int_A \frac{W(x,y)}{F(x,y)}\,\mbox{d}x\mbox{d}y$

gdzie domena integracji jest mniejsza niż , , ale dodatkowo ograniczone przez . Ponieważ i są gładkie, a $A$ $x=[-1,1]$ $y=[-1,1]$ $F(x,y)>0$ $F$ $W$ $W \ne 0$ na granicach, późniejsza relacja oznacza, że całka może być pojedyncza na granicach. Całka jest jednak skończona. Do tej pory obliczam tę całkę za pomocą zagnieżdżonej integracji numerycznej. To się udaje, ale powoli. Poszukuję bardziej odpowiedniej (szybszej) metody rozwiązania całki, być może metody Monte-Carlo. Potrzebuję jednak takiego, który nie stawia punktów na granicy nie-sześciennej dziedziny A i prawidłowo przyjmuje granicę niewłaściwej całki. Czy transformacja integralna może pomóc w tym ogólnym wyrażeniu? Zauważ, że mogę rozwiązać dla jako funkcję a nawet obliczyć dla kilku specjalnych funkcji wagowych $F(x,y)$ $y$ $x$ $I$ . $W(x,y)$

monte-carlo quadrature

— wysoki
źródło

F (x, y) \neq 0

$F(x,y) \ne 0$

A

$A$

F (x, y)

$F(x,y)$

Algorytm GSL QAGS: gnu.org/software/gsl/manual/html_node/… . Dzięki za zmiany (nie widziałem opcji pisania równań)!

— highsciguy

Odpowiedzi:

Zastrzeżenie: Napisałem pracę doktorską na temat kwadratury adaptacyjnej, więc ta odpowiedź będzie poważnie stronnicza w stosunku do mojej pracy.

QAGS GSL to stary QUADPACK integrator i nie jest on całkowicie niezawodny, szczególnie w obecności osobliwości. Zwykle powoduje to, że użytkownicy żądają o wiele więcej cyfr dokładności, niż faktycznie potrzebują, przez co integracja jest dość droga.

Jeśli używasz GSL, możesz wypróbować mój własny kod, CQUAD , opisany w tym artykule . Jest przeznaczony do radzenia sobie z osobliwościami, zarówno na krawędziach przedziału, jak i w domenie. Zauważ, że oszacowanie błędu jest dość solidne, więc pytaj tylko o tyle cyfr, ile faktycznie potrzebujesz.

Jeśli chodzi o integrację z Monte-Carlo, zależy to od tego, jakiej dokładności szukasz. Nie jestem również pewien, jak dobrze będzie działać w pobliżu osobliwości.

— Pedro
źródło

Na pewno się temu przyjrzę, ponieważ najłatwiej będzie go wdrożyć. W rzeczywistości zauważyłem, że procedura QAGS nie była superstabilna dla tego problemu.

— highsciguy

Czy istnieje sposób, aby wpłynąć na występowanie „GSL_EDIVERGE”? Wydaje się, że pojawia się dla niektórych parametrów.

— highsciguy

@highsciguy: Algorytm zwraca GSL_EDIVERGE, gdy uważa, że całka nie jest skończona. Gdybyś mógł podać mi przykład, dla którego zawodzi, mógłbym przyjrzeć się temu bliżej.

— Pedro

Trochę trudno jest wyodrębnić prostą procedurę, ponieważ jest ona osadzona w ogólnym kodzie dla całek n-dimensynowych. Zobaczę ... Ale dla ustalonego y, 1 / sqrt (F (x, y)) powinien zachowywać się jak 1 / sqrt (x), ponieważ x zbliża się do zer F (x, y) od F (x, y) można następnie zapisać jako wielomian w x. Ale może być tak, że zachowanie 1 / sqrt (x) zaczyna się późno. Może być również tak, że liczbowa prekursja całki nie jest zbyt dobra.

— highsciguy

@highsciguy: Tak, to zły pomysł. Większość reguł kwadraturowych zakłada, że całka ma pewien stopień gładkości, a jeśli ustawisz ją na zero z dowolnego punktu, wprowadzisz nieciągłość. Osiągniesz znacznie lepsze wyniki, jeśli użyjesz rzeczywistego odstępu!

— Pedro

Metody Monte Carlo zasadniczo nie mogą konkurować z kwadraturą adaptacyjną, chyba że masz całkę o dużych wymiarach, w której nie możesz sobie pozwolić na kombinatoryczną eksplozję punktów kwadratury z wymiarem.

$\int_{[0,1]^n} f(x)\; d^nx$ $n$ $M$ $M^n$ $k$ $N=(kM)^n$ $k$ $(2k-1)$ $e = {\cal O}(h^5)={\cal O}(M^{-(2k-1)})$

e = O (N^{- (2 k - 1) / n}) .

$e = {\cal O}(N^{-(2k-1)/n}).$

e = O (N^{- 1 / 2})

$e = {\cal O}(N^{-1/2})$

k > n / 4 + 1 / 2

$k>n/4+1/2$

$k\le 8$ $n=30$ $M=1$ $N=8^{30}$ punkty integracji, o wiele więcej niż kiedykolwiek w ocenie życia. Innymi słowy, tak długo, jak można ocenić wystarczającą liczbę punktów integracji, kwadratura na poddziałach domeny integracji jest zawsze bardziej wydajnym podejściem. Są to przypadki, w których masz całkę wielowymiarową, dla której nie możesz już oszacować punktów integracji nawet na jednym podziale, że ludzie używają metod Monte Carlo, pomimo ich gorszego porządku zbieżności.

— Wolfgang Bangerth
źródło

Wypróbuj zagnieżdżoną kwadraturę podwójnie wykładniczą (zobacz implementacje Ooura ). Technika ta wykorzystuje zmienną transformację, która sprawia, że transformowany całka zachowuje się bardzo płynnie na granicach i jest bardzo skuteczny w obsłudze osobliwości na granicy. Na jego stronie internetowej znajduje się również bardzo dobra lista odniesień do kwadratury DE.

— GertVdE
źródło