ECUACIÓN DEL CALOR

Trabajo realizado por estudiantes
Título	Ecuación del calor. Grupo GRwM
Asignatura	EDP
Curso	2023-24
Autores	Arturo Barrena y Mario Ríos
Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura

1 . Introducción

La visualización de soluciones de la ecuación del calor en una dimensión constituye una herramienta poderosa para desentrañar los intricados patrones de propagación térmica en sistemas unidimensionales. En este estudio, nos sumergimos en la riqueza matemática que rodea la representación gráfica de soluciones, explorando la diversidad de comportamientos térmicos que emanan de este modelo fundamental. A través de la creación visual, buscamos capturar la esencia dinámica de la ecuación del calor y ofrecer una perspectiva única sobre cómo la temperatura evoluciona en función del tiempo y la posición.

La ecuación del calor en una dimensión sirve como un punto de partida teórico sólido para comprender la difusión térmica en estructuras lineales. Al trazar diferentes soluciones, nos sumergimos en el tejido mismo de la propagación térmica, revelando detalles intrincados que, de otra manera, podrían escapar a una descripción puramente analítica. Este enfoque gráfico no solo ilustra la dinámica temporal de la temperatura, sino que también destaca la influencia de condiciones iniciales y parámetros en la evolución térmica a lo largo del espacio unidimensional.

2 . Planteamiento del sistema

Se considera una varilla metálica que ocupa el intervalo [0, 1] y que se encuentra aislada por su superficie lateral, de manera que la conducción de calor sólo se produce en la dirección longitudinal. La temperatura inicial de la varilla es 0 grados. En el extremo izquierdo se consigue mantener la temperatura a 0 grados mientras que en el derecho la temperatura es siempre de 1 grado. De acuerdo con todas estas indicaciones, el sistema de EDP que modeliza el comportamiento de la temperatura es el siguiente:

[math]\left \{ \begin{array}{ll} \frac{\partial u}{\partial t}- \frac{\kappa}{c} \frac{\partial ^2 u}{\partial x^2}=0, & \quad 0 \lt x \lt 1, t \gt 0, \\ u(x,0)=0, & \quad 0 \lt x \lt 1, \\ u(0,t)=0, & \quad t \gt 0, \\ u(1,t)=1 & \quad t \gt 0, \end{array} \right. [/math]

donde [math]u(x,t)[/math] es la temperatura en la posición [math]x[/math] en el tiempo [math]t[/math], [math]\kappa\gt0[/math] es la constante de conductividad térmica y [math]c[/math] es el calor específico del material, el cual también se considera constante. Cabe destacar que por simplificar el problema se suelen tomar los valores [math]\kappa[/math] y [math]c[/math] como constantes tal y como se ha hecho en este caso, ya que en bastantes casos de interés es razonable ignorar su variación. Si simplificamos aún más, suponiendo que ambas constantes son igual a 1, se tiene el siguiente sistema:

[math]\left \{ \begin{array}{ll} \frac{\partial u}{\partial t}-\frac{\partial ^2 u}{\partial x^2}=0, & \quad 0 \lt x \lt 1, t \gt 0, \\ u(x,0)=0, & \quad 0 \lt x \lt 1, \\ u(0,t)=0, & \quad t \gt 0, \\ u(1,t)=1 & \quad t \gt 0. \end{array} \right. [/math]

3 . Deducción de la Solución estacionaria

Antes de tratar de resolver el sistema, intentemos conjeturar qué podría suceder. Dado que la temperatura en el lado derecho es mayor que la del lado izquierdo según las condiciones frontera, el calor comenzará a fluir desde el extremo más caliente, aumentando la temperatura dentro de la barra y causando una salida de calor hacia el lado frío. Por otro lado, el aumento interior de la temperatura hace que el flujo de entrada de calor disminuya con el tiempo (ya que el flujo de entrada es exactamente la derivada de la temperatura con respecto a x cambiada de signo), mientras que el flujo de salida aumenta. Esperamos que tarde o temprano los dos flujos se equilibren entre sí y que la temperatura eventualmente alcance una distribución de estado estacionario. Por tanto, es razonable suponer que si se hace el límite cuando el tiempo tiende a infinito de la función temperatura [math]u[/math], se obtendrá una función que tan sólo depende de [math]x[/math], ya que la temperatura debe haber llegado al estado estacionario mencionado en algún momento. Por tanto, bajo esta suposición, la siguiente función [math]v[/math] está bien definida (depende sólo de x):

[math] v(x)=\lim_{{t \to \infty}} u(x,t). [/math]

Además, cumple el sistema de EDO:

[math]\left \{ \begin{array}{ll} \frac{\partial ^2 v}{\partial x^2}=0 & \quad 0 \lt x \lt 1, \\ v(0)=0, \\ v(1)=1. \\ \end{array} \right. [/math]

Si resolvemos este problema de valor inicial obtendremos:

[math]v(x)=x. [/math]

Es decir, que una vez que la temperatura se estabilice y deje de variar con el tiempo, seguirá el comportamiento de esta función lineal [math]v[/math]. Dicho resultado será muy importante ya que nos permitirá homogeneizar el sistema planteado en la anterior sección.

[math][/math]

4 . Resolución de la parte no estacionaria

Para comenzar a resolver el sistema de EDP inicial, homogeneizamos las condiciones frontera utilizando la solución estacionaria que hemos obtenido en el anterior apartado, planteando un problema equivalente con condiciones de frontera homogéneas. Para ello, es conveniente introducir la función:

[math] w(x,t)=u(x,t)-v(x) [/math]

De esta forma, la solución se compondrá de dos partes, la parte estacionaria y la parte no estacionaria, es decir que la podemos escribir de la siguiente forma:

[math]u(x,t)= v(x) + w(x,t) [/math]

donde [math]w(x,t)[/math] corresponde a la parte no estacionaria de la solución que es la abordaremos en esta sección. Veamos como podemos obtener las nuevas condiciones frontera y condiciones iniciales:

[math] w(x,0)=u(x,0)-v(x)=-x \\ w(0,t)=u(0,t)-v(0)= 0 - 0 = 0\\ w(1,t)=u(1,t)-v(1)= 1 - 1 = 0. [/math]

La función [math]w[/math] introducida nos permite llegar al sistema homogéneo, que será de la forma:

[math]\left \{ \begin{array}{ll} \frac{\partial w}{\partial t}-\frac{\partial ^2 w}{\partial x^2}=0 & \quad 0 \lt x \lt 1, t \gt 0, \\ w(x,0)=-x, & \quad 0 \lt x \lt 1, \\ w(0,t)=0, & \quad t \gt 0, \\ w(1,t)=0 & \quad t \gt 0. \end{array} \right. [/math]

Ahora procedemos a la resolución del sistema por medio del método de separación de variables, es decir suponiendo que la solución es de la forma [math]w(x,t)=X(x)T(t)[/math]. De esta suposición, se obtienen dos problemas de valor inicial uno asociado a la parte espacial de la solución X(x) y el otro asociado a la parte temporal T(t). El primero de ellos será:

[math] \left \{ \begin{array}{ll} X'' + \lambda X=0\\ X(0)= 0, X(1)= 0 \\ \end{array} \right. [/math]

Para este problema tendremos como solución cualquiera de las siguientes auto-funciones [math]X_k=B\sin{(k \pi x)} [/math] para [math]k \in N [/math], con [math]B[/math] constante. Por otro lado, para el problema [math] T'+ \lambda T= 0 [/math] obtendremos como solución cualquiera de las siguientes auto-funciones [math] T_k= Ce^{-k^2\pi^2t} [/math] para [math]k \in N [/math], con [math]C[/math] constante. Sin embargo, el producto entre cualquiera de estas auto-funciones nos da una función solución de la ecuación diferencial y las condiciones frontera, que sin embargo, no cumple con la condición inicial. Para tratar de conseguir una solución que cumpla también la condición inicial proponemos la siguiente solución resultante de hacer la suma infinita del producto de las auto-funciones y hacer el desarrollo en serie de Fourier de la condición inicial del sistema:

[math] w(x,t)=\sum_{k=1}^{\infty} \frac{2(-1)^{k}}{k \pi}\sin{(k \pi x)} e^{-k^2\pi^2t} [/math].

De modo que la expresión para la solución final será:

[math] u(x,t)=x + \sum_{k=1}^{\infty} \frac{2(-1)^{k}}{k \pi}\sin{(k \pi x)} e^{-k^2\pi^2t} [/math].

5 . Dibujo para [math]t\in[0,1][/math]

En esta sección veremos visualmente cuál es la forma de la solución obtenida anteriormente. Para ello, aproximaremos la forma de la solución considerando los 10 primeros términos de la serie. Además, tan sólo la dibujaremos en el intervalo temporal [math][0,1][/math], ya que este breve intervalo de tiempo es suficientemente grande como para llegar a observar claramente cómo la solución se acerca a su estado estacionario. Representarla para tiempos mayores que este no tiene sentido, teniendo en cuenta que la solución seguirá manteniéndose estacionaria con el tiempo tal y cómo ya sabemos.

Una vez dicho esto, dibujamos la aproximación de la solución u de forma tridimensional:

Aproximación de la función u, dibujada para [math](x,t)\in[0,1]\times[0,1].[/math]

Se puede ver claramente cómo a medida que avanza el tiempo la solución va llegando a un estado estacionario. Esto se ve ya que la solución comienza a dejar de variar con el tiempo. Además, este estado estacionario viene claramente dado por la solución estacionaria obtenida anteriormente [math]v(x)=x[/math]. En la gráfica, se puede ver cómo a medida que pasa el tiempo la solución se va aproximando a dicha solución estacionaria la cual al estar representada tridimensionalmente tiene forma de plano. Además, podemos observar las oscilaciones cuando [math]t=0[/math] resultantes de hacer la serie de Fourier en la condición inicial, aproximando mediante senos.

9 . Cambio de condiciones de frontera

Ahora cambiamos las condiciones en la frontera, en lugar de suponer que la temperatura es 0 grados en el extremo derecho, suponemos que la barra está aislada térmicamente en dicho extremo. De este modo obtendremos el siguiente sistema

[math]\left \{ \begin{array}{ll} \frac{\partial u}{\partial t}- \frac{\partial ^2 u}{\partial x^2}=0 & \quad 0 \lt x \lt 1, 0 \lt t \lt T, \\ u(0,t)=0, & \quad 0 \lt t \lt T, \\ u_{x}(1,t)=0 & \quad 0 \lt t \lt T, \\ u(x,0)=max \{0,1-4|x-1/2|\}, & \quad 0 \lt x \lt 1. \end{array} \right. [/math]

Se puede observar que el sistema ya está homogeneizado ya que las condiciones de frontera son nulas, por tanto solo haría falta resolver el problema por el método de separación de variables:

[math] u(x,t)=\sum_{k=1}^{\infty} c_{k}\sin{((\frac{\pi}{2}+k\pi)x)} e^{-(\frac{\pi}{2}+k\pi)^2t} [/math].

Véase que no podemos trabajar con la base trigonométrica usual [math]\{ 1, sen(\pi kx), cos((\pi kx)\}_{k=1}^{\infty} [/math] ya que no nos permite describir adecuadamente la solución. En su lugar, vamos a trabajar con la base [math]\{ \frac{1}{2}, sen((\frac{\pi}{2}+k\pi)x), cos((\frac{\pi}{2}+k\pi)x) \}_{k=1}^{\infty} [/math] que ha sido modificada de tal forma que podamos representar adecuadamente la solución .Extendiendo de forma impar la función [math]max \{0,1-4|x-1/2|\}[/math], vamos a tener la base [math]\{sen((\frac{\pi}{2}+k\pi)x)\}_{k=1}^{\infty} [/math] pues forma un conjunto completo en [math][0,1][/math] y es ortogonal. Es lógico que la base este compuesta por senos ya que al extender de forma impar es la manera más coherente de representar dicha función. Verificamos que la base es ortonormal:

[math] \lt sen((\frac{\pi}{2}+n\pi)x), sen((\frac{\pi}{2}+m\pi)x)\gt_{L^2(0,1)} = \int_{0}^{1}sen((\frac{\pi}{2}+n\pi)x) \cdot sen((\frac{\pi}{2}+m\pi)x) dx = 0[/math]

Gracias a esta comprobación podemos ver que la base es ortogonal y completa por tanto podemos calcular los [math] c_{k} [/math] correspondientes a la solución que son los coeficientes de Fourier en esta nueva base. Planteamos la expresión de los [math] c_{k} [/math]:

[math] c_{k} = \int_{0}^{1}sen((\frac{\pi}{2}+n\pi)x) \cdot f(x) dx [/math]

Si realizamos la integral obtenemos:

[math] c_{k} = (8 \pi n cos(\pi n) - 8 sen(\pi n))/(\pi^2 n^2) + 2[/math]

También aportaremos una aproximación de estos coeficientes usando el método numérico del trapecio, para poder representar adecuadamente la solución.

10 . Principio del Máximo

El teorema del principio del máximo nos permite asegurar la unicidad del sistema. Dicho teorema nos dice lo siguiente:

Sea [math] u \in C^{2,1} (Q_T) \cap C(\overline{Q_T})[/math] que verifica [math]u_t - \Delta u\leq 0 [/math] en [math]Q_T[/math].Entonces u alcanza su máximo en la frontera parabólica:

[math]\max \limits_{(x,t) \in \overline{Q_T}} u = \max \limits_{(x,t) \in \partial _P Q_T} u [/math]

Nuestra solución correspondiente al apartado anterior cumplen todas las hipótesis del teorema del máximo, en primer lugar [math] u \in C^{2,1} (Q_T) \cap C(\overline{Q_T})[/math] además por el apartado anterior sabemos que la u que es la solución correspondiente al sistema cumple la expresión [math]u_t - \Delta u\leq 0 [/math]. Además la condición inicial converge uniformemente de manera que la condición que necesitamos para la continuidad se cumple a la perfección. (HAY QUE JUSTIFICAR POR QUE LA CONDICION INICIAL CONVERGE UNIFORMEMENTE).

De esta forma si consideramos dos soluciones u y v de tal forma que definimos w=u-v que satisface:

[math]\left \{ \begin{array}{ll} w_t - \Delta w =0 & \quad x \in (0,1), T \gt 0 \\ w(0,t)=0, & \quad x \in (0,1) , \\ w(x,T)=0, & \quad x \in \{0,1\}, T \gt 0. \end{array} \right. [/math]

Aplicando el principio del máximo obtenemos: [math] \min \limits_{(x,t) \in \partial _P Q_T}w \leq w \leq \max \limits_{(x,t) \in \partial _P Q_T}w [/math] de aqui deducimos que [math]w \equiv 0[/math] por tanto podemos concluir que u=v, quedando mostrado que la solución del sistema es única.