La définition du système ODE rigide

17

Considérons un IVP pour le système ODE , y ( x 0 ) = y 0 . Le plus souvent, ce problème est considéré comme rigide lorsque la matrice de Jacobi fy=F(X,y)y(X0)=y0a à lafois desvaleurs propres avec une très grande partie réelle négative et des valeurs propres avec une très petite partie réelle négative (je ne considère que le cas stable).Fy(X0,y0)

En revanche, dans le cas d'une seule équation, par exemple l'équation de Prothero-Robinson , elle est dite rigide lorsque λ - 1 .y=λy+g+λgλ-1

Il y a donc deux questions:

  1. Pourquoi de petites valeurs propres sont incluses dans la définition de la rigidité pour les systèmes ODE? Je crois que la présence uniquement de très grandes pièces réelles négatives est assez suffisante pour que le système soit rigide, car cela nous oblige à utiliser de petits pas de temps pour des méthodes explicites.

  2. Oui, je sais que les problèmes de raideur les plus courants (résultant par exemple des PDE paraboliques) ont des valeurs propres grandes et petites. Donc, la deuxième question: existe-t-il un bon exemple naturel de grand système rigide sans très petites valeurs propres (ou alternativement avec un rapport doux )?λmax/λmin


OK, modifions la question. Considérons deux systèmes ODE linéaires bidimensionnels: d'abord avec les valeurs propres {-1000000, -0.00000001} et ensuite avec {-1000000, -999999}. Quant à moi, les deux sont raides. Mais si l'on considère la définition du rapport de rigidité, le deuxième système ne l'est pas. La question principale: pourquoi le rapport de rigidité est-il envisagé?

Et la deuxième partie de la question est toujours importante, paraphrase: je cherche un grand système ODE "naturel" avec de grandes valeurs propres négatives et un rapport de rigidité doux (pas supérieur, disons, à 100).

faleichik
la source
2
Bienvenue sur scicomp.se. Vos questions trouvent une réponse complète sur wikipedia: en.m.wikipedia.org/wiki/Stiff_equation
David Ketcheson
Je pense qu'entre le commentaire de @DavidKetcheson et les différentes sources que j'ai citées, vous verrez que le rapport de rigidité n'est qu'une indication. Ce n'est pas parfait; c'est pourquoi ce n'est pas dans la définition. Il se trouve que c'est une caractéristique de nombreux systèmes rigides, mais pas de tous. Et en ce qui concerne la deuxième partie, je pense que vous aurez du mal à la trouver à moins qu'elle n'ait une structure spéciale ou ne se présente dans une application. Je vous ai donné un exemple d'une telle application où le rapport de rigidité n'est pas toujours élevé, et je vous encourage à consulter le livre de Hairer et Wanner.
Geoff Oxberry
1
@David: Je ne peux pas être d'accord avec toi. Prenons par exemple le problème unidimensionnel y '= - 50 (y-cos x). La "valeur propre" est -50. On ne peut pas résoudre ce problème avec Euler explicite avec des pas supérieurs à 2/50. Si nous remplaçons -50 par -50000, la restriction sur le pas de temps devient 2/50000. Quelles «unités» pouvons-nous choisir ici pour surmonter cette barrière?
faleichik
2
@faleichik La partie de votre exemple fixe l'échelle de temps du "collecteur lent" (qui est probablement l'échelle de temps qui vous intéresse, bien qu'il soit concevable que vous soyez intéressé par des échelles de temps beaucoup plus courtes). Je ne crois pas qu'il soit possible de définir la rigidité sans choisir une échelle de temps d'observation (peut-être implicitement en énonçant des propriétés que vous souhaitez conserver sur des périodes plus longues). Le rapport de rigidité ne quantifie que la séparation des échelles entre les échelles de temps les plus rapides et les plus lentes du système autonome . cosX
Jed Brown
1
Il y a une nouvelle et meilleure réponse à cette question dans cet article .
David Ketcheson

Réponses:

10

La rigidité implique une certaine séparation des écailles. En général, si vous êtes intéressé par la phase du mode le plus rapide du système, vous devez le résoudre et le système n'est pas rigide. Mais souvent, vous vous intéressez à la dynamique à long terme d'un "collecteur lent" plutôt qu'à la vitesse précise à laquelle une solution du collecteur lent s'approche d'elle.

Les réactions chimiques et les écoulements réactifs sont des exemples courants de systèmes rigides. L' oscillateur van der Pol est un problème de référence commun pour les intégrateurs ODE qui a un paramètre de rigidité accordable.

Un océan est un autre exemple qui peut être utile à visualiser. Les tsunamis (ondes de gravité de surface) se déplacent à la vitesse d'un avion et produisent une structure de vagues complexe, mais se dissipent sur de longues échelles de temps et sont pour la plupart sans conséquence sur la dynamique à long terme de l'océan. Les tourbillons, ou d'autre part, voyagent environ 100 fois plus lentement à des vitesses assez piétonnes, mais provoquent des températures de mélange et de transport, de la salinité et des traceurs biogéochimiques qui sont pertinents. Mais la même physique qui propage une onde de gravité de surface prend également en charge un tourbillon (une structure de quasi-équilibre), donc la vitesse des tourbillons, le chemin sous Coriolis et le taux de dissipation dépendent de la vitesse de l'onde de gravité. Cela représente une opportunité pour un schéma d'intégration temporelle conçu pour les systèmes rigides de franchir l'échelle de temps de l'onde de gravité et de résoudre uniquement les échelles de temps dynamiques pertinentes. VoirMousseau, Knoll et Reisner (2002) pour une discussion de ce problème avec une comparaison des schémas de fractionnement et d'intégration du temps totalement implicite.

Connexes: Quand faut-il utiliser des méthodes implicites dans l'intégration des EDP hyperboliques?

Notez que les processus diffusifs sont généralement considérés comme rigides car l'échelle de temps la plus rapide dans le système discret dépend du maillage, avec une échelle de , mais l'échelle de temps de la physique pertinente est indépendante du maillage. En fait, les échelles de temps les plus rapides pour un maillage donné représentent une relaxation spatialement locale à la variété plus lente sur laquelle évoluent des échelles spatiales plus longues, de sorte que les méthodes implicites peuvent être très précises même dans des normes fortes même si elles ne résolvent pas les échelles les plus rapides.(ΔX)2

Jed Brown
la source
10

Partie 1

Les petites valeurs propres ne sont pas incluses dans la définition de la rigidité pour les systèmes ODE (problème de valeur initiale). Il n'y a pas de définition satisfaisante de la rigidité à ma connaissance, mais les meilleures définitions que j'ai rencontrées sont:

Si une méthode numérique avec une région finie de stabilité absolue, appliquée à un système avec des conditions initiales, est forcée d'utiliser dans un certain intervalle d'intégration une longueur de pas qui est excessivement petite par rapport à la régularité de la solution exacte dans cet intervalle , alors le système est dit rigide dans cet intervalle. (Lambert, JD (1992), Méthodes numériques pour les systèmes différentiels ordinaires , New York: Wiley.)

Un IVP [problème de valeur initiale] est rigide dans un certain intervalle[0,b]

Les équations rigides sont des équations où certaines méthodes implicites, en particulier BDF, fonctionnent mieux, généralement énormément mieux, que les méthodes explicites. (CF Curtiss & JO Hirschfelder (1952): Intégration d'équations rigides. PNAS, vol. 38, p. 235-243)

L'article de Wikipedia sur les équations rigides continue d'attribuer les «déclarations» suivantes à Lambert:

  1. Un système à coefficient constant linéaire est rigide si toutes ses valeurs propres ont une partie réelle négative et que le rapport de rigidité est grand.

  2. La rigidité se produit lorsque les exigences de stabilité, plutôt que celles de précision, contraignent la longueur du pas. [Notez que cette «observation» est essentiellement la définition d'Ascher et Petzold.]

  3. La rigidité se produit lorsque certains composants de la solution se désintègrent beaucoup plus rapidement que d'autres.

Chacune de ces observations a des contre-exemples (bien qu'il soit vrai que je n'ai pas pu en produire une du haut de ma tête).

Partie 2

Le meilleur exemple que je pourrais trouver serait probablement d'intégrer tout type de grand système de réaction de combustion dans la cinétique chimique dans des conditions qui provoquent l'inflammation. Le système d'équations sera rigide jusqu'à l'allumage, puis il ne sera plus rigide car le système a passé un transitoire initial. Le rapport de la plus grande à la plus petite valeur propre ne doit pas être grand, sauf autour de l'événement d'allumage, bien que de tels systèmes tendent à confondre les intégrateurs rigides, sauf si vous définissez des tolérances d'intégration extrêmement strictes.

Le livre de Hairer et Wanner donne également plusieurs autres exemples dans sa première section (Partie IV, section 1) qui illustrent de nombreux autres exemples d'équations rigides. (Wanner, G., Hairer, E., Résolution d'équations différentielles ordinaires II: Problèmes rigides et différentiels-algébriques (2002), Springer.)

Enfin, il convient de souligner l'observation de CW Gear:

Bien qu'il soit courant de parler de « équations différentielles raides, » une équation en soi est pas raide, un problème de valeur initiale particulière pour cette équation peut être raide, dans certaines régions, mais la taille de ces régions dépendent des valeurs initiales et la tolérance aux erreurs. (CW Gear (1982): Détection et traitement automatiques des équations différentielles ordinaires oscillatoires et / ou rigides. Dans: Intégration numérique des équations différentielles, Notes de cours en mathématiques., Vol. 968, p. 190-206.)

Geoff Oxberry
la source
Cher Geoff, merci pour la tolérance :-) Je voulais garder ma question simple, mais j'ai fini par être considéré comme inexpérimenté. En fait, je connais toutes ces définitions, mais.
faleichik
1. Les petites valeurs propres agissent implicitement dans la définition du rapport de rigidité: il est grand lorsque le démonateur est petit. 2. Pour le cas linéaire unidimensionnel, le rapport de rigidité est toujours un, même pour les équations rigides. 3. Avez-vous une référence pour le problème de cinétique chimique que vous avez suggéré? Et 4. Je vais essayer de clarifier la question dans les commentaires.
faleichik
2
Vous pouvez trouver les mécanismes chimiques au format CHEMKIN ici . Les problèmes sont suffisamment importants pour que des fichiers d'entrée soient nécessaires et les équations sont configurées automatiquement à l'aide d'un package de chimie. Je suggère d'utiliser les fichiers d'entrée en conjonction avec le package de chimie Cantera et la suite de solveurs ODE / DAE SUNDIALS , qui sont tous deux open source. Vous pouvez ensuite résoudre ces problèmes en C ++ ou MATLAB.
Geoff Oxberry
Personnellement, je prends la phrase Curtiss-Hirschfelder comme ma définition de travail de la rigidité; si RK explicite ou Adams prend trop de temps pour résoudre votre problème, alors il est probablement rigide.
JM
2

En fait, Jed Brown a réglé la question pour moi. Ce que je fais maintenant, c'est simplement de mettre ses mots en contexte.

  1. Les deux systèmes ODE linéaires 2D d'en haut sont rigides (c'est-à-dire difficiles à résoudre avec des méthodes explicites) sur des intervalles de temps relativement grands (par exemple [0,1]).

  2. Les systèmes linéaires avec un rapport de rigidité élevé peuvent être considérés comme "plus rigides" car il est très probable que l' on doive les intégrer sur un grand intervalle de temps. Cela est dû à des composants lents correspondant aux plus petites valeurs propres: la solution tend lentement vers l'état d'équilibre, et cet état d'équilibre est généralement important à atteindre.

  3. En revanche, l'intégration de systèmes à faible rapport de rigidité sur de grands intervalles n'est pas intéressante: dans ce cas l'état stationnaire est atteint très rapidement et on peut juste l'extrapoler.

Merci à tous pour cette discussion!

faleichik
la source
1

L'amplitude absolue des valeurs propres (dans un problème linéaire et autonome) n'a à elle seule aucune signification; c'est un artefact des unités dans lesquelles vous choisissez d'exprimer le problème.

La chaîne de commentaires devient incontrôlable, donc j'en fais une réponse. Je ne vais pas répondre à la question complète; comme je l'ai dit, voir wikipedia ou les autres réponses ici. Je réponds juste au bout qui dit

Considérons deux systèmes ODE linéaires bidimensionnels: d'abord avec les valeurs propres {-1000000, -0.00000001} et ensuite avec {-1000000, -999999}. Quant à moi, les deux sont raides. Mais si l'on considère la définition du rapport de rigidité, le deuxième système ne l'est pas. La question principale: pourquoi le rapport de rigidité est-il envisagé?

Bon, considérons un exemple du deuxième cas:

y1(t)=-1000000y1(t)
y2(t)=-999999y2(t)

t=1000000t

y1(t)=-y1(t)
y2(t)=-0,999999y2(t)

Remarque 1: J'ai choisi un système diagonal pour le rendre totalement évident, mais si vous l'essayez avec un autre système avec ces valeurs propres, vous verrez le même effet, car multiplier une matrice par une constante multiplie ses valeurs propres par la même constante.

|λ|1

David Ketcheson
la source
David, vous n'avez pas considéré l'intervalle d'intégration. Soit [0,1] dans le premier cas. En supposant que les limites de stabilité d'Euler soient explicites, le pas maximum autorisé est de 2/1000000. Nous devons donc faire au moins 500 000 pas. Lorsque vous modifiez le temps, la taille de pas maximale augmente à 2, mais tout l'intervalle d'intégration devient 1 000 000 et nous atteignons à nouveau le minimum de 500 000 pas.
faleichik
@faleichik Oui, vous l'avez maintenant. La rigidité n'a pas à voir avec l'amplitude absolue des valeurs propres mais avec leur taille par rapport à votre échelle de temps d'intérêt, comme Jed l'a noté ci-dessus.
David Ketcheson