(Δx, Δt) = (0.125, 0.0078125)
Métodos Numéricos Avanzados. 2023-2024
\[ \newcommand{\diff}{\mathrm{d}} \newcommand{\bm}[1]{\boldsymbol{#1}} \]
La ecuación de calor corresponde a \(F = -\nabla \rho\)
\[ \frac{\partial \rho}{\partial t} = \Delta \rho \]
También se puede escribir como \(F = - \rho \nabla \ln \rho\), y esto es interesante en algunos contextos. No lo haremos así.
Pongamos como dato inicial
\[ \rho_0 = \min\Big(1,\max(2-4|x|,0)\Big) \]
Puesto que \(F = -\nabla \rho\), una elección natural es \[ F_{i+\frac 1 2} = -\frac{\rho_{i+\frac 1 2} - \rho_{i-\frac 1 2}}{\Delta x} \] Resulta el método \[ \frac{d \rho_i}{d t} = - \frac{2 \rho_i (t) - \rho_{i+1} (t) - \rho_{i-1} (t) }{(\Delta x)^2} \]
Tenemos tres opciones sencillas:
Periódicas, igual que antes \(\rho_0 = \rho_N\) y \(\rho_{N+1} = \rho_1\).
Neumann: fijamos el flujo igual a cero \(F_{\frac 1 2} = F_{N+\frac 1 2} = 0\). Es decir, formalmente \(\rho_0 = \rho_{N+1} = 0\) y trabajar con todas las incógnitas interiores.
Dirichlet: fijamos \(\rho(t,-1) = \rho(t,1) = 0\). Podemos elegir \(\rho_1 = \rho_{N+1} = 0\) y trabajar con un punto menos, o ajustar la malla.
\[ \rho_i^{n+1} = \rho_i^n - \frac{\Delta t}{(\Delta x)^2} \left( 2 \rho_i^n - {\rho_{i+1}^n - \rho_{i-1}^n } \right) \]
No hay una condición CFL estrictamente hablando.
Si queremos que sea estable en el sentido \(|\rho_i^{n+1}| \le \max_j |\rho_j^n|\) basta que
\[ 2 \frac{\Delta t}{(\Delta x)^2} \le 1 \]
Esto obliga a que \(\Delta t\) se muy pequeño.
T = 2
N = 16;
x = range(-1.0,1.0,N+1); x = vec(x[1:end-1])
Δx = x[2] - x[1]
C = 0.5;
Δt = C*(Δx)^2;
ρ = ρ0(x);
@gif for n=1:ceil(T/Δt)
plot(x, ρ, ylim=[0.0,1.5],
title=@sprintf("t = %0.3f", float((n-1)*Δt)),
thickness_scaling = 1.5,linewidth=3,label="")
ρn = copy(ρ)
ρ[1] = ρn[1] - C*( 2*ρn[1] - ρn[N] - ρn[2] )
ρ[2:end-1] = ρn[2:end-1] - C*( 2*ρn[2:end-1] - ρn[1:end-2] - ρn[3:end] )
ρ[N] = ρn[N] - C*( 2*ρn[N] - ρn[1] - ρn[N-1] )
end
\[ \rho_i^{n+1} + \frac{\Delta t}{(\Delta x)^2} \left( 2 \rho_i^{n+1} - {\rho_{i+1}^{n+1} - \rho_{i-1}^{n+1} } \right) = \rho_i^n \]
Es incondicionalmente estable. Funciona incluso con un \(\Delta t\) muy grande.
Construímos \(A\) tal que \(A \bm \rho^{n+1} = \rho^n\).