Optimale ODE-Methode für eine feste Anzahl von RHS-Bewertungen

14

In der Praxis ist die Laufzeit der numerischen Lösung eines IVP x ( t 0 ) = x 0 wird oft von der Auswertungsdauer der rechten Seite (RHS)dominiert. Nehmen wir daher an, dass alle anderen Operationen sofort ausgeführt werden (dh ohne Rechenaufwand). Wenn die Gesamtlaufzeit zum Lösen des IVP begrenzt ist, ist dies gleichbedeutend mit der Begrenzung der Anzahl der Auswertungen vonauf einige.

x˙(t)=f(t,x(t)) for t[t0,t1]
x(t0)=x0
ffNN

Uns interessiert nur der Endwert .x(t1)

Ich suche nach theoretischen und praktischen Ergebnissen, die mir bei der Auswahl der besten ODE-Methode in einem solchen Umfeld helfen.

Wenn zum Beispiel ist, können wir den IVP mit zwei expliziten Euler-Schritten der Breite ( t 1 - t 0 ) / 2 oder einem Schritt der Breite t 1 - t 0 mit der Mittelpunktmethode lösen . Mir ist nicht sofort klar, welches man bevorzugt. Für größere N kann man natürlich auch über Mehrschrittmethoden, iterierte Runge-Kutta-Schemata usw. nachdenken.N=2(t1t0)/2t1t0N

Was ich suche, sind Ergebnisse ähnlich denen, die zum Beispiel für Quadraturregeln existieren: Wir können Gewichte { w i } und zugehörige Punkte { x i } so auswählen, dass die Quadraturregel n i = 1 w i ist g ( x i ) ist für alle Polynome g genau, so dass d e g ( g ) 2 n - 1 ist .n{wi}{xi}i=1nwig(xi)gdeg(g)2n1

Aus diesem Grund suche ich nach einer Ober- oder Untergrenze für die globale Genauigkeit von ODE-Methoden bei einer begrenzten Anzahl zulässiger Bewertungen der RHS . Es ist in Ordnung, wenn die Grenzen nur für einige RHS-Klassen gelten oder zusätzliche Einschränkungen für die Lösung x darstellen (genau wie das Ergebnis für die Quadraturregel, die nur für Polynome bis zu einem bestimmten Grad gilt).fx

EDIT: Einige Hintergrundinformationen: Dies ist für harte Echtzeitanwendungen, dh das Ergebnis muss vor einem bekannten Termin verfügbar sein. Daher die Begrenzung der Anzahl der RHS-Bewertungen N als dominierender Kostenfaktor. Typischerweise sind unsere Probleme steif und vergleichsweise klein.x(t1)N

EDIT2: Leider habe ich nicht die genauen Timing-Anforderungen, aber es ist sicher anzunehmen, dass eher klein sein wird (definitiv <100, wahrscheinlich näher an 10). In Anbetracht der Echtzeitanforderungen müssen wir einen Kompromiss zwischen der Genauigkeit der Modelle (wobei bessere Modelle zu längeren Ausführungszeiten der RHS und damit zu einem niedrigeren N führen ) und der Genauigkeit der ODE-Methode (wobei bessere Methoden höhere Werte erfordern) finden Werte von N ).NNN

Florian Brucker
quelle
Die übliche Entsprechung von Runge-Kutta-Verfahren mit Newton-Cotes-Verfahren mit festem Schritt gilt für den Fall, dass das RK-Verfahren auf das IVP angewendet wird, ; Das Anwenden der klassischen Methode vierter Ordnung auf diesen IVP ist beispielsweise gleichbedeutend mit dem Durchführen von Simpsons Regel für f ( x ) . y=f(x)f(x)
JM
@JM: Mir ist das bewusst. Ich wollte die Quadraturregeln nur als Beispiel für die Charakterisierung der Genauigkeit einer numerischen Methode für einen bestimmten Satz von Eingaben verwenden, wenn die Anzahl der Funktionsbewertungen begrenzt ist. Ansonsten interessieren mich "echte" ODEs, also solche, die sich nicht auf die Standardintegration reduzieren.
Florian Brucker
1
Das wird immer interessanter. Nun hat die Zahl für sich genommen keine Bedeutung. Es könnte hilfreich sein, λ N / T zu kennen , wobei T die Länge des Integrationsintervalls und λ die Lipschitz-Konstante von f in Bezug auf x ist . Dies wird uns sagen, wie steif das Problem wirklich ist. Vorausgesetzt, es ist steif, ist ein wahrscheinlicher Kandidat die BDF-Methode 2. Ordnung. NλN/TTλfx
David Ketcheson
@DavidKetcheson: Mich interessiert eher der allgemeine Ansatz zur Auswahl einer geeigneten Methode für ein bestimmtes Problem als die optimale Methode für ein bestimmtes Problem. Wir haben eine größere Anzahl von Modellen, die sich in Bezug auf Steifigkeit und Timing stark unterscheiden.
Florian Brucker
Sie sagen, dass die Bewertung von sehr teuer ist. Können Sie überhaupt einen Jakobianer berechnen? Was ist mit einer Annäherung, die die Hauptsteifigkeit korrigieren kann? Sie sind nicht in guter Verfassung, wenn Ihr Problem sehr schwergängig ist und Sie keine Möglichkeit haben, es zu beheben. f
Jed Brown

Antworten:

7

Ich denke, ein wichtiger Hinweis zur Beantwortung Ihrer Frage ist dieses Papier von Hosea und Shampine . Jetzt werde ich einige Hintergrundinformationen geben.

Im Allgemeinen kann die Schrittgröße, die Sie bei der numerischen Integration eines IVP verwenden können, durch Stabilität oder Genauigkeit eingeschränkt sein. Wenn Sie den hinsichtlich der Stabilität besten Löser auswählen möchten, müssen Sie den Bereich der absoluten Stabilität berücksichtigen . Bei einer einstufigen Methode ist dies

S={zC:|P(z)|1}.

Hier ist die Stabilitätsfunktion des Verfahrens; siehe zB den Text von Hairer et. al. Eine notwendige Bedingung für die Stabilität ist, dass λ h S ist, wobei λ über den Eigenwerten des Jacobians von f liegt und h die Schrittgröße ist. Dies ist nicht immer eine ausreichende Bedingung für nichtlineare Probleme, ist jedoch normalerweise eine gute Faustregel und wird in der Praxis angewendet.P(z)λhSλfh

Eine ausführliche Beschreibung des Problems, (explizite) Methoden zu finden, die große stabile Schrittgrößen ermöglichen, finden Sie in meinem Artikel über Stabilitätspolynome und in diesem Artikel über die Optimierung von Runge-Kutta-Methoden für Simulationen komprimierbarer Flüssigkeiten .

Stabilität ist das relevante Problem, wenn Sie feststellen, dass die größte stabile Schrittgröße bereits eine ausreichende Genauigkeit bietet. Andererseits kann die Schrittweite durch Ihre Genauigkeitsanforderungen eingeschränkt sein. In der Regel wird eine lokale Fehlerkontrolle durchgeführt. Die Lösung wird mit zwei Methoden berechnet, und ihre Differenz wird als Schätzung des Fehlers in der weniger genauen verwendet. Die Schrittweite wird adaptiv gewählt, um die vorgeschriebene Toleranz möglichst genau zu erreichen.

Zwei theoretische Maßnahmen sind der Schlüssel zur Vorhersage der Genauigkeitseffizienz. Die erste ist die Genauigkeitsreihenfolge des Verfahrens, die die Rate beschreibt, mit der sich der Fehler Null nähert, wenn die Schrittgröße verringert wird. Der zweite ist der Genauigkeits-Effizienz-Index (siehe das oben im ersten Satz aufgeführte Paper von Hosea und Shampine), der die in den Fehlerbegriffen auftretenden Konstanten berücksichtigt und Vergleiche zwischen Methoden derselben Reihenfolge ermöglicht.

Genauigkeit und Stabilitätseffizienz einer Vielzahl von Methoden können mit NodePy auf einfache und automatisierte Weise berechnet werden (Haftungsausschluss: NodePy wurde von mir entwickelt).

David Ketcheson
quelle
Vielen Dank. Das Papier von Hosea und Shampine ist in der Tat sehr interessant. Kennen Sie ähnliche Ergebnisse für steife Probleme? Mir ist bewusst, dass man normalerweise implizite Methoden für diese verwendet, aber diese sind nicht von vornherein an die Anzahl der RHS-Bewertungen gebunden, so dass sie in meinem Fall von geringem Nutzen sind.
Florian Brucker
Ich kenne so etwas nicht für steife Probleme, aber ich vermute, dass etwas existiert. Wie Sie sagen, ist die Frage bei der Verwendung impliziter Methoden subtiler. Ein Ansatz könnte darin bestehen, Rosenbrock-Methoden zu verwenden, die steife Probleme gut bewältigen, aber eine feste Anzahl von RHS-Bewertungen aufweisen.
David Ketcheson
6

In dieser Richtung gibt es nicht viele Ergebnisse, da es schwieriger ist, als nur die Genauigkeit zu korrigieren, da Sie aus Stabilitätsgründen häufig Zeitschritte auswählen müssen, die kleiner sind, als Sie für die gewünschte Genauigkeit benötigen würden. Die Ergebnisse werden also zwischen steifen und nicht steifen Fällen aufgeteilt. Im ersteren Fall werden die Anforderungen an Zeitschritte und RHS-Bewertungen im Allgemeinen nicht von der Genauigkeit bestimmt, und im letzteren Fall sind sie es.

Ich werde mich auf explizite Methoden konzentrieren, da der implizite Fall weitaus weniger offensichtlich ist, wie viele RHS-Bewertungen Sie verwenden müssen. Dies hängt ganz davon ab, wie Sie sich für die Lösung des resultierenden Systems entscheiden.

Für nicht steife Systeme:

Es gibt Stufenbeschränkungen für explizite Runge-Kutta-Methoden, für die angegeben ist, wie viele Stufen (RHS-Bewertungen) erforderlich sind, um eine bestimmte Genauigkeitsreihenfolge zu erreichen. Nach der vierten Ordnung übersteigt die Anzahl der Stufen die Genauigkeitsordnung und die Ungleichheit wächst weiter. Das große ODE-Buch von Butcher: http://books.google.com/books/about/Numerical_Methods_for_Ordinary_Different.html?id=opd2NkBmMxsC

leistet gute Arbeit bei der Erklärung einiger dieser 'Nichtexistenz'-Beweise.

Ihr Beispiel für eine Quadraturregel führt entweder zu einer mehrstufigen Methode wie Adams-bashforth oder zu sogenannten spektralverzögerten Korrekturmethoden. Für adams-bashforth benötigen Sie nur eine RHS-Bewertung pro Schritt, aber da die Stabilitätsbereiche für diese Methoden im Allgemeinen so klein sind, erledigen Sie im Allgemeinen den gleichen Arbeitsaufwand hinsichtlich der RHS-Bewertungen wie für eine Runge-Kutta-Methode mit denselben bestellen.

Hier ist ein Artikel zur spektralen verzögerten Korrektur:

https://www.google.com/search?q=spectral+deferred+correction&aq=f&oq=spectral+deferred+correction&aqs=chrome.0.57j0l2j62.3336j0&sourceid=chrome&ie=UTF-8

Ich bin mir nicht sicher, wie sich diese Integrationsmethoden gegenüber expliziten Standardmethoden verhalten. Sie erfordern häufig viel mehr Speicher, um Lösungszustände an Quadraturknoten zu speichern, und deshalb habe ich sie selbst nie verwendet.

Für steife Systeme:

S2S2S1S

Reid. Atcheson
quelle
2
Es kann vorkommen, dass die Verwendung einer Methode mit variablen Schritten (oder sogar einer Methode mit variabler Reihenfolge) effizienter ist, als hartnäckig an einer Methode mit festen Schritten festzuhalten. Man könnte zum Beispiel in Betracht ziehen, eine extrapolative Methode wie Bulirsch-Stoer zu verwenden: Führen Sie in einigen Schritten einige Auswertungen durch und bauen Sie dann (angeblich) genauere Schätzungen aus den Ergebnissen dieser Schritte auf.
JM
Wahr. Tatsächlich entsprechen viele der optimalen Methoden in gewissem Sinne einer variablen Schrittversion eines anderen Zeitschrittmachers. Runge-Kutta-Chebshev kann zum Beispiel als Vorwärts-Euler angesehen werden, wobei die variablen Zeitschritte Chebyshev-Punkte sind.
Reid.Atcheson
@JM: Genau. Aber gibt es eine Möglichkeit, die Genauigkeit dieser Ansätze in Bezug auf die Anzahl der RHS-Bewertungen zu beurteilen, abgesehen von numerischen Experimenten (die angesichts der hohen Anzahl möglicher Ansätze sehr aufwändig wären)?
Florian Brucker
@Florian, im Allgemeinen nicht. Sie haben von Lorenz 'Gleichungen gehört, nehme ich an?
JM
1
@JM: Ja :) Deshalb habe ich das Quadratur-Beispiel erwähnt, bei dem die Genauigkeit in einer Teilmenge (Polynomen) des ursprünglichen Problemraums gemessen wird. Ich würde mich über Ergebnisse freuen, die nur für eine bestimmte Teilmenge von Problemen funktionieren.
Florian Brucker
3

1014f(x)

Es gibt natürlich Ausnahmen (sehr große Systeme, sehr steife Systeme), aber ein allgemeines Gefühl in der Community ist, dass die Frage des Entwurfs von ODE-Solvern für "Standard" -Systeme gelöst ist. Aus diesem Grund denke ich, dass die von Ihnen gestellte Frage nicht sehr interessant ist - sie befasst sich mit einer Komponente des ODE-Solver-Designs, die nicht mehr von Bedeutung ist. Dies kann auch den Mangel an Literatur zu diesem Thema erklären.

Wolfgang Bangerth
quelle
+1. Wann immer jemand nach effizienten ODE-Lösern fragt, gehe ich davon aus, dass er an riesigen Systemen von ODEs interessiert ist, die aus PDE-Halbdiskretisierungen oder großen n-Körper-Problemen stammen.
David Ketcheson
Kannst du mir bitte erklären, wie das mit meiner Frage zusammenhängt? Ich sehe den Zusammenhang nicht, da ich mich für den Fall interessiere, dass die Bewertung f(x)nicht kostenlos, sondern so teuer ist, dass die Anzahl der Bewertungen begrenzt ist.
Florian Brucker
@ DavidKetcheson: Dies ist hier nicht der Fall. Es ist eher so, dass wir sehr strenge Timing-Anforderungen (harte Echtzeit) für schwache Hardware (eingebettete Geräte) haben. Die ODE-Systeme selbst sind vergleichsweise klein.
Florian Brucker
NNNN
NNN<1000
1

O(dim3)O(dim2)

Der erste Punkt ist also, sicherzustellen, dass Ihre RHS wirklich teurer ist als die zugrunde liegende lineare Algebra.

Der zweite Punkt: Aus der Literatur ist bekannt, dass Löser, die auf "teuren" Methoden (dh expliziten RK-Methoden) basieren, manchmal schneller sind als "billigere" Methoden (explizite Mehrschrittmethoden).

Zusammenfassend denke ich, dass Sie nicht nur die RHS-Bewertungen berücksichtigen sollten.

Faleichik
quelle
N