Diferencia entre revisiones de «Ecuación del calor(Grupo Eau De Parfum(EDP)))»

Revisión del 20:14 6 mar 2024

1 Introducción

2 Conocimientos previos

Antes de comenzar será de gran ayuda introducir una serie de conceptos que utilizaremos a lo largo este trabajo:

• Ley del calor de Newton: describe la tasa de pérdida o ganancia de calor de un objeto en relación con el ambiente circundante.La ley se expresa comúnmente como:

[math] \frac{d T(t)}{d t} = - k (T(t) - T_{\mathrm{a}}) [/math],

siendo [math]T(t) [/math]la temperatura del objeto en un tiempo dado, [math]T_{\mathrm{a}}[/math] la temperatura ambiente y [math]k[/math] la constante de proporcionalidad.

• Principio de conservación de la energía calorífica: establece que la variación de energía calorífica sobre un cuerpo V se debe al balance entre el calor que entra y sale del cuerpo más una producción externa.

• Ley de Fourier:intaura que el calor fluye desde regiones de alta temperatura a regiones de baja temperatura, y la cantidad de flujo de calor depende de la diferencia de temperatura y las propiedades del material.La expresión matemática de la Ley de Fourier es:

[math] \mathbf{q} = - k {\nabla} T [/math],

donde [math] \mathbf{q}[/math] es el vector de flujo de calor por unidad de superficie y [math] {\nabla} T [/math] es el gradiente del campo de temperatura en el interior del material.

• Principio del máximo: Sea [math]u \in C^{2,1}(Q_T)\cap (\overline{Q_T})[/math] tal que [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]\partial _p Q_T[/math]:

[math]\max_{\overline{Q_T}} u = \max_{\partial _p Q_T} u[/math]

3 Primer problema

Consideramos una varilla metálica de longitud 1 m , que se encuentra aislada por su superficie lateral, de manera que la conducción de calor únicamente se produce en la dirección longitudinal. La temperatura inicial de la varilla es 0ºC. En el extremo izquierdo se consigue mantener la temperatura a 0 ºC mientras que en el derecho la temperatura es siempre de 1 ºC. Además tanto la conductividad térmica como el calor específico toman un valor constante de 1.Por lo que nuestro problema a resolver queda:

Este problema se resuelve mediante el método de variación de las constantes, pero un factor clave para su implementación es que las condiciones de contorno estén homogeneizadas, es decir, que sean igual a cero. Esto se conseguirá mediante un cambio de variable, en el que es imprescindible la llamada solución estacionaria. La solución estacionaria de nuestro problema es aquella que no varía con el tiempo [math] v=\lim_{t \to \infty} u(x,t)[/math]y por tanto, su derivada respecto al tiempo es nula. Nuestro sistema queda ahora:

Resolviendo obtenemos que [math] v(x)=x [/math] es nuestra solución estacionaria, cuya expresión gráfica es:

import matplotlib.pyplot as plt
import numpy as np
# Define el rango de x entre 0 y 1
x = np.linspace(0, 1, 100)  
# Define la función  estacionaria ü(x) = x
ü_x = x
# Crea la gráfica
plt.plot(x, v_x)
# Etiquetas de los ejes
plt.xlabel('x')
plt.ylabel('v(x)')
# Título de la gráfica
plt.title('Gráfica de la solución estacionaria v(x)=x')
# Muestra la gráfica
plt.show()

Una vez hallada la solución estacionaria, podemos homogeneizar nuestro problema mediante el cambio de variable [math]w(x,t)= u(x,t)-v(x)[/math].Por tanto, sustituyendo en el sistema inicial obtenemos:

Nuestro sistema es ya homogéneo,por lo que aplicamos el método de separación de variables (que consiste en expresar la solución como el producto de una función que depende únicamente de la variable espacial x y una función dependiente de la variable temporal) y obtenemos como solución homogénea:

Y deshaciendo el cambio, nuestra solución al problema original es por tanto:

Mostramos gráficamente la solución, observando que se cumplen las condiciones de frontera y además según avanza el tiempo, la solución deja de depender del tiempo y alcanza un estado estacionario, determinado justamente por la solución estacionaria obtenida anteriormente, [math] \ddot u(x)=x [/math]. Cabe destacar que en [math] t=0 [/math] obtenemos las oscilaciones correspondientes a la serie de Fourier de la condición inicial, que se aproxima mediante senos.

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# Parámetros
L = 1       # Longitud de la varilla
T = 1      # Tiempo final
D=1       # Coeficiente de conductividad térmica
N = 10      # Número de términos en la serie
x_points = 100  # Puntos en la dirección x
t_points = 100  # Puntos en la dirección t

# Espacios
x = np.linspace(0, L, x_points)
t = np.linspace(0, T, t_points)
X, T = np.meshgrid(x, t)

# Solución usando los primeros N términos de la serie
def u(x, t, N):
    sum = 0
    for k in range(1, N+1):
        lambda_k = k**2 * np.pi**2 
        #Calculasmos los coeficientes C_k
        C_k = 2 * (-1)**(k) / (k * np.pi)  
        sum += C_k * np.sin(k*np.pi * x) * np.exp(-lambda_k * t) 
    return x+sum

# Calcular la solución en la malla
Z = u(X, T, N)

# Graficar
fig = plt.figure(figsize=(12,8))
ax = fig.add_subplot(111, projection='3d')
surf=ax.plot_surface(X, T, Z, cmap='viridis')

ax.set_xlabel('Posición x')
ax.set_ylabel('Tiempo t')
ax.set_zlabel('Temperatura T(x,t)')
ax.set_title('Evolución de la Temperatura en la Varilla')

fig.colorbar(surf, shrink=0.5, aspect=5)
plt.show()

4 Cambio de condiciones

Supongamos ahora que la temperatura en el extremo derecho es [math]0 Cº[/math] y, además, inicialmente la temperatura de la barra viene dada por la función [math]u(x,0) = max\{0, 1−4|x−1/2|\}[/math]. Entonces, obtenemos este sistema:

Supongamos que, en tiempo infinito, \( t \rightarrow \infty \), se llega a una solución estacionaria y que, por tanto, ya no varía en tiempo. De modo, que tomando esta suposición llegamos al siguiente sistema con \( u(x,t) \xrightarrow{\text{t tiende a }\infty} v(x) \).

Cuya solución es la constante 0. De modo que el sistema homogeneizado, coincide con el original, pues \( w(x,t)=u(x,t)-v(x) =u(x,t)-0=u(x,t)\) Antes de exigir la condición inicial, la solución coincide con la del primer apartado, pues nos encontramos ante la misma ecuación diferencial, salvo por la función \( u(x,0)\), la cual no toma parte hasta el final de la resolución de esta suponiendo separación de variables. De modo que

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

Para hallar el valor de los coeficientes [math] C_k[/math], buscamos la serie de Fourier de la función [math] max\{0, 1−4|x−1/2|\}[/math] extendida de forma impar. Esto es posible, ya que

[math] \int^1_0 (max\{0, 1−4|x−1/2|\})^2 dx=0.1667 \lt \infty[/math]

El cálculo de los coeficientes se ha realizado numéricamente aproximando por la fórmula del trapecio. Los coeficientes se pueden consultar ejecutando el siguiente código en Python

import numpy as np 
def f(x):
    return np.maximum(0, 1 - 4 * np.abs(x - 1/2))

def fourier_coefficients(f, n, x):
    coefficients = np.zeros(n)
    for k in range(1,n+1): #empiezo en k=1 o en k=0??
        coefficients[k-1] = 2 * np.trapz(f(x) * np.sin(k * np.pi * x), x)
    return coefficients

fourier_coefficients(f, 10, np.linspace(0, 1, 1000))

De modo que ya se ha calculado la solución de este nuevo problema. Visualmente responde a la siguiente gráfica, resultado de ejecutar este código:

import sympy as sp
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# Parámetros
L = 1       # Longitud de la varilla
T = 1       # Tiempo final
D = 1       # Coeficiente de conductividad térmica
N = 10      # Número de términos en la serie
# Espacios
x = np.linspace(0, L, 1000) # Puntos en la dirección x
t = np.linspace(0, T, 1000) # Puntos en la dirección t
X, T = np.meshgrid(x, t)

def f(x):
    return np.maximum(0, 1 - 4 * np.abs(x - 1/2))

def fourier_coefficients(f, n, x):
    coefficients = np.zeros(n)
    for k in range(1,n+1): 
        coefficients[k-1] = 2 * np.trapz(f(x) * np.sin(k * np.pi * x), x)
    return coefficients.tolist()  # Convertir a lista
coef=fourier_coefficients(f, N, x)

# Solución usando los primeros N términos de la serie
def u(x, t, n):
    sum = 0
    for k in range(1, n+1):
        lambda_k = k* np.pi
        sum += coef[k-1] * np.sin(lambda_k * x) * np.exp(-lambda_k**2 * t)
    return sum

Z = u(X, T, N)

# Graficar

fig = plt.figure(figsize=(12, 8))  # Cambio aquí para hacer la figura más grande
ax = fig.add_subplot(111, projection='3d')
surf=ax.plot_surface(X, T, Z,cmap='viridis')

ax.set_xlabel('Posición x')
ax.set_ylabel('Tiempo t')
ax.set_zlabel('Temperatura T(x,t)')
ax.set_title('Evolución de la Temperatura en la Varilla')
fig.colorbar(surf, shrink=0.5, aspect=5)

plt.show()

​ Vemos como la temperatura sufre una gran variación en tiempos cercanos a 0, y rápidamente tiende a su estado estacionario, [math]v(x)=0[/math]. Para poder visualizarlo mejor mostraremos una grafica restringida a [math]t \in [0,0.05][/math] Al principio, u es nulo fuera del origen. Una vez que [math]t\gt0[/math] u se hace positivo en todas partes, lo que significa que la temperatura se difunde casi instantáneamente a lo largo del eje x, resultando en una velocidad de propagación infinita. \ref Terminamos el estudio de esta ecuación de difusión con el estudio del flujo, la inversa de la variación de la variable x de la solución. Para ello derivamos la solución que hemos obtenido. En términos de [math]C_k[/math], puesto que no son dependientes de x, obtenemos: [math]\frac{u(x,t)}{dx}=\sum_{k=1}^{\infty} C_k e^{(-k^2\pi ^2 t )}\cos(k \pi x)\frac{1}{k \pi}[/math] como flujo de la ecuación homogénea. Adjuntamos el código con el cual se ha calculado la expresión de forma numérica y su respectiva representación gráfica