A definição de sistema rígido de EDO

17

Considere um IVP para o sistema ODE , . Geralmente esse problema é considerado rígido quando a matriz de Jacobi \ frac {\ parcial f} {\ parcial y} (x_0, y_0) tem ambos os valores próprios com parte real negativa muito grande e valores próprios com parte real negativa muito pequena (considero apenas o estável caso).y ( x 0 ) = y 0 fy=f(x,y)y(x0 0)=y0 0fy(x0 0,y0 0)

Por outro lado, no caso de apenas uma equação, por exemplo, a equação de Prothero-Robinson y=λy+g+λg , é chamada de rígida quando λ-1 1 .

Portanto, existem duas perguntas:

  1. Por que pequenos autovalores estão incluídos na definição de rigidez para sistemas ODE? Acredito que a presença de apenas partes reais negativas muito grandes é suficiente para o sistema ficar rígido, porque isso nos faz usar pequenos intervalos de tempo para métodos explícitos.

  2. Sim, eu sei que os problemas rígidos mais comuns (por exemplo, decorrentes de EDPs parabólicas) têm valores próprios grandes e pequenos. Portanto, a segunda pergunta: existe um bom exemplo natural de sistema rígido grande sem valores próprios muito pequenos (ou alternativamente com razão moderada λmax/λmin )?


OK, vamos modificar a pergunta. Considere dois sistemas ODE lineares bidimensionais: primeiro com valores próprios {-1000000, -0,00000001} e segundo com {-1000000, -999999}. Quanto a mim, os dois são rígidos. Mas se considerarmos a definição da taxa de rigidez, o segundo sistema não é. A principal questão: por que razão da rigidez é considerada?

E a segunda parte da pergunta ainda é importante, vamos parafrasear: estou procurando um sistema ODE grande e "natural" com grandes autovalores negativos e uma taxa de rigidez leve (não superior a, digamos, 100).

faleichik
fonte
2
Bem-vindo ao scicomp.se. Suas perguntas são respondidas completamente na wikipedia: en.m.wikipedia.org/wiki/Stiff_equation
David Ketcheson
Penso que, entre o comentário de @DavidKetcheson e as várias fontes que citei, você verá que a taxa de rigidez é apenas uma diretriz. Não é perfeito; é por isso que não está na definição. Acontece ser uma característica de muitos, mas não todos, sistemas rígidos. E quanto à segunda parte, acho que você terá dificuldade em encontrá-la, a menos que ela tenha uma estrutura especial ou surja em um aplicativo. Dei a você um exemplo de aplicação em que a taxa de rigidez nem sempre é grande e encorajo você a olhar para o livro de Hairer e Wanner.
precisa saber é o seguinte
11
@ David: Eu não posso concordar com você. Tomemos, por exemplo, problema unidimensional y '= - 50 (y-cos x). O "valor próprio" é -50. Não se pode resolver esse problema com Euler explícito com tamanhos de etapas maiores que 2/50. Se substituirmos -50 por -50000, a restrição no timestep se tornará 2/50000. Que "unidades" podemos escolher aqui para superar essa barreira?
faleichik
2
@faleichik A parte do seu exemplo fixa a escala de tempo do "coletor lento" (que provavelmente é a escala de tempo em que você está interessado, embora seja possível que você esteja interessado em escalas de tempo muito mais curtas). Não acredito que seja possível definir rigidez sem escolher uma escala de tempo observacional (talvez implicitamente afirmando propriedades que você deseja conservar por períodos mais longos). A taxa de rigidez quantifica apenas a separação de escalas entre as escalas de tempo mais rápidas e mais lentas do sistema autônomo . porquex
Jed Brown
11
Há uma nova e melhor resposta a essa pergunta neste documento .
David Ketcheson

Respostas:

10

A rigidez envolve alguma separação de escalas. Em geral, se você estiver interessado na fase do modo mais rápido do sistema, precisará resolvê-lo e o sistema não será rígido. Mas, freqüentemente, você está interessado na dinâmica de longo prazo de um "coletor lento", e não na taxa precisa com a qual uma solução do coletor lento se aproxima dele.

Reações químicas e fluxos de reação são exemplos comuns de sistemas rígidos. O oscilador van der Pol é um problema de referência comum para integradores de ODE que possui um parâmetro de rigidez ajustável.

Um oceano é outro exemplo que talvez seja útil para visualizar. Os tsunamis (ondas de gravidade da superfície) viajam na velocidade de um avião e produzem estruturas de ondas complexas, mas se dissipam em escalas de tempo longas e são principalmente irrelevantes para a dinâmica de longo prazo do oceano. Os redemoinhos, por outro lado, viajam cerca de 100 vezes mais devagar a velocidades bastante pedestres, mas causam mistura e temperatura de transporte, salinidade e traçadores biogeoquímicos relevantes. Mas a mesma física que propaga uma onda de gravidade superficial também suporta um redemoinho (uma estrutura de quase-equilíbrio), de modo que a velocidade do redemoinho, o caminho sob Coriolis e a taxa de dissipação dependem da velocidade da onda de gravidade. Isso representa uma oportunidade para um esquema de integração de tempo projetado para sistemas rígidos percorrer a escala de tempo da onda de gravidade e resolver apenas as escalas de tempo dinâmicas relevantes. VejoMousseau, Knoll e Reisner (2002) para discussão desse problema com uma comparação de esquemas de integração de tempo totalmente divididos e totalmente implícitos.

Relacionado: Quando métodos implícitos devem ser usados ​​na integração de EDPs hiperbólicas?

Observe que os processos difusivos são geralmente considerados rígidos porque a escala de tempo mais rápida no sistema discreto é dependente de malha, escalando com , mas a escala de tempo da física relevante é independente de malha. De fato, as escalas de tempo mais rápidas para uma determinada malha representam relaxamento espacialmente local para o coletor mais lento, no qual evoluem escalas espaciais mais longas, de modo que os métodos implícitos podem ser muito precisos, mesmo em normas fortes, apesar de não resolver as escalas mais rápidas.(Δx)2

Jed Brown
fonte
10

Parte 1

Valores próprios pequenos não são incluídos na definição de rigidez para sistemas ODE (problema de valor inicial). Não sei qual é a definição de rigidez satisfatória, mas as melhores definições que encontrei são:

Se um método numérico com uma região finita de estabilidade absoluta, aplicado a um sistema com quaisquer condições iniciais, for forçado a usar em um determinado intervalo de integração um comprimento de passo que seja excessivamente pequeno em relação à suavidade da solução exata naquele intervalo , então o sistema é considerado rígido nesse intervalo. (Lambert, JD (1992), Numerical Methods for Ordinary Differential Systems , Nova York: Wiley.)

Um IVP [problema de valor inicial] fica rígido em algum intervalo se o tamanho da etapa necessário para manter a estabilidade do método Euler para a frente for muito menor que o tamanho da etapa necessária para representar a solução com precisão. (Ascher, UM e Petzold, LP (1998), Métodos de Computador para Equações Diferenciais Ordinárias e Equações Diferenciais-Algébricas , Filadélfia: SIAM.)[0 0,b]

Equações rígidas são equações em que certos métodos implícitos, em particular o BDF, têm melhor desempenho, geralmente tremendamente melhor do que os explícitos. (CF Curtiss e JO Hirschfelder (1952): Integração de equações rígidas. PNAS, vol. 38, p. 235-243)

O artigo da Wikipedia sobre equações rígidas continua atribuindo as seguintes "declarações" a Lambert:

  1. Um sistema de coeficiente linear constante é rígido se todos os seus autovalores tiverem parte real negativa e a taxa de rigidez for grande.

  2. A rigidez ocorre quando os requisitos de estabilidade, e não os de precisão, restringem o comprimento da etapa. [Observe que essa "observação" é essencialmente a definição de Ascher e Petzold.]

  3. A rigidez ocorre quando alguns componentes da solução se deterioram muito mais rapidamente que outros.

Cada uma dessas observações tem contraexemplos (embora, reconhecidamente, eu não tenha conseguido produzir um em cima da minha cabeça).

Parte 2

Provavelmente, o melhor exemplo que eu poderia ter seria integrar qualquer tipo de grande sistema de reação de combustão na cinética química em condições que resultem em ignição. O sistema de equações ficará rígido até a ignição e depois não será mais rígido porque o sistema passou por um transiente inicial. A proporção entre o maior e o menor valor próprio não deve ser grande, exceto em torno do evento de ignição, embora esses sistemas tendam a confundir integradores rígidos, a menos que você defina tolerâncias de integração extremamente estritas.

O livro de Hairer e Wanner também fornece vários outros exemplos em sua primeira seção (Parte IV, seção 1) que ilustram muitos outros exemplos de equações rígidas. (Wanner, G., Hairer, E., Resolvendo Equações Diferenciais Ordinárias II: Problemas Rígidos e Diferenciais-Algébricos (2002), Springer.)

Por fim, vale ressaltar a observação do CW Gear:

Embora seja comum falar sobre "equações diferenciais rígidas", uma equação por si só não é rígida, um problema de valor inicial específico para essa equação pode ser rígido em algumas regiões, mas o tamanho dessas regiões depende dos valores iniciais e do valor inicial. tolerância a erros. (CW Gear (1982): Detecção e tratamento automáticos de equações diferenciais ordinárias oscilatórias e / ou rígidas. Em: Integração numérica de equações diferenciais, Notas de aula em Math., Vol. 968, p. 190-206.)

Geoff Oxberry
fonte
Caro Geoff, obrigado pela tolerância :-) Eu queria manter minha pergunta simples, mas acabei sendo considerada inexperiente. Na verdade, eu sei todas essas definições, mas.
faleichik
1. Autovalores pequenos atuam implicitamente na definição da razão de rigidez: é grande quando o demoninador é pequeno. 2. No caso linear unidimensional, a taxa de rigidez é sempre uma, mesmo para equações rígidas. 3. Você tem alguma referência para o problema de cinética química que sugeriu? E 4. Vou tentar esclarecer a pergunta nos comentários.
faleichik
2
Você pode encontrar mecanismos químicos no formato CHEMKIN aqui . Os problemas são grandes o suficiente para que os arquivos de entrada sejam necessários e as equações são configuradas automaticamente usando um pacote de química. Sugiro usar os arquivos de entrada em conjunto com o pacote químico Cantera e o SUNDIALS , o pacote de soluções ODE / DAE , ambos de código aberto. Você pode resolver esses problemas em C ++ ou MATLAB.
precisa saber é o seguinte
Pessoalmente, tomo a sentença Curtiss-Hirschfelder como minha definição de rigidez; se RK ou Adams explícitos estiverem demorando muito para resolver seu problema, provavelmente será difícil.
JM
2

De fato, Jed Brown esclareceu a questão para mim. O que estou fazendo agora é apenas colocar as palavras dele no contexto.

  1. Os dois sistemas ODE lineares 2d acima são rígidos (isto é, difíceis de resolver com métodos explícitos) em intervalos de tempo relativamente altos (por exemplo, [0,1]).

  2. Os sistemas lineares com alta taxa de rigidez podem ser considerados "mais rígidos" porque provavelmente é necessário integrá-los em um grande intervalo de tempo. Isso ocorre devido aos componentes lentos que correspondem aos menores autovalores: a solução tende lentamente ao estado estacionário, e esse estado estacionário é geralmente importante de ser alcançado.

  3. Por outro lado, não é interessante a integração de sistemas com pequena taxa de rigidez em grandes intervalos: nesse caso, o estado estacionário é atingido muito rapidamente e podemos extrapolá-lo.

Obrigado a todos por esta discussão!

faleichik
fonte
1

A magnitude absoluta dos valores próprios (em um problema linear e autônomo) por si só não tem nenhum significado; é um artefato das unidades em que você escolhe expressar o problema.

A cadeia de comentários está ficando fora de controle, então estou fazendo disso uma resposta. Eu não vou responder a pergunta completa; como eu disse, veja a Wikipedia ou as outras respostas aqui. Só estou respondendo o que diz

Considere dois sistemas ODE lineares bidimensionais: primeiro com valores próprios {-1000000, -0,00000001} e segundo com {-1000000, -999999}. Quanto a mim, os dois são rígidos. Mas se considerarmos a definição da taxa de rigidez, o segundo sistema não é. A principal questão: por que razão da rigidez é considerada?

Ok, vamos considerar um exemplo do segundo caso:

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

t=1000000t

y1 1(t)=-y1 1(t)
y2(t)=-0.999999y2(t)

Nota 1: Escolhi um sistema diagonal para torná-lo totalmente óbvio, mas se você tentar com outro sistema com esses autovalores, verá o mesmo efeito, pois multiplicar uma matriz por uma constante multiplica seus autovalores pela mesma constante.

|λ|1 1

David Ketcheson
fonte
David, você não considerou o intervalo de integração. Seja [0,1] no primeiro caso. Supondo limitações explícitas de estabilidade do Euler, a etapa máxima permitida é 2/1000000. Então, precisamos fazer pelo menos 500 000 etapas. Quando você escala o tempo, o tamanho máximo da etapa aumenta para 2, mas todo o intervalo de integração se torna 1 000 000 e atingimos o mínimo de 500 000 etapas novamente.
faleichik
@faleichik Sim, agora você já conseguiu. A rigidez não tem a ver com a magnitude absoluta dos valores próprios, mas com seu tamanho em relação à sua escala de tempo de interesse, como Jed observou acima.
David Ketcheson