Heat Equation in Python

Heat Equation in Python |  In this article, we will discuss the heat equation and how we can solve it using Python programming. So, firstly we must know what the heat equation is?

What is Heat Equation?

The heat equation is a partial differential equation that describes how the temperature of material changes over time. It is given by:- ∂ T / ∂ t = α ∇ ² T

Where T is the temperature, t is time, α is the thermal diffusivity, and ∇² is the Laplacian operator. This equation states that the rate of change of temperature with respect to time is proportional to the second derivative of temperature with respect to space.

In Python, the heat equation can be solved numerically using various methods, such as the finite difference method, the finite element method, or the spectral method. The basic idea is to discretize the temperature field into a grid of points, and then iteratively update the temperature at each point based on the temperature at neighboring points and the time step. The resulting numerical solution can be visualized using various plotting tools in Python, such as Matplotlib.

1D Heat Equation In Python

Here’s an example implementation of the 1D heat equation in Python using the finite difference method:

import numpy as np
import matplotlib.pyplot as plt

# Define parameters
L = 1.0  # Length of rod
nx = 100  # Number of grid points
dx = L / (nx - 1)  # Grid spacing
nt = 1000  # Number of time steps
dt = 0.001  # Time step size
alpha = 1.0  # Thermal diffusivity
T0 = np.zeros(nx)  # Initial temperature profile
T0[0] = 1.0  # Set left boundary temperature to 1.0

# Set up a coefficient matrix
A = np.zeros((nx, nx))
for i in range(1, nx - 1):
    A[i, i - 1] = -alpha * dt / dx**2
    A[i, i] = 1 + 2 * alpha * dt / dx**2
    A[i, i + 1] = -alpha * dt / dx**2

# Iterate over time steps
T = T0.copy()
for n in range(nt):
    T = np.linalg.solve(A, T)

# Plot final temperature profile
x = np.linspace(0.0, L, nx)
plt.plot(x, T)
plt.xlabel("x")
plt.ylabel("Temperature")
plt.show()

This code solves the 1D heat equation on a rod using the finite difference method, with a time step size of 0.001 and a thermal diffusivity of 1.0. The initial temperature profile is set to zero everywhere except at the left boundary, which is set to 1.0. The code then iterates over time steps to solve for the temperature profile at each time. The final temperature profile is then plotted using Matplotlib’s ‘plot’ function. Note that the boundary conditions are implemented implicitly in the coefficient matrix ‘A’.

2D Heat Equation In Python

Here’s an example implementation of the 2D heat equation in Python using the finite difference method:-

import numpy as np
import matplotlib.pyplot as plt

# Define parameters
Lx = 1.0   # Length of x-axis
Ly = 1.0   # Length of y-axis
nx = 50    # Number of grid points in the x-direction
ny = 50    # Number of grid points in the y-direction
dx = Lx / (nx - 1)  # Grid spacing in the x-direction
dy = Ly / (ny - 1)  # Grid spacing in the y-direction
nt = 1000  # Number of time steps
dt = 0.001  # Time step size
alpha = 1.0  # Thermal diffusivity
T0 = np.zeros((nx, ny))  # Initial temperature profile
T0[0, :] = 1.0  # Set left boundary temperature to 1.0

# Create grid
x = np.linspace(0.0, Lx, nx)
y = np.linspace(0.0, Ly, ny)
X, Y = np.meshgrid(x, y)

# Set up a coefficient matrix
Ax = np.zeros((nx, nx))
Ay = np.zeros((ny, ny))
for i in range(1, nx - 1):
    Ax[i, i - 1] = -alpha * dt / dx**2
    Ax[i, i] = 1 + 2 * alpha * dt / dx**2
    Ax[i, i + 1] = -alpha * dt / dx**2
for i in range(1, ny - 1):
    Ay[i, i - 1] = -alpha * dt / dy**2
    Ay[i, i] = 1 + 2 * alpha * dt / dy**2
    Ay[i, i + 1] = -alpha * dt / dy**2

# Iterate over time steps
T = T0.copy()
for n in range(nt):
    # Update x-direction
    for i in range(ny):
        T[:, i] = np.linalg.solve(Ax, T[:, i])
    # Update y-direction
    for i in range(nx):
        T[i, :] = np.linalg.solve(Ay, T[i, :])

# Plot final temperature profile
plt.pcolormesh(X, Y, T)
plt.xlabel('x')
plt.ylabel('y')
plt.colorbar()
plt.show()

This code solves the 2D heat equation on a square domain using the finite difference method, with a time step size of 0.001 and a thermal diffusivity of 1.0. The initial temperature profile is set to zero everywhere except at the left boundary, which is set to 1.0. The code then iterates over time steps to solve for the temperature profile at each time. The final temperature profile is then plotted using Matplotlib’s ‘pcolormesh’ function, which creates a 2D colored plot of the temperature field. Note that the boundary conditions are implemented implicitly in the coefficient matrices ‘Ax’ and ‘Ay’.

If you enjoyed this post, share it with your friends. Do you want to share more information about the topic discussed above or do you find anything incorrect? Let us know in the comments. Thank you!