Ver el código de «Ecuación del calor CCP»

{{ TrabajoED 
| Ecuación del calor CCP
| [[Categoría:EDP|EDP]] EDP
|[[Categoría:EDP25/26|Ecuación del calor CCP]] Curso 2025-26
| Claudia Pozurama Pliego<br>Covadonga Díaz García<br>Paula Blanco Díez
}}

[[Archivo:Poster calor CCP (1).png|center|800px]]


A continuación se muestran los códigos utilizados en el proyecto.

<syntaxhighlight lang="python">
#DATO INICIAL...............................................................
import numpy as np
import matplotlib.pyplot as plt

# Definimos los tramos
x1 = np.linspace(0, 1/3, 300)
x2 = np.linspace(1/3, 2/3, 300)
x3 = np.linspace(2/3, 1, 300)

# Figura
plt.figure(figsize=(9,5))

# Tramos de la función
plt.plot(x1, 10*np.ones_like(x1), color='red', linewidth=3)
plt.plot(x2, np.zeros_like(x2), color='blue', linewidth=3)
plt.plot(x3, 10*np.ones_like(x3), color='red', linewidth=3)

# Líneas de discontinuidad
plt.axvline(1/3, color='skyblue', linestyle='--', linewidth=1)
plt.axvline(2/3, color='skyblue', linestyle='--', linewidth=1)

# Texto
plt.text(0.15, 10.4, 'Zona caliente', color='red', fontsize=10)
plt.text(0.72, 10.4, 'Zona caliente', color='red', fontsize=10)
plt.text(0.45, 0.6, 'Zona fría', color='blue', fontsize=10)

# Ejes y título
plt.xlabel('$x$')
plt.ylabel('Temperatura')
plt.title('Dato inicial: distribución de temperatura')

# Límites
plt.xlim(0,1)
plt.ylim(-1,11)

# Cuadrícula
plt.grid(alpha=0.3)

plt.tight_layout()
plt.savefig('Dato inicial.png', dpi=300)
plt.show()

#EVOLUCIÓN TEMPERATURA...............................................................
import numpy as np
import matplotlib.pyplot as plt

# Mallado
x = np.linspace(0, 1, 1000)

# Dato inicial
def u0(x):
    return np.where((x >= 1/3) & (x <= 2/3), 0, 10)

# Restamos estado estacionario
def v0(x):
    return u0(x) - 10

# Número de términos de Fourier
N = 100

# Coeficientes (trapecio)
b = np.zeros(N)
for n in range(1, N+1):
    integrando = v0(x) * np.sin(n*np.pi*x)
    b[n-1] = 2 * np.trapz(integrando, x)

# Solución aproximada
def u(x, t):
    if t == 0:
        return u0(x)
    suma = np.zeros_like(x)
    for n in range(1, N+1):
        suma += b[n-1] * np.exp(-n**2 * np.pi**2 * t) * np.sin(n*np.pi*x)
    return 10 + suma

# Figura
plt.figure(figsize=(9,5))

# t=0 (DISCONTINUA)
x1 = np.linspace(0, 1/3, 300)
x2 = np.linspace(1/3, 2/3, 300)
x3 = np.linspace(2/3, 1, 300)

plt.plot(x1, 10*np.ones_like(x1), color='black', linewidth=3, label='t=0')
plt.plot(x2, np.zeros_like(x2), color='black', linewidth=3)
plt.plot(x3, 10*np.ones_like(x3), color='black', linewidth=3)

# t>0 (SUAVE)
times = [0.001, 0.01, 0.05, 0.1]
colors = ['blue', 'green', 'orange', 'red']

for t, c in zip(times, colors):
    plt.plot(x, u(x, t), color=c, linewidth=2.5, label=f't={t}')

# Líneas verticales
plt.axvline(1/3, linestyle='--', alpha=0.4)
plt.axvline(2/3, linestyle='--', alpha=0.4)

# Texto
plt.text(0.1, 10.4, 'Zona caliente', fontsize=10)
plt.text(0.38, 0.7, 'Zona fría inicial', fontsize=10)

# Ejes
plt.xlabel('x')
plt.ylabel('Temperatura')
plt.title('Evolución de la temperatura en la barra')

plt.xlim(0,1)
plt.ylim(-1,11)

plt.grid(alpha=0.3)
plt.legend()

plt.tight_layout()
plt.savefig('Evolucion temperatura.png', dpi=300)
plt.show()

# ZOOM EN LA ZONA CENTRAL ..............................................................

plt.figure(figsize=(6,4))

mask = (x >= 0.4) & (x <= 0.6)
x_zoom = x[mask]

# t=0
plt.plot(x_zoom, u0(x_zoom), color='black', linewidth=3, label='t=0')

# t pequeño
t_zoom = 0.001
plt.plot(x_zoom, u(x_zoom, t_zoom), color='blue', linewidth=2.5, label='t=0.001')

# Líneas verticales
plt.axvline(1/3, linestyle='--', alpha=0.3)
plt.axvline(2/3, linestyle='--', alpha=0.3)

# CLAVE: zoom fuerte en eje Y
plt.ylim(-0.005, 0.05)

# Etiquetas
plt.title('Zoom en la zona central')
plt.xlabel('x')
plt.ylabel('Temperatura')

plt.grid(alpha=0.3)
plt.legend()

plt.tight_layout()
plt.savefig('zoom_centro.png', dpi=300)
plt.show()

#ECUACIÓN DEL CALOR ......................................................................
import numpy as np
import matplotlib.pyplot as plt

# Mallado
x = np.linspace(0, 1, 1000)

# =======================
# MODELO 1: ECUACIÓN DEL CALOR
# =======================

def u0(x):
    return np.where((x >= 1/3) & (x <= 2/3), 0, 10)

def v0(x):
    return u0(x) - 10

N = 100

b = np.zeros(N)
for n in range(1, N+1):
    integrando = v0(x) * np.sin(n*np.pi*x)
    b[n-1] = 2 * np.trapz(integrando, x)

def u(x, t):
    if t == 0:
        return u0(x)
    suma = np.zeros_like(x)
    for n in range(1, N+1):
        suma += b[n-1] * np.exp(-n**2 * np.pi**2 * t) * np.sin(n*np.pi*x)
    return suma + 10

# =======================
# MODELO 2: VELOCIDAD FINITA
# =======================

def u_finite_speed(x, t, c=0.18, delta=0.025):
    left_front = 1/3 - c*t
    right_front = 2/3 + c*t

    u = np.full_like(x, 10.0)

    inside = (x >= left_front) & (x <= right_front)
    u[inside] = 0.0

    left_transition = (x > left_front - delta) & (x < left_front)
    u[left_transition] = 10 * (left_front - x[left_transition]) / delta

    right_transition = (x > right_front) & (x < right_front + delta)
    u[right_transition] = 10 * (x[right_transition] - right_front) / delta

    return np.clip(u, 0, 10)

# =======================
# PLOT DOBLE
# =======================

fig, axs = plt.subplots(1, 2, figsize=(12,5))

# Colores consistentes
colors = ['tab:blue', 'tab:orange', 'tab:green', 'tab:red']
times = [0, 0.01, 0.03, 0.06]

# --------- IZQUIERDA (ecuación del calor)
ax = axs[0]

# t = 0 por tramos
x1 = np.linspace(0, 1/3, 300)
x2 = np.linspace(1/3, 2/3, 300)
x3 = np.linspace(2/3, 1, 300)

ax.plot(x1, 10*np.ones_like(x1), color=colors[0], linewidth=2, label='t = 0')
ax.plot(x2, np.zeros_like(x2), color=colors[0], linewidth=2)
ax.plot(x3, 10*np.ones_like(x3), color=colors[0], linewidth=2)

# t > 0
for t, c in zip(times[1:], colors[1:]):
    ax.plot(x, u(x, t), color=c, linewidth=2, label=f't = {t}')

ax.axvline(1/3, linestyle='--', alpha=0.3)
ax.axvline(2/3, linestyle='--', alpha=0.3)

ax.set_title("Ecuación del calor (difusión instantánea)")
ax.set_xlabel("x")
ax.set_ylabel("Temperatura")
ax.grid(alpha=0.3)
ax.legend()

# --------- DERECHA (velocidad finita)
ax = axs[1]

for t, c in zip(times, colors):
    ax.plot(x, u_finite_speed(x, t), color=c, linewidth=2, label=f't = {t}')

ax.axvline(1/3, linestyle='--', alpha=0.3)
ax.axvline(2/3, linestyle='--', alpha=0.3)

ax.set_title("Modelo esquemático (velocidad finita)")
ax.set_xlabel("x")
ax.set_ylabel("Temperatura")
ax.grid(alpha=0.3)
ax.legend()

# Ajuste final
plt.tight_layout()
plt.savefig('ec_calor.png', dpi=300)
plt.show()

#BOLLO ...........................................................
import numpy as np
import matplotlib.pyplot as plt

x = np.linspace(0, 1, 1000)

def bollo(x, t):
    center = 0.5
    dist = np.abs(x - center)
    
    c = 0.25
    front = c * t
    
    u = np.zeros_like(x)
    
    # zonas que ya han sido alcanzadas por el calor
    heated = dist >= (0.5 - front)
    
    # crecimiento suave desde los bordes
    u[heated] = 10 * (1 - np.exp(-5 * (dist[heated] - (0.5 - front))))
    
    return np.clip(u, 0, 10)

times = [0, 0.5, 1, 2]

plt.figure(figsize=(6,4))

for t in times:
    plt.plot(x, bollo(x, t), label=f"t={t}")

plt.xlabel("posición en el bollo")
plt.ylabel("temperatura")
plt.title("Calentamiento progresivo de un bollo")
plt.legend()
plt.grid(alpha=0.3)
plt.savefig('cal_prog.png', dpi=300)
plt.show()
</syntaxhighlight>