simulate heat equation in python

main.py
import numpy as np
import matplotlib.pyplot as plt

# Parameters
L = 1.0  # Length of the rod
T = 0.1  # Total simulation time
k = 1.0  # Thermal diffusivity
dx = 0.01  # Grid spacing
dt = 0.0001  # Time step size
r = k * dt / dx**2  # Stability parameter

# Grid
x = np.arange(0, L + dx, dx)
N = len(x)

# Initial condition
T0 = np.zeros(N)
T0[int(0.4*N):int(0.6*N)] = 1.0

# Main loop
T = T0.copy()
t = 0
while t < T:
    # Compute new temperature
    T[1:-1] += r * (T[2:] - 2*T[1:-1] + T[:-2])

    # Boundary conditions
    T[0] = 0
    T[-1] = 0

    t += dt

# Plot results
plt.plot(x, T0, label='Initial')
plt.plot(x, T, label='Final')
plt.legend()
plt.show()
668 chars
38 lines