update exercise 2: enhance heat equation simulation in Jupyter notebook, improve matrix construction and add temperature profile visualization

This commit is contained in:
Félix MARQUET
2025-04-01 11:58:58 +02:00
parent 2b0bdf89cf
commit fe9af73a23
2 changed files with 119 additions and 35 deletions

View File

@@ -1,18 +1,83 @@
# Eq en régime permanent
# d2T / dx2 = 0
# T(x =0) = T0
# T(x = 1) = T1
# T(i)i E [0, N+1]
# x(i+1) - x(i) = deltaX = 1/N
# d2T/dx2 = (T(i-1) - 2*T(i) + T(i+1))/deltaX**2
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.cm as cm
# Question 1: Résoudre cette équation pour les conditions aux bord
# d2T / dx2 => dT / dx = A => T(x) = Ax + b
# Ici b = T0 d'où T(x) = Ax + T0
# T1 - T0 / 1 - 0 = A = T1 - T0 donc T(x) = (T1 - T0)x + T0
# Question 2: Ecrire l'équation de la chaleur discrétisée pour i E [1, N-2], puis pour i = 0 et i = N
# d2T/dx2 = 0
# T(i+1) - 2*T(i) + T(i-1) = 0
# Pour i = 0
# -2T0 + T(1) = -T0-
# Pour i = N
# Constantes et paramètres
N = 1000
delta_x = 1 / N
delta_x2 = delta_x**2
D = 1
delta_t = 0.9 * (delta_x2 / 2 * D)
T0 = 0
T1 = 0
if delta_t < (delta_x2 / 2 * D):
M = np.zeros((N + 1, N + 1))
# Remplissage de la matrice M
for i in range(N + 1):
if i == 0:
M[i, i] = 1 - 2 * delta_t * D / delta_x2
M[i, i + 1] = delta_t * D / delta_x2
elif i == N:
M[i, i - 1] = delta_t * D / delta_x2
M[i, i] = 1 - 2 * delta_t * D / delta_x2
else:
M[i, i - 1] = delta_t * D / delta_x2
M[i, i] = 1 - 2 * delta_t * D / delta_x2
M[i, i + 1] = delta_t * D / delta_x2
# Vecteur spatial pour tracer
x = np.linspace(0, 1, N + 1)
# Condition initiale T₀(x) = sin(πx)
T = np.sin(np.pi * x)
T[0] = T0
T[N] = T1
# Temps de simulation
temps_final = 0.1
nb_iterations = int(temps_final / delta_t)
# Nombre de courbes souhaité
nb_courbes = 100000
# Sauvegarder plus de profils à différents temps
profils = [T.copy()]
iterations_sauvegarde = np.linspace(0, nb_iterations, nb_courbes, dtype=int)
temps_sauvegarde = [0] + [i * delta_t for i in iterations_sauvegarde[1:]]
# Simulation de l'évolution temporelle
for n in range(nb_iterations):
T = M @ T
if n in iterations_sauvegarde:
profils.append(T.copy())
# Tracer les profils de température avec un dégradé de couleurs
plt.figure(figsize=(12, 8))
colormap = cm.viridis
# Affichage des courbes avec un dégradé de couleurs
for i, t in enumerate(temps_sauvegarde):
if i < len(profils):
couleur = colormap(i / len(temps_sauvegarde))
plt.plot(x, profils[i], color=couleur, alpha=0.7)
# Afficher seulement quelques légendes pour éviter l'encombrement
indices_legende = np.linspace(0, len(temps_sauvegarde) - 1, 10, dtype=int)
for i in indices_legende:
if i < len(profils):
couleur = colormap(i / len(temps_sauvegarde))
plt.plot([], [], color=couleur, label=f"t = {temps_sauvegarde[i]:.4f}s")
plt.xlabel("Position x")
plt.ylabel("Température T")
plt.title("Évolution de la température au cours du temps (50 courbes)")
plt.legend(loc="upper right")
plt.grid(True)
plt.colorbar(plt.cm.ScalarMappable(cmap=colormap), label="Temps (s)", ax=plt.gca())
plt.show()
else:
print("Constante incorrect")
print(f"delta_t = {delta_t}")
print(f"(delta_x2 / 2 * D) = {(delta_x2 / 2 * D)}")

View File

@@ -250,22 +250,30 @@
]
},
{
"metadata": {},
"metadata": {
"jupyter": {
"is_executing": true
},
"ExecuteTime": {
"start_time": "2025-04-01T09:53:17.236001Z"
}
},
"cell_type": "code",
"source": [
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"import matplotlib.cm as cm\n",
"\n",
"# Constantes et paramètres\n",
"N = 10\n",
"delta_t = (delta_x2 / 2 * D)\n",
"delta_x2 = 1/N\n",
"N = 1000\n",
"delta_x = 1/N\n",
"delta_x2 = delta_x**2\n",
"D = 1\n",
"delta_t = 0.9 * (delta_x2 / 2 * D)\n",
"T0 = 0\n",
"T1 = 0\n",
"\n",
"if delta_t < (delta_x2 / 2 * D):\n",
"\n",
" M = np.zeros((N + 1, N + 1))\n",
"\n",
" # Remplissage de la matrice M\n",
@@ -281,9 +289,6 @@
" M[i, i] = (1 - 2 * delta_t * D / delta_x2)\n",
" M[i, i + 1] = delta_t * D / delta_x2\n",
"\n",
" # Affichage de la matrice M\n",
" print(f\"Matrice M: {M}\")\n",
"\n",
" # Vecteur spatial pour tracer\n",
" x = np.linspace(0, 1, N+1)\n",
"\n",
@@ -296,30 +301,44 @@
" temps_final = 0.1\n",
" nb_iterations = int(temps_final / delta_t)\n",
"\n",
" # Sauvegarder quelques profils à différents temps\n",
" # Nombre de courbes souhaité\n",
" nb_courbes = 100000\n",
"\n",
" # Sauvegarder plus de profils à différents temps\n",
" profils = [T.copy()]\n",
" temps_sauvegarde = [0, 0.01, 0.02, 0.05, 0.1]\n",
" indices_sauvegarde = [int(t/delta_t) for t in temps_sauvegarde[1:]]\n",
" iterations_sauvegarde = np.linspace(0, nb_iterations, nb_courbes, dtype=int)\n",
" temps_sauvegarde = [0] + [i * delta_t for i in iterations_sauvegarde[1:]]\n",
"\n",
" # Simulation de l'évolution temporelle\n",
" for n in range(nb_iterations):\n",
" T = M @ T\n",
" if n+1 in indices_sauvegarde:\n",
" if n in iterations_sauvegarde:\n",
" profils.append(T.copy())\n",
"\n",
" print(f\"Solution profils: {profils}\")\n",
" # Tracer les profils de température avec un dégradé de couleurs\n",
" plt.figure(figsize=(12, 8))\n",
" colormap = cm.viridis\n",
"\n",
" # Tracer les profils de température\n",
" plt.figure(figsize=(10, 6))\n",
" # Affichage des courbes avec un dégradé de couleurs\n",
" for i, t in enumerate(temps_sauvegarde):\n",
" if i < len(profils):\n",
" plt.plot(x, profils[i], label=f't = {t}s')\n",
" couleur = colormap(i / len(temps_sauvegarde))\n",
" plt.plot(x, profils[i], color=couleur, alpha=0.7)\n",
"\n",
" # Afficher seulement quelques légendes pour éviter l'encombrement\n",
" indices_legende = np.linspace(0, len(temps_sauvegarde)-1, 10, dtype=int)\n",
" for i in indices_legende:\n",
" if i < len(profils):\n",
" couleur = colormap(i / len(temps_sauvegarde))\n",
" plt.plot([], [], color=couleur, label=f't = {temps_sauvegarde[i]:.4f}s')\n",
"\n",
" plt.xlabel('Position x')\n",
" plt.ylabel('Température T')\n",
" plt.title('Évolution de la température au cours du temps')\n",
" plt.legend()\n",
" plt.title('Évolution de la température au cours du temps (50 courbes)')\n",
" plt.legend(loc='upper right')\n",
" plt.grid(True)\n",
" plt.colorbar(plt.cm.ScalarMappable(cmap=colormap),\n",
" label='Temps (s)', ax=plt.gca())\n",
" plt.show()\n",
"else:\n",
" print(\"Constante incorrect\")\n",