Physical Context

We aim to analyze how heat disperses in a long metallic rod of length L. One end is positioned at x=0, while the opposite end is at x=L. Since the rod is significantly longer than it is wide, we can simplify our analysis to account only for variations in x and t. Given the specific heat capacity is known, determining the function T(x, t) will provide insight into the heat distribution.

The rod is insulated along its length, allowing heat transfer solely through its ends. Thus, the temperature profile is influenced by three factors:

1. The initial temperature distribution T(x, 0) — known as the initial condition.
2. The temperatures at the rod's ends, T(0, t) and T(L, t) — these are the boundary conditions.
3. A governing law that dictates heat transfer within the rod. The heat equation serves as a mathematical expression of this physical law.

When tasked with solving a PDE under a set of initial and boundary conditions, we refer to it as an initial boundary value problem (IBVP). To keep our analysis straightforward, we will focus on the IBVP for the heat equation with boundary conditions T(0, t) = T(L, t) = 0°C, known as homogeneous boundary conditions.

Deriving the Heat Equation

The heat equation can be derived from the conservation of energy principle: the rate of change of heat at a point in the bar correlates with the net heat flow into that point, adhering to the continuity equation. If we denote Q as the heat at each location and V as the vector field representing heat flow, we can express this relationship.

According to the Second Law of Thermodynamics, when two identical bodies are in thermal contact, heat flows from the hotter to the cooler body at a rate proportional to their temperature difference. Thus, V is proportional to the negative temperature gradient, leading to the equation V = -k∇T. In a one-dimensional scenario, this simplifies to V = (-k∂T/∂x)x, where x is the unit vector in the +x-direction. Additionally, since Q = ρcT, we derive the heat equation by substituting these variables.

Solving the Heat Equation

To progress, we must demonstrate the existence of a unique solution for the heat equation under physically relevant initial and boundary conditions. While a formal proof exceeds our current scope, we can use empirical reasoning. The laws of thermodynamics assert that regardless of the initial temperature distribution along the bar, the system will inevitably reach thermal equilibrium, governed by the heat equation.

This indicates that if f(x, t) and g(x, t) are distinct functions satisfying the same IBVP for the heat equation, they will share the same form. Given the linear nature of the heat equation, any linear combination of solutions will also be a solution.

Consider a function that satisfies the heat equation:

Where n is a positive integer greater than zero. This function also meets the boundary conditions, leading us to conclude that the general solution can be expressed as:

The challenge lies in determining the coefficients Aₙ to ensure this general solution aligns with the initial conditions. This process is known as a Fourier sine series expansion for the initial conditions, with Aₙ referred to as the Fourier coefficients.

Calculating the Fourier Coefficients

The initial condition T(x, 0) is a piecewise continuous function on the interval [0, L], reaching zero at the boundaries. Notably, the set of such functions forms a vector space under addition and scalar multiplication, which we denote as 𝕊ᴸ.

This vector space can be equipped with an inner product, allowing us to express any function f∈𝕊ᴸ as a linear combination of basis functions. For integers m, n > 0, the functions of the form sin(nπx/L) are orthonormal. Consequently, any function f∈𝕊ᴸ can be represented as:

The coefficients for this representation are provided by the Euler integral. To illustrate, let's derive the Fourier coefficients for a unit sawtooth pulse:

Clearly, saw(x)∈𝕊¹, yielding Fourier coefficients as follows:

The Fourier series expansion for the sawtooth wave is:

This animation depicts how the Fourier series converges to the sawtooth function as the number of sine terms increases.

The Gibbs phenomenon, an inherent error in the Fourier series of a discontinuous function, causes an overestimation of the function's value at points of discontinuity. Though the Gibbs phenomenon cannot be entirely eliminated, increasing the number of terms in the Fourier series minimizes the error to the point of discontinuity.

Thus, solving the homogeneous IBVP for the heat equation involves using the Euler integral to determine the Fourier coefficients.

Example: Uniform Temperature Distribution

Consider an insulated, one-meter-long bar with k/ρc = 0.1 m²/s (a simplified assumption for illustration) that starts at a uniform temperature of 100°C. Upon attaching cooling elements set to 0°C at both ends at t=0, the initial and boundary conditions are:

The Fourier coefficients calculated for this scenario are:

Let’s verify that this generates a Fourier series consistent with the initial and boundary conditions.

Now, let's visualize the solution:

And here’s a three-dimensional representation:

Example: Temperature "Spike"

Now imagine the bar begins at 0°C except for a 10-centimeter section in the middle, where the temperature is 100°C. With diffusivity set at 0.0075 m²/s, the Fourier coefficients can be calculated in a more general form:

The solution, in this case, is:

Although the bar will eventually stabilize to a uniform temperature of 0°C, it takes time for the temperature to equilibrate, so we’ll conclude the animation after five seconds. Here’s the 3D plot:

Example: Random Temperature Distribution

In another scenario, assume the rod, made of copper with a diffusivity of 1.11×10⁻⁴ m²/s, exhibits a random temperature distribution measured at 10-centimeter intervals as recorded in the following table:

For such cases, numerical integration is typically preferable for finding Fourier coefficients, as demonstrated using MATLAB. If interested, the code is available here on Pastebin.

Here’s an animation illustrating the solution:

After 240 seconds, the temperature changes very slowly, prompting us to conclude the simulation at this point. The most dynamic behavior occurs within the first 60 seconds, which is reflected in the accompanying 3D plot.

Concluding Thoughts

This concludes our exploration of the heat equation's fundamental principles. With a solid understanding of how to solve the heat equation for basic scenarios, you’re equipped to tackle more complex problems in future discussions.

A detailed explanation of the heat equation with practical examples.

An in-depth derivation of the heat equation from calculus, ideal for learners.

Copyright

All images, animations, and code in this article are original works. Feel free to use them elsewhere, provided you give appropriate credit. The MATLAB code linked should function on most modern machines, but I cannot guarantee universal compatibility.

UPDATE (July 10, 2019): The second article in this series is now available. Click here to read it.