Section 9.1 : The Heat Equation
Before we get into actually solving partial differential equations and before we even start discussing the method of separation of variables we want to spend a little bit of time talking about the two main partial differential equations that we’ll be solving later on in the chapter. We’ll look at the first one in this section and the second one in the next section.
The first partial differential equation that we’ll be looking at once we get started with solving will be the heat equation, which governs the temperature distribution in an object. We are going to give several forms of the heat equation for reference purposes, but we will only be really solving one of them.
We will start out by considering the temperature in a 1-D bar of length \(L\). What this means is that we are going to assume that the bar starts off at \(x = 0\) and ends when we reach \(x = L\). We are also going to so assume that at any location, \(x\) the temperature will be constant at every point in the cross section at that \(x\). In other words, temperature will only vary in \(x\) and we can hence consider the bar to be a 1-D bar. Note that with this assumption the actual shape of the cross section (i.e. circular, rectangular, etc.) doesn’t matter.
Note that the 1-D assumption is actually not all that bad of an assumption as it might seem at first glance. If we assume that the lateral surface of the bar is perfectly insulated (i.e. no heat can flow through the lateral surface) then the only way heat can enter or leave the bar as at either end. This means that heat can only flow from left to right or right to left and thus creating a 1-D temperature distribution.
The assumption of the lateral surfaces being perfectly insulated is of course impossible, but it is possible to put enough insulation on the lateral surfaces that there will be very little heat flow through them and so, at least for a time, we can consider the lateral surfaces to be perfectly insulated.
Okay, let’s now get some definitions out of the way before we write down the first form of the heat equation.
\[\begin{align*}&u\left( {x,t} \right) = {\mbox{Temperature at any point }}x{\mbox{ and any time }}t\\ & c\left( x \right) = {\mbox{Specific Heat}}\\ & \rho \left( x \right) = {\mbox{Mass Density}}\\ & \varphi \left( {x,t} \right) = {\mbox{Heat Flux}}\\ & Q\left( {x,t} \right) = {\mbox{Heat energy generated per unit volume per unit time}}\end{align*}\]We should probably make a couple of comments about some of these quantities before proceeding.
The specific heat, \(c\left( x \right) > 0\), of a material is the amount of heat energy that it takes to raise one unit of mass of the material by one unit of temperature. As indicated we are going to assume, at least initially, that the specific heat may not be uniform throughout the bar. Note as well that in practice the specific heat depends upon the temperature. However, this will generally only be an issue for large temperature differences (which in turn depends on the material the bar is made out of) and so we’re going to assume for the purposes of this discussion that the temperature differences are not large enough to affect our solution.
The mass density, \(\rho \left( x \right)\), is the mass per unit volume of the material. As with the specific heat we’re going to initially assume that the mass density may not be uniform throughout the bar.
The heat flux, \(\varphi \left( {x,t} \right)\), is the amount of thermal energy that flows to the right per unit surface area per unit time. The “flows to the right” bit simply tells us that if \(\varphi \left( {x,t} \right) > 0\) for some \(x\) and \(t\) then the heat is flowing to the right at that point and time. Likewise, if \(\varphi \left( {x,t} \right) < 0\) then the heat will be flowing to the left at that point and time.
The final quantity we defined above is \(Q\left( {x,t} \right)\) and this is used to represent any external sources or sinks (i.e. heat energy taken out of the system) of heat energy. If \(Q\left( {x,t} \right) > 0\) then heat energy is being added to the system at that location and time and if \(Q\left( {x,t} \right) < 0\) then heat energy is being removed from the system at that location and time.
With these quantities the heat equation is,
\[\begin{equation}c\left( x \right)\rho \left( x \right)\frac{{\partial u}}{{\partial t}} = - \frac{{\partial \varphi }}{{\partial x}} + Q\left( {x,t} \right)\label{eq:eq1}\end{equation}\]While this is a nice form of the heat equation it is not actually something we can solve. In this form there are two unknown functions, \(u\) and \(\varphi \), and so we need to get rid of one of them. With Fourier’s law we can easily remove the heat flux from this equation.
Fourier’s law states that,
\[\varphi \left( {x,t} \right) = - {K_0}\left( x \right)\frac{{\partial u}}{{\partial x}}\]where \({K_0}\left( x \right) > 0\) is the thermal conductivity of the material and measures the ability of a given material to conduct heat. The better a material can conduct heat the larger \({K_0}\left( x \right)\) will be. As noted the thermal conductivity can vary with the location in the bar. Also, much like the specific heat the thermal conductivity can vary with temperature, but we will assume that the total temperature change is not so great that this will be an issue and so we will assume for the purposes here that the thermal conductivity will not vary with temperature.
Fourier’s law does a very good job of modeling what we know to be true about heat flow. First, we know that if the temperature in a region is constant, i.e.\(\frac{{\partial u}}{{\partial x}} = 0\), then there is no heat flow.
Next, we know that if there is a temperature difference in a region we know the heat will flow from the hot portion to the cold portion of the region. For example, if it is hotter to the right then we know that the heat should flow to the left. When it is hotter to the right then we also know that \(\frac{{\partial u}}{{\partial x}} > 0\) (i.e. the temperature increases as we move to the right) and so we’ll have \(\varphi < 0\) and so the heat will flow to the left as it should. Likewise, if \(\frac{{\partial u}}{{\partial x}} < 0\) (i.e. it is hotter to the left) then we’ll have \[\varphi > 0\] and heat will flow to the right as it should.
Finally, the greater the temperature difference in a region (i.e. the larger \(\frac{{\partial u}}{{\partial x}}\) is) then the greater the heat flow.
So, if we plug Fourier’s law into \(\eqref{eq:eq1}\), we get the following form of the heat equation,
\[\begin{equation}c\left( x \right)\rho \left( x \right)\frac{{\partial u}}{{\partial t}} = \frac{\partial }{{\partial x}}\left( {{K_0}\left( x \right)\frac{{\partial u}}{{\partial x}}} \right) + Q\left( {x,t} \right)\label{eq:eq2}\end{equation}\]Note that we factored the minus sign out of the derivative to cancel against the minus sign that was already there. We cannot however, factor the thermal conductivity out of the derivative since it is a function of \(x\) and the derivative is with respect to \(x\).
Solving \(\eqref{eq:eq2}\) is quite difficult due to the non uniform nature of the thermal properties and the mass density. So, let’s now assume that these properties are all constant, i.e.,
\[c\left( x \right) = c\hspace{0.25in}\rho \left( x \right) = \rho \hspace{0.25in}{K_0}\left( x \right) = {K_0}\]where \(c\), \(\rho \) and \({K_0}\) are now all fixed quantities. In this case we generally say that the material in the bar is uniform. Under these assumptions the heat equation becomes,
\[\begin{equation}c\rho \frac{{\partial u}}{{\partial t}} = {K_0}\frac{{{\partial ^2}u}}{{\partial {x^2}}} + Q\left( {x,t} \right)\label{eq:eq3}\end{equation}\]For a final simplification to the heat equation let’s divide both sides by \(c\rho \) and define the thermal diffusivity to be,
\[k = \frac{{{K_0}}}{{c\rho }}\] \[\begin{equation}\frac{{\partial u}}{{\partial t}} = k\frac{{{\partial ^2}u}}{{\partial {x^2}}} + \frac{{Q\left( {x,t} \right)}}{{c\rho }}\label{eq:eq4}\end{equation}\]To most people this is what they mean when they talk about the heat equation and in fact it will be the equation that we’ll be solving. Well, actually we’ll be solving \(\eqref{eq:eq4}\) with no external sources, i.e.\(Q\left( {x,t} \right) = 0\), but we’ll be considering this form when we start discussing separation of variables in a couple of sections. We’ll only drop the sources term when we actually start solving the heat equation.
Now that we’ve got the 1-D heat equation taken care of we need to move into the initial and boundary conditions we’ll also need in order to solve the problem. If you go back to any of our solutions of ordinary differential equations that we’ve done in previous sections you can see that the number of conditions required always matched the highest order of the derivative in the equation.
In partial differential equations the same idea holds except now we have to pay attention to the variable we’re differentiating with respect to as well. So, for the heat equation we’ve got a first order time derivative and so we’ll need one initial condition and a second order spatial derivative and so we’ll need two boundary conditions.
The initial condition that we’ll use here is,
\[u\left( {x,0} \right) = f\left( x \right)\]and we don’t really need to say much about it here other than to note that this just tells us what the initial temperature distribution in the bar is.
The boundary conditions will tell us something about what the temperature and/or heat flow is doing at the boundaries of the bar. There are four of them that are fairly common boundary conditions.
The first type of boundary conditions that we can have would be the prescribed temperature boundary conditions, also called Dirichlet conditions. The prescribed temperature boundary conditions are,
\[u\left( {0,t} \right) = {g_1}\left( t \right)\hspace{0.25in}\hspace{0.25in}\hspace{0.25in}u\left( {L,t} \right) = {g_2}\left( t \right)\]The next type of boundary conditions are prescribed heat flux, also called Neumann conditions. Using Fourier’s law these can be written as,
\[ - {K_0}\left( 0 \right)\frac{{\partial u}}{{\partial x}}\left( {0,t} \right) = {\varphi _1}\left( t \right)\hspace{0.25in}\hspace{0.25in}\hspace{0.25in} - {K_0}\left( L \right)\frac{{\partial u}}{{\partial x}}\left( {L,t} \right) = {\varphi _2}\left( t \right)\]If either of the boundaries are perfectly insulated, i.e. there is no heat flow out of them then these boundary conditions reduce to,
\[\frac{{\partial u}}{{\partial x}}\left( {0,t} \right) = 0\hspace{0.25in}\hspace{0.25in}\hspace{0.25in}\frac{{\partial u}}{{\partial x}}\left( {L,t} \right) = 0\]and note that we will often just call these particular boundary conditions insulated boundaries and drop the “perfectly” part.
The third type of boundary conditions use Newton’s law of cooling and are sometimes called Robins conditions. These are usually used when the bar is in a moving fluid and note we can consider air to be a fluid for this purpose.
Here are the equations for this kind of boundary condition.
\[ - {K_0}\left( 0 \right)\frac{{\partial u}}{{\partial x}}\left( {0,t} \right) = - H\left[ {u\left( {0,t} \right) - {g_1}\left( t \right)} \right]\hspace{0.25in} - {K_0}\left( L \right)\frac{{\partial u}}{{\partial x}}\left( {L,t} \right) = H\left[ {u\left( {L,t} \right) - {g_2}\left( t \right)} \right]\]where \(H\) is a positive quantity that is experimentally determined and \({g_1}\left( t \right)\) and \({g_2}\left( t \right)\) give the temperature of the surrounding fluid at the respective boundaries.
Note that the two conditions do vary slightly depending on which boundary we are at. At \(x = 0\) we have a minus sign on the right side while we don’t at \(x = L\). To see why this is let’s first assume that at \(x = 0\) we have \(u\left( {0,t} \right) > {g_1}\left( t \right)\). In other words, the bar is hotter than the surrounding fluid and so at \(x = 0\) the heat flow (as given by the left side of the equation) must be to the left, or negative since the heat will flow from the hotter bar into the cooler surrounding liquid. If the heat flow is negative then we need to have a minus sign on the right side of the equation to make sure that it has the proper sign.
If the bar is cooler than the surrounding fluid at \(x = 0\), i.e. \(u\left( {0,t} \right) < {g_1}\left( t \right)\) we can make a similar argument to justify the minus sign. We’ll leave it to you to verify this.
If we now look at the other end, \(x = L\), and again assume that the bar is hotter than the surrounding fluid or, \(u\left( {L,t} \right) > {g_2}\left( t \right)\). In this case the heat flow must be to the right, or be positive, and so in this case we can’t have a minus sign. Finally, we’ll again leave it to you to verify that we can’t have the minus sign at \(x = L\) is the bar is cooler than the surrounding fluid as well.
Note that we are not actually going to be looking at any of these kinds of boundary conditions here. These types of boundary conditions tend to lead to boundary value problems such as Example 5 in the Eigenvalues and Eigenfunctions section of the previous chapter. As we saw in that example it is often very difficult to get our hands on the eigenvalues and as we’ll eventually see we will need them.
It is important to note at this point that we can also mix and match these boundary conditions so to speak. There is nothing wrong with having a prescribed temperature at one boundary and a prescribed flux at the other boundary for example so don’t always expect the same boundary condition to show up at both ends. This warning is more important that it might seem at this point because once we get into solving the heat equation we are going to have the same kind of condition on each end to simplify the problem somewhat.
The final type of boundary conditions that we’ll need here are periodic boundary conditions. Periodic boundary conditions are,
\[u\left( { - L,t} \right) = u\left( {L,t} \right)\hspace{0.25in}\hspace{0.25in}\frac{{\partial u}}{{\partial x}}\left( { - L,t} \right) = \frac{{\partial u}}{{\partial x}}\left( {L,t} \right)\]Note that for these kinds of boundary conditions the left boundary tends to be \(x = - L\) instead of \(x = 0\) as we were using in the previous types of boundary conditions. The periodic boundary conditions will arise very naturally from a couple of particular geometries that we’ll be looking at down the road.
We will now close out this section with a quick look at the 2-D and 3-D version of the heat equation. However, before we jump into that we need to introduce a little bit of notation first.
The del operator is defined to be,
\[\nabla = \frac{\partial }{{\partial x}}\vec i + \frac{\partial }{{\partial y}}\vec j\hspace{0.25in}\hspace{0.25in}\hspace{0.25in}\nabla = \frac{\partial }{{\partial x}}\vec i + \frac{\partial }{{\partial y}}\vec j + \frac{\partial }{{\partial z}}\vec k\]depending on whether we are in 2 or 3 dimensions. Think of the del operator as a function that takes functions as arguments (instead of numbers as we’re used to). Whatever function we “plug” into the operator gets put into the partial derivatives.
So, for example in 3-D we would have,
\[\nabla f = \frac{{\partial f}}{{\partial x}}\vec i + \frac{{\partial f}}{{\partial y}}\vec j + \frac{{\partial f}}{{\partial z}}\vec k\]This of course is also the gradient of the function \(f\left( {x,y,z} \right)\).
The del operator also allows us to quickly write down the divergence of a function. So, again using 3-D as an example the divergence of \(f\left( {x,y,z} \right)\) can be written as the dot product of the del operator and the function. Or,
\[\nabla \centerdot f = \frac{{\partial f}}{{\partial x}} + \frac{{\partial f}}{{\partial y}} + \frac{{\partial f}}{{\partial z}}\]Finally, we will also see the following show up in our work,
\[\nabla \centerdot \left( {\nabla f} \right) = \frac{\partial }{{\partial x}}\left( {\frac{{\partial f}}{{\partial x}}} \right) + \frac{\partial }{{\partial y}}\left( {\frac{{\partial f}}{{\partial y}}} \right) + \frac{\partial }{{\partial z}}\left( {\frac{{\partial f}}{{\partial k}}} \right) = \frac{{{\partial ^2}f}}{{\partial {x^2}}} + \frac{{{\partial ^2}f}}{{\partial {y^2}}} + \frac{{{\partial ^2}f}}{{\partial {z^2}}}\] \[{\nabla ^2}f = \frac{{{\partial ^2}f}}{{\partial {x^2}}} + \frac{{{\partial ^2}f}}{{\partial {y^2}}} + \frac{{{\partial ^2}f}}{{\partial {z^2}}}\]and is called the Laplacian. The 2-D version of course simply doesn’t have the third term.
Okay, we can now look into the 2-D and 3-D version of the heat equation and where ever the del operator and or Laplacian appears assume that it is the appropriate dimensional version.
The higher dimensional version of \(\eqref{eq:eq1}\) is,
\[\begin{equation}c\rho \frac{{\partial u}}{{\partial t}} = - \nabla \centerdot \,\varphi + Q\label{eq:eq5}\end{equation}\]and note that the specific heat, \(c\), and mass density, \(\rho \), are may not be uniform and so may be functions of the spatial variables. Likewise, the external sources term, \(Q\), may also be a function of both the spatial variables and time.
Next, the higher dimensional version of Fourier’s law is,
\[\varphi = - {K_0}\nabla u\]where the thermal conductivity, \({K_0}\), is again assumed to be a function of the spatial variables.
If we plug this into \(\eqref{eq:eq5}\) we get the heat equation for a non uniform bar (i.e. the thermal properties may be functions of the spatial variables) with external sources/sinks,
\[c\rho \frac{{\partial u}}{{\partial t}} = \nabla \centerdot \,\left( {{K_0}\nabla u} \right) + Q\]If we now assume that the specific heat, mass density and thermal conductivity are constant (i.e. the bar is uniform) the heat equation becomes,
\[\begin{equation}\frac{{\partial u}}{{\partial t}} = k{\nabla ^2}u + \frac{Q}{{cp}}\label{eq:eq6}\end{equation}\]where we divided both sides by \(c\rho \) to get the thermal diffusivity, \(k\) in front of the Laplacian.
The initial condition for the 2-D or 3-D heat equation is,
\[u\left( {x,y,t} \right) = f\left( {x,y} \right)\hspace{0.25in}\hspace{0.25in}{\mbox{or}}\hspace{0.25in}\hspace{0.25in}u\left( {x,y,z,t} \right) = f\left( {x,y,z} \right)\]depending upon the dimension we’re in.
The prescribed temperature boundary condition becomes,
\[u\left( {x,y,t} \right) = T\left( {x,y,t} \right)\hspace{0.25in}\hspace{0.25in}{\mbox{or}}\hspace{0.25in}\hspace{0.25in}u\left( {x,y,z,t} \right) = T\left( {x,y,z,t} \right)\]where \(\left( {x,y} \right)\) or \(\left( {x,y,z} \right)\), depending upon the dimension we’re in, will range over the portion of the boundary in which we are prescribing the temperature.
The prescribed heat flux condition becomes,
\[ - {K_0}\nabla u\,\centerdot \,\vec n = \varphi \left( t \right)\]where the left side is only being evaluated at points along the boundary and \(\vec n\) is the outward unit normal on the surface.
Newton’s law of cooling will become,
\[ - {K_0}\nabla u\,\centerdot \,\vec n = H\left( {u - {u_B}} \right)\]where \(H\) is a positive quantity that is experimentally determine, \({u_B}\) is the temperature of the fluid at the boundary and again it is assumed that this is only being evaluated at points along the boundary.
We don’t have periodic boundary conditions here as they will only arise from specific 1‑D geometries.
We should probably also acknowledge at this point that we’ll not actually be solving \(\eqref{eq:eq6}\) at any point, but we will be solving a special case of it in the Laplace’s Equation section.