Utilize o Método das Linhas para resolver a equação de reação-difusão transiente com uma cinética de segunda ordem:
\frac{\partial C_{A}}{\partial t}=D \frac{\partial^{2} C_{A}}{\partial x^{2}}-k C_{A}^{2}
Sujeita à condição de contorno em x=L :
\left.\frac{\partial C_{A}}{\partial x}\right|_{x=L}=0
Enquanto em x=0 temos C_{A}(0, t)=1. A concentração inicial é C_{A}(x, 0)=0.
a) Se discretizarmos espaço com o método das diferenças finitas e usarmos Euler Explícito para a parte temporal, quais são as equações necessárias para calcular a concentrações C_{i}^{(k)} em cada nó i no passo de integração k ?
b) Se ao invés disso usarmos Euler Implícito, qual é o vetor residual e o Jacobiano do sistema não-linear?

Passo 1

O que o problema pede e dados fornecidos.
  • Equação PDE: \frac{\partial C_{A}}{\partial t}=D \frac{\partial^{2} C_{A}}{\partial x^{2}}-k C_{A}^{2}
  • Condições de contorno e inicial: \left.\frac{\partial C_{A}}{\partial x}\right|_{x=L}=0, \;C_{A}(0,t)=1, \;C_{A}(x,0)=0.
  • Vou discretizar o espaço em nós igualmente espaçados i=0,1,\dots,N com \Delta x = L/N. Assim x_i = i\Delta x. O nó i=0 corresponde a x=0 (Dirichlet) e i=N a x=L (Neumann).
  • Denoto por \;C_i^{(k)} a aproximação de C_A(x_i,t_k), com passo temporal \Delta t e t_k = k\Delta t.
  • Decidimos que o vetor de incógnitas por passo temporal será \;[C_1^{(k)},C_2^{(k)},\dots,C_N^{(k)}]^T. O valor \;C_0^{(k)} é conhecido (Dirichlet = 1) e será usado como dado.

Passo 2

Discretização espacial (diferenças finitas) — explicando a fórmula e o tratamento da Neumann.
  • Para o termo de segunda derivada usamos a diferença central clássica, que para um nó interior i dá (erro O(\Delta x^2)):
\left.\frac{\partial^{2}C}{\partial x^{2}}\right|_{x=x_i}\approx \frac{C_{i-1}-2C_i+C_{i+1}}{\Delta x^{2}}.
Explicação rápida: essa fórmula vem da soma das expansões de Taylor de C_{i\pm1} e cancelamento dos termos ímpares. É a escolha natural porque é simétrica e de segunda ordem em \Delta x.
  • Condição de Neumann \left.\frac{\partial C}{\partial x}\right|_{x=L}=0 no nó i=N: usamos um ponto fantasma C_{N+1} e a aproximação central da derivada:
\left.\frac{\partial C}{\partial x}\right|_{x=L}\approx\frac{C_{N+1}-C_{N-1}}{2\Delta x}=0\quad\Rightarrow\quad C_{N+1}=C_{N-1}.
Substituindo esse C_{N+1} na fórmula da segunda derivada em i=N obtemos (fórmula usando somente pontos reais):
\left.\frac{\partial^{2}C}{\partial x^{2}}\right|_{i=N}\approx \frac{C_{N-1}-2C_N+C_{N+1}}{\Delta x^{2}}=\frac{C_{N-1}-2C_N+C_{N-1}}{\Delta x^{2}}=\frac{2\,(C_{N-1}-C_N)}{\Delta x^{2}}.
  • Dirichlet em i=0: \;C_0^{(k)}=1 para todo k.

Passo 3

Euler Explícito (parte temporal) — derivação e fórmulas finais (resposta a a)).
  • Aproximação temporal explícita (Euler explícito):
\left.\frac{\partial C}{\partial t}\right|_{t=t_k}\approx\frac{C_i^{(k+1)}-C_i^{(k)}}{\Delta t}.
  • Substituindo a discretização espacial e a reação no esquema explícito temos, para nós interiores i=1,2,\dots,N-1:
\frac{C_i^{(k+1)}-C_i^{(k)}}{\Delta t} = D\frac{C_{i-1}^{(k)}-2C_i^{(k)}+C_{i+1}^{(k)}}{\Delta x^{2}}-k\big(C_i^{(k)}\big)^{2}.
Isolando C_i^{(k+1)} (forma prática de cálculo):
\boxed{\,C_i^{(k+1)} = C_i^{(k)} + \Delta t\left[ D\frac{C_{i-1}^{(k)}-2C_i^{(k)}+C_{i+1}^{(k)}}{\Delta x^{2}} - k\big(C_i^{(k)}\big)^{2}\right],\quad i=1,\dots,N-1.\,}
  • No nó i=N (Neumann) usamos a fórmula da segunda derivada com o ponto fantasma já eliminado:
\frac{C_N^{(k+1)}-C_N^{(k)}}{\Delta t} = D\frac{2\,(C_{N-1}^{(k)}-C_N^{(k)})}{\Delta x^{2}}-k\big(C_N^{(k)}\big)^{2}.
Isolando:
\boxed{\,C_N^{(k+1)} = C_N^{(k)} + \Delta t\left[ D\frac{2\,(C_{N-1}^{(k)}-C_N^{(k)})}{\Delta x^{2}} - k\big(C_N^{(k)}\big)^{2}\right].\,}
  • Lembrando que \;C_0^{(k+1)}=1 (condição de contorno fixa).
  • Observação prática (didática): esse esquema é simples de implementar nó-a-nó, mas é condicionalmente estável; para parte difusiva linear a condição CFL (aproximação) é \Delta t \lesssim \Delta x^{2}/(2D) — com reação a condição pode ficar ainda mais restritiva.

Passo 4

Euler Implícito — definição do residual para o sistema não-linear (início da resposta b)).
  • Euler implícito (backward Euler) aproxima \partial C/\partial t por diferença backward:
\frac{C_i^{(k)}-C_i^{(k-1)}}{\Delta t} \approx \left.\frac{\partial C}{\partial t}\right|_{t=t_k}.
Aqui as quantidades espaciais aparecem avaliadas no novo passo k, o que leva a um sistema não-linear por causa do termo -k(C_i^{(k)})^{2}.
  • Definimos o vetor incógnita no passo k como \mathbf{C}^{(k)}=[C_1^{(k)},\dots,C_N^{(k)}]^T. O vetor conhecido do passo anterior é \mathbf{C}^{(k-1)}.
  • Definimos o residual componente a componente R_i(\mathbf{C}^{(k)}) tal que a solução do sistema implícito satisfaz R_i(\mathbf{C}^{(k)})=0. Para i=1,\dots,N-1 (nós interiores):
\begin{aligned} R_i(\mathbf{C}^{(k)}) &= C_i^{(k)}-C_i^{(k-1)} -\Delta t\left[ D\frac{C_{i-1}^{(k)}-2C_i^{(k)}+C_{i+1}^{(k)}}{\Delta x^{2}} - k\big(C_i^{(k)}\big)^{2}\right]. \end{aligned}
  • Para i=N (Neumann), usando a fórmula \partial^2 C/\partial x^2\approx 2(C_{N-1}-C_N)/\Delta x^2, o residual é:
\begin{aligned} R_N(\mathbf{C}^{(k)}) &= C_N^{(k)}-C_N^{(k-1)} -\Delta t\left[ D\frac{2\,(C_{N-1}^{(k)}-C_N^{(k)})}{\Delta x^{2}} - k\big(C_N^{(k)}\big)^{2}\right]. \end{aligned}
  • Observação: C_0^{(k)}=1 é conhecida e, se aparecer em R_1, deve ser colocada no termo explícito (não é incógnita). Assim o sistema não-linear tem dimensão N (ou N incógnitas listadas).

Passo 5

Jacobiano do sistema não-linear (continuação da resposta b)) — cálculo dos elementos e estrutura.
  • O Jacobiano \mathbf{J} tem entradas J_{ij}=\partial R_i/\partial C_j avaliadas no chute atual (ou na solução, quando convergir). Vamos derivar as entradas para nós interiores i=1,\dots,N-1: usando a forma expandida do residual (para clareza):
\begin{aligned} R_i &= C_i^{(k)}-C_i^{(k-1)} -\Delta t\left(\frac{D}{\Delta x^{2}}C_{i-1}^{(k)} - \frac{2D}{\Delta x^{2}}C_i^{(k)} + \frac{D}{\Delta x^{2}}C_{i+1}^{(k)} - k (C_i^{(k)})^{2}\right) \\\ &= C_i^{(k)}-C_i^{(k-1)} -\Delta t\frac{D}{\Delta x^{2}}C_{i-1}^{(k)} + 2\Delta t\frac{D}{\Delta x^{2}}C_i^{(k)} -\Delta t\frac{D}{\Delta x^{2}}C_{i+1}^{(k)} + \Delta t\,k (C_i^{(k)})^{2}. \end{aligned}
  • Derivadas parciais (para i=1,\dots,N-1):
    • Para j=i-1: \frac{\partial R_i}{\partial C_{i-1}} = -\Delta t\,\frac{D}{\Delta x^{2}}.
    • Para j=i: diferencie termos que contêm C_i: \frac{\partial R_i}{\partial C_i} = 1 + 2\Delta t\frac{D}{\Delta x^{2}} + 2\Delta t\,k\,C_i^{(k)}.
    • Para j=i+1: \frac{\partial R_i}{\partial C_{i+1}} = -\Delta t\,\frac{D}{\Delta x^{2}}.
    • Todas as outras entradas na linha i são zero. Logo a estrutura é tri-diagonal nas linhas interiores.
  • Para a última linha i=N (Neumann) usamos a forma específica do residual; expandindo e diferenciando obtemos:
    • Para j=N-1: \frac{\partial R_N}{\partial C_{N-1}} = -\Delta t\,\frac{2D}{\Delta x^{2}}.
    • Para j=N: \frac{\partial R_N}{\partial C_{N}} = 1 + 2\Delta t\frac{D}{\Delta x^{2}} + 2\Delta t\,k\,C_N^{(k)}.
    • Não há coluna N+1 porque eliminamos o ponto fantasma.
  • Resumo da estrutura do Jacobiano (assumindo incógnitas C_1,\dots,C_N): é tri-diagonal com elementos (para i=1,\dots,N-1):
\boxed{\begin{aligned} J_{i,i-1} &= -\Delta t\,\frac{D}{\Delta x^{2}},\\[4pt] J_{i,i} &= 1 + 2\Delta t\,\frac{D}{\Delta x^{2}} + 2\Delta t\,k\,C_i^{(k)},\\[4pt] J_{i,i+1} &= -\Delta t\,\frac{D}{\Delta x^{2}}. \end{aligned}}
E na última linha i=N:
\boxed{\begin{aligned} J_{N,N-1} &= -\Delta t\,\frac{2D}{\Delta x^{2}},\\[4pt] J_{N,N} &= 1 + 2\Delta t\,\frac{D}{\Delta x^{2}} + 2\Delta t\,k\,C_N^{(k)}. \end{aligned}}
  • Observação prática (didática): ao montar \mathbf{R} e \mathbf{J} para uma iteração de Newton, avalie C_i^{(k)} no chute atual para obter os termos 2\Delta t\,k\,C_i^{(k)} no diagonal; resolva \mathbf{J}\,\delta = -\mathbf{R} para a correção \delta e atualize \mathbf{C}\leftarrow\mathbf{C}+\delta.]

Resposta

a) Para i=1,\dots,N-1:
C_i^{(k+1)} = C_i^{(k)} + \Delta t\left[ D\frac{C_{i-1}^{(k)}-2C_i^{(k)}+C_{i+1}^{(k)}}{\Delta x^{2}} - k\big(C_i^{(k)}\big)^{2}\right].
Para i=N:
C_N^{(k+1)} = C_N^{(k)} + \Delta t\left[ D\frac{2\,(C_{N-1}^{(k)}-C_N^{(k)})}{\Delta x^{2}} - k\big(C_N^{(k)}\big)^{2}\right].
b) Residual ( incógnitas \mathbf{C}^{(k)}=[C_1^{(k)},\dots,C_N^{(k)}]^T ): para i=1,\dots,N-1
R_i(\mathbf{C}^{(k)})=C_i^{(k)}-C_i^{(k-1)}-\Delta t\!\left[ D\frac{C_{i-1}^{(k)}-2C_i^{(k)}+C_{i+1}^{(k)}}{\Delta x^{2}} - k\big(C_i^{(k)}\big)^{2}\right],
e para i=N
R_N(\mathbf{C}^{(k)})=C_N^{(k)}-C_N^{(k-1)}-\Delta t\!\left[ D\frac{2\,(C_{N-1}^{(k)}-C_N^{(k)})}{\Delta x^{2}} - k\big(C_N^{(k)}\big)^{2}\right].
Jacobiano \mathbf{J} (tri-diagonal): para i=1,\dots,N-1
J_{i,i-1}=-\Delta t\,\frac{D}{\Delta x^{2}},\qquad J_{i,i}=1+2\Delta t\,\frac{D}{\Delta x^{2}}+2\Delta t\,k\,C_i^{(k)},\qquad J_{i,i+1}=-\Delta t\,\frac{D}{\Delta x^{2}},
e na linha i=N
J_{N,N-1}=-\Delta t\,\frac{2D}{\Delta x^{2}},\qquad J_{N,N}=1+2\Delta t\,\frac{D}{\Delta x^{2}}+2\Delta t\,k\,C_N^{(k)}.

Exercícios relacionados