Alternativas à análise de estabilidade de von neumann para métodos de diferenças finitas

Estou trabalhando na resolução das equações de poroelasticidade unidimensional acopladas (modelo de biot), dadas como:

$$-(\lambda+ 2\mu) \frac{\partial^2 u}{\partial x^2} + \frac{\partial p}{\partial x} = 0$$ $$\frac{\partial}{\partial t} \left[ \gamma p + \frac{\partial u}{\partial x}\right] -\frac{\kappa}{\eta}\left[\frac{\partial^2 p}{\partial x^2}\right] =q(x,t)$$ no domínio e com as condições de contorno: $\Omega=(0,1)$

$p=0, (\lambda + 2\mu)\frac{\partial u}{\partial x}=-u_0$ em e u = 0, \ frac {\ parcial p} {\ parcial x} = 0 em x = 1 .u = 0 , ∂ $x=0$x=$u=0, \frac{\partial p}{\partial x} = 0$$x=1$

Eu discretizei essas equações usando um esquema de diferença finita centrada:

γp t + 1 i -p t

$$(\lambda + 2\mu)\frac{u^{t+1}_{i+1}-2u^{t+1}_i+u^{t+1}_{i-1}}{\Delta x^2} + \frac{p^{t+1}_{i+1}-p^{t+1}_{i-1}}{2\Delta x} = 0$$ $$\gamma \frac{p^{t+1}_i-p^t_i}{\Delta t} + \frac{u^{t+1}_{i+1}-u^{t+1}_{i-1}}{2\Delta x \Delta t} -\left[ \frac{u^t_{i+1}-u^t_{i-1}}{2\Delta x \Delta t}\right] - \frac{\kappa}{\eta}\left[ \frac{p^{t+1}_{i+1} -2p^{t+1}_i + p^{t+1}_{i-1}}{\Delta x^2}\right]=q^{t+1}_i$$

Atualmente, estou elaborando os detalhes da convergência do esquema analisando sua consistência e estabilidade. A parte da consistência parece bastante direta para mim, mas já estou prevendo algumas dificuldades com a análise de estabilidade. Primeiro de tudo, existem duas variáveis ​​e duas equações. Em segundo lugar, há também um termo derivado espaço-temporal misto na segunda equação. Eu estou familiarizado com a análise de estabilidade de von neumann e posso ver que será muito difícil estabelecer estabilidade com esse método. Existem alternativas à análise von neumann que eu possa usar?

Se você substituir, pelo menos por sua análise, por , poderá escrever seu sistema como onde todas as constantes são definidas como e onde o subscrito se refere à discretização do espaço, tanto das variáveis ​​quanto dos operadores diferenciais. Seu esquema é então obtido aproximando via Euler implícito.$\frac{\partial u}{\partial x}$$u_x$

$$\begin{bmatrix} 0 & 0 \\ I & I \end{bmatrix} \frac{d}{dt} \begin{bmatrix} p_h(t) \\ u_{x,h}(t) \end{bmatrix}+\begin{bmatrix} -\partial_h & \partial_h \\ -\Delta_h & 0 \end{bmatrix}\begin{bmatrix} p_h(t) \\ u_{x,h}(t) \end{bmatrix}= \begin{bmatrix} q_h(t) \\ 0 \end{bmatrix} \quad (*)$$ $1$${}_h$$\frac{d}{dt}$

Agora a estrutura diferencial-algébrica (DAE) é evidente. Para as variáveis ​​existem equações diferenciais (no tempo) e algébricas.

Se você puder mostrar que é invertível, cf. esta pré - impressão [pág. 3] e a edição abaixo, em que o DAE é do índice 1 ou sabe-se que Euler implícito e sem estranheza é convergente, consulte o Teorema 5.12 neste livro . (Isenção de responsabilidade: este livro não está disponível gratuitamente e foi escrito pelo meu supervisor de doutorado)$\begin{bmatrix} -\partial_h & \partial_h \\ I & I \end{bmatrix}$

Com essa abordagem, você pode contornar a análise de estabilidade.

Para uma prova direta da estabilidade de , tentaria usar a Equação para aplicar a análise de estabilidade de von Neumann usando as funções próprias de e investigando o efeito de nas funções próprias.$L^2$$(*)$$\Delta_h$$\partial_h$

No entanto , se a estabilidade não puder ser estabelecida para , isso não significa que seu esquema não é convergente - devido à substituição de . De um modo geral, pode-se esperar estabilidade para esquemas que se aproximam das variáveis ​​reais, e não para esquemas que se aproximam de suas derivadas.$(*)$$u\leftarrow u_x$

APÊNDICE: Diz-se que um DAE é o índice 1, se puder ser transformado em uma EDO sem diferenciar as equações.

Digamos, o DAE está no formato A invertibilidade de implica que há uma transformação variável que eventualmente troca as colunas dos coeficientes para que com invertível (propriedade de classificação completa de ) e invertível (o complemento de Schur).

$$\begin{bmatrix} E_1 \\ 0 \end{bmatrix}\dot y + \begin{bmatrix} A_1 \\ A_2 \end{bmatrix} y = f.$$ $\begin{bmatrix} E_1 \\ A_2 \end{bmatrix}$$\tilde y \to y$$\begin{bmatrix} E_1 \\ A_2 \end{bmatrix} \to \begin{bmatrix} \tilde E_{11} & \tilde E_{12} \\ \tilde A_{21} & \tilde A_{22} \end{bmatrix}$$\tilde A_{22}$$A_2$$\tilde A_{11} - \tilde E_{12} \tilde A_{22}^{-1} \tilde A_{21}$

Para o sistema isso significa que a parte algébrica definida com pode ser usada para resolver uma parte de . Então, pode-se eliminar da parte diferencial (a segunda linha de bloco em ), para obter um ODE para as demais variáveis.$(*)$$A_2:=[-\partial_h ~ \partial_h]$$\tilde y_2$$(p_h,u_{x,h})$$\frac{d}{dt}\tilde y_2$$(*)$

Não estou familiarizado com as equações apresentadas aqui, mas lembro-me de aprender outro método para verificar a estabilidade de um esquema numérico em meus cursos. É conhecido como análise de equações modificadas.

Aqui está uma boa referência para isso,

http://193.146.160.29/gtb/sod/usu/$UBUG/repositorio/10291890_Warming.pdf

Na referência acima, a conexão entre a teoria da estabilidade baseada na análise de equações modificadas e a análise de estabilidade de Von Neumann é estabelecida.

Depois de um pouco de pesquisa on-line, me deparei com as seguintes referências,

Este artigo discute a modelagem por diferenças finitas das equações poroelásticas de Biot em frequências sísmicas. Também possui uma seção sobre estabilidade do esquema numérico.

Este artigo apresenta uma estratégia de solução para dissociar o sistema acoplado e verificar a estabilidade do esquema numérico.