Heat Conduction Cheat Sheet

I am dumping some equations here I need now and then! The sections about 3-dimensional temperature waves summarize what is described at length in the second part of this post.

Temperature waves are interesting for simulating yearly and daily oscillations in the temperature below the surface of the earth or near wall/floor of our ice/water tank. Stationary solutions are interesting to assess the heat transport between pipes and the medium they are immersed in, like the tubes making up the heat exchanger in the tank or the solar/air collector.

Heat equation – conservation of energy
Solutions for oscillating sources – temperature waves
1D solution
Helpers’ for the 3D case (spherical)
‘Helpers’ for the 2D case (cylindrical)
3D solution
Comparison of surface energy densities: 1D versus 3D
Stationary solutions
1D – plane
3D – sphere
2D – cylinder: pipe.
Comparison of overall heat flow: 1D versus 2D

Heat equation – conservation of energy [Top]

Energy is conserved locally. It cannot be destroyed or created, but it it is also not possible to remove energy in one spot and make it reappear in a different spot. The energy density η in a volume element can only change because energy flows out of this volume, at a flow density j (energy per area and time).

$\frac{\partial \eta}{\partial t} + \frac{\partial \vec{j}}{\partial\vec{r}} = 0$

In case of heat energy, the sensible heat energy ‘contained’ in a volume element is the volume times mass density ρ [kg/m³] times specific heat c [J/kgK] times the temperature difference in K (from a ‘zero point’). The flow of heat energy is proportional to the temperature gradient (with constant λ – heat conductivity [J/mK], and heat flows from hot to colder spots.

$\rho c \frac{\partial T}{\partial t} + \frac{\partial}{\partial\vec{r}} (- \lambda \frac{\partial T}{\partial\vec{r}}) = 0$

Re-arranging and assuming that the three properties ρ, c, and λ are constant in space and time, they can be combined into a single property called thermal diffusivity D

$D = \frac{\lambda}{\rho c}$

$\frac{\partial T}{\partial t} = D \frac{\partial}{\partial\vec{r}} \frac{\partial T}{\partial\vec{r}} = D \Delta T$

In one dimensions – e.g. heat conduction to/from an infinite plane – the equation is

$\frac{\partial T}{\partial t} = D \frac{d^{2} T}{d x^{2}}$

1D solution – temperature waves in one dimension [Top]

I covered it already here in detail. I’m using complex solutions as some manipulations are easier to do with the exponential functions than with trigonometric functions, keeping in mind we are finally interested in the real part.

Boundary condition – oscillating temperature at the surface; e.g. surface temperature of the earth in a year. Angular frequency ω is 2π over period T (e.g.: one year)

$T(t,0) = T_0 e^{i \omega t}$

Ansatz: Temperature wave, temperature oscillating with ω in time and with to-be-determined complex β in space.

$T(t,x) = T_0 e^{i \omega t + \beta x}$

Plugging into 1D heat equation, you get β as a function of ω and the properties of the material:

$i \omega = D \beta^2$

$\beta = \pm \sqrt{\frac{i \omega}{D}} = \pm \sqrt{i} \sqrt{\frac{\omega}{D}} = \pm (1 + i){\sqrt 2} \sqrt{\frac{\omega}{D}} = \pm (1 + i) \sqrt{\frac{\omega}{2D}}$

The temperature should better decay with increasing x – only the solution with a negative sense makes sense, then $T(\infty) = T_0$ . The temperature well below the surface, e.g. deep in the earth, is the same as the yearly average of the air temperature (neglecting the true geothermal energy and related energy flow and linear temperature gradient).

Solution – temperature as function of space and time:

$T(t,x) = T_0 e^{i \omega t - (1 + i) \sqrt{\frac{\omega}{2D}} x} = T_0 e^{i (\omega t - \sqrt{\frac{\omega}{2D}} x)} e^{-\sqrt{\frac{\omega}{2D}} x}$

Introducing parameter k:

$\sqrt{\frac{\omega}{2D}} = k$

Concise version of the solution function:

$T(t,x) = T_0 e^{i (\omega t - kx)} e^{-kx}$

Strip off the real part:

$Re(T(t,x)) = T_0 cos(\omega t - kx) e^{-kx}$

Relations connecting the important wave parameters:

$\tau = \frac {2 \pi}{\omega}$

$\lambda = \frac {2 \pi}{k}$

‘Helpers’ for the 3D case (spherical) [Top]

Basic stuff

$r = \sqrt{x^2 + y^2 + z^2}$

$\frac{\partial r}{\partial \vec{r}} = (\frac{\partial}{\partial x},\frac{\partial}{\partial y},\frac{\partial}{\partial z})\sqrt{x^2 + y^2 + z^2} = \frac{\vec{r}}{r}$

$\frac{\partial \vec{r}}{\partial \vec{r}} = (\frac{\partial}{\partial x},\frac{\partial}{\partial y},\frac{\partial z}{\partial z})(x,y,z) = 3$

$\Delta T = (\frac{\partial^2}{\partial x^2} + \frac{\partial^2}{\partial y^2} + \frac{\partial^2}{\partial z^2})T(r)$

Inserting, to obtain a nicely looking Laplacian in spherical symmetry

$\Delta T = \frac{\partial}{\partial\vec{r}} \frac{\partial}{\partial\vec{r}} T(\sqrt{x^2 + y^2 + z^2}) = \frac{\partial}{\partial\vec{r}} \frac{\partial r}{\partial\vec{r}} (\frac{dT}{dr}) = \frac{\partial}{\partial\vec{r}} (\frac{\vec{r}}{r} \frac{dT}{dr})$

$= \frac{3}{r} \frac{dT}{dr} - \frac{1}{r^2} \frac{\partial r}{\partial\vec{r}} \vec{r} \frac{dT}{dr} + \frac{\vec{r}}{r} \frac{\vec{r}}{r} \frac{d^2 T}{dr^2}$

$= \frac{3}{r} \frac{dT}{dr} - \frac{1}{r} \frac{dT}{dr}+ \frac{d^2 T}{dr^2} = \frac{2}{r} \frac{dT}{dr} + \frac{d^2 T}{dr^2}$

$= \frac{1}{r}(\frac{dT}{dr} + \frac{dT}{dr} + r \frac{d^2T}{dr^2}) = \frac{1}{r} \frac{d}{dr} (T + r \frac{dT}{dr}) = \frac{1}{r} \frac{d^2}{dr^2}(rT)$

‘Helpers’ for the 2D case (cylindrical) [Top]

Basic stuff

$r = \sqrt{x^2 + y^2}$

$\frac{\partial r}{\partial \vec{r}} = (\frac{\partial}{\partial x}, \frac{\partial}{\partial y})\sqrt{x^2 + y^2 } = \frac{\vec{r}}{r}$

$\frac{\partial \vec{r}}{\partial \vec{r}} = (\frac{\partial}{\partial x}, \frac{\partial}{\partial y})(x,y) = 2$

$\Delta T = (\frac{\partial^2}{\partial x^2} + \frac{\partial^2}{\partial y^2})T(r)$

Inserting, to obtain a nicely looking Laplacian in cylindrical symmetry

$\Delta T = \frac{\partial}{\partial\vec{r}} \frac{\partial}{\partial\vec{r}} T(\sqrt{x^2 + y^2}) = \frac{\partial}{\partial\vec{r}} \frac{\partial r}{\partial\vec{r}} (\frac{dT}{dr})$

$= \frac{\partial}{\partial\vec{r}} (\frac{\vec{r}}{r} \frac{dT}{dr}) = \frac{2}{r} \frac{dT}{dr} - \frac{1}{r^2} \frac{\partial r}{\partial\vec{r}} \vec{r} \frac{dT}{dr} + \frac{\vec{r}}{r} \frac{\vec{r}}{r} \frac{d^2 T}{dr^2}$

$= \frac{2}{r} \frac{dT}{dr} - \frac{1}{r} \frac{dT}{dr}+ \frac{d^2 T}{dr^2} = \frac{1}{r} \frac{dT}{dr} + \frac{d^2 T}{dr^2} = \frac{1}{r} \frac{d}{dr} (r \frac{dT}{dr})$

3D solution – temperature waves in three dimensions [Top]

Boundary condition – oscillating temperature at the surface of a sphere with radius R

$T(t,R) = T_R e^{i \omega t}$

Ansatz – a wave with amplitude decrease as 1/r. Why try 1/r? Because energy flow density is the gradient of temperature, and energy flow density would better decrease as 1/m² .

$T(t,r) = \frac{A}{r} e^{i \omega t + \beta r}$

Plugging in, getting β

$i\omega \frac{A}{r} e^{i \omega t + \beta r} = D \Delta T = \frac{D}{r} \frac{d^2}{dr^2}(rT)$

$= \frac{D}{r} \frac{d^2}{dr^2}(Ae^{i \omega t + \beta r}) = \frac{AD}{r} \beta^2 e^{i \omega t + \beta r}$

$i\omega = D \beta^2$

Same β as in 1D case, using the decaying solution

$T(t,r) = \frac{A}{r} e^{i \omega t + \beta r} = \frac{A}{r} e^{i (\omega t - kr)} e^{-kr}$

Inserting boundary condition

$T(t,R) = \frac{A}{R} e^{i \omega t + \beta R} = T_R e^{i \omega t}$

$\frac{A}{R} e^{\beta R} = T_R \Rightarrow A = T_R R e^{-\beta R}$

$T(t,r) = \frac{T_R R}{r} e^{-\beta R} e^{i\omega t + \beta r)} = \frac{T_R R}{r} e^{i\omega t + \beta(r-R)}$

$= \frac{T_R R}{r} e^{i(\omega t - k (r-R))}e^{-k(r-R))}$

The ‘amplitude’ A is complex as β is complex. Getting the real part – this is what you would compare with measurements:

$Re (T(t,r)) = \frac{T_R R}{r} cos(\omega t - k (r-R))e^{-k(r-R))}$

Comparison of surface energy densities: 1D versus 3D temperature waves [Top]

This is to estimate the magnitude of the error you introduce when solving an actually 3D problem in only one dimension; replacing the curved (spherical) surface by a plane.

One dimension – energy flow density is just a number:

$(t,x) = - \kappa \frac{dT}{dx} = - \kappa \beta T_0 e^{i \omega t + \beta x}$

Real part of this, at the surface (x=0)

$Re(j(t,0)) = - \kappa T_0 Re(\beta e^{i \omega t}) = - Re((-k -ik) \kappa T_0 e^{i \omega t})$

$= \kappa T_0 k (cos(\omega t) - sin(\omega t)) = \kappa T_0 k \sqrt{2} (cos(\omega t)\frac{1}{\sqrt{2}} - sin(\omega t))\frac{1}{\sqrt{2}})$

$= \kappa T_0 k \sqrt{2} (cos(\omega t)\cos(\frac{\pi}{4} - sin(\omega t))\sin(\frac{\pi}{4}) = \kappa T_0 k \sqrt{2} cos(\omega t + \frac{\pi}{4})$

How should this be compared to the 3D case? The time average (e.g. yearly) average is zero, to one could compare the average value for half period, when the cosine is positive or negative (‘summer’ or ‘winter’ average). But then, you can as well compare the amplitudes.

Introducing new parameters

$l = \frac{1}{k}$

$j_{amp} = \frac{\kappa T_0}{l}$

3D case: Energy flow density is a vector

$\vec{j}(t,\vec{r}) = -\kappa \frac{\partial T}{\partial \vec{r}} = -\kappa \frac{\partial}{\partial \vec{r}} \frac{T_R R}{r} e^{i\omega t + \beta(r-R)}$

$= -\kappa T_R R e^{i\omega t} [-\frac{1}{r^2} \frac{\vec{r}}{r} e^{\beta(r-R)} + \frac{1}{r} \beta \frac{\vec{r}}{r} e^{\beta(r-R)} ]$

$= \kappa T_R R e^{i\omega t} e^{\beta(r-R)} \frac{\vec{r}}{r} [\frac{1}{r^2} - \frac{\beta}{r} ]$

$= \frac{\vec{r}}{r} \kappa \frac{T_R R}{r} e^{-k(r-R)} e^{i(\omega t - k(r-R))} [\frac{1}{r} + k + ik]$

The vector points radially of course, its absolute value is

$j(t,r)= \kappa \frac{T_R R}{r} e^{-k(r-R)} e^{i(\omega t - k(r-R))} [\frac{1}{r} + k + ik]$

At the surface of the sphere the ‘ugly part’ is zero as

$\vec{r} = \vec{R}$

$r = R$

$k(r-R) = 0$

Real part:

$Re(j(t,r)) = \kappa T_R Re (e^{i(\omega t} [\frac{1}{R} + k + ik] )$

$= \kappa T_R [(\frac{1}{R} + k) cos(\omega t) - k sin(\omega t) ]$

$= \kappa T_R [k \sqrt{2} cos(\omega t + \frac{\pi}{4}) + \frac{1}{R} cos(\omega t)]$

Here, I was playing with somewhat realistic parameters for the properties of the conducting material. If the sphere has a radius of a few meters, you can ‘compensate for the curvature’ by tweaking parameters and obtain a 1D solution in the same order of magnitude.

Temporal change – there is a ‘base’ phase different between temperature and energy flow of (about) π/4 which is also changed by introducing curvature. I varied ρ,c, and λ with the goal to make the j curves overlap as much as possible. It is sufficient and most effective to change specific heat only. If the surface is curved, energy ‘spreads out more’. So to make it ‘as fast as’ the 3D wave you need to compensate by a giving it a higher D.

I did not bother to shift the temperature to, say, 10°C as a yearly average. But this is just a linear shift tat will not change anything else – 0°C is arbitrary.

1D stationary solution – plane [Top]

Stationary means, that nothing changes with time. The time derivative is zero, and so is the (spatial) curvature:

$\frac{\partial T}{\partial t} = 0 = D \frac{d^{2} T}{d x^{2}}$

The solution is a straight line, and you need to know the temperature at two different points. Indicating the surface x=0 again with 0 and the endpoint x_E with E, and using the definition of j in terms of temperature gradient and distance from the surface (x_E – 0 = Δx).

$|j(x = 0)| = \lambda |\frac{dT}{dx}| = \lambda \frac{|T_E - T_0|}{x_E} = \lambda \frac{|T_E - T_0|}{\Delta x}$ $

3D stationary solution- sphere [Top]

The time derivative is zero, so the Laplacian is zero:

$\frac{\partial T}{\partial t} = 0 = \Delta T(t, r) = \frac{1}{r} \frac{d^2}{dr^2}(rT)$

Ansatz, guessing something simple

$T(r) = \frac{A + Br}{r} = \frac{A}{r} + B$

Boundary conditions, as for the 1D case:

$T(R_0) = T_0$

$T(R_E) = T_E$

Plugging in – getting functions for all r:

$T(r) = \frac{1}{R_0 - R_E} [R_E T_E(\frac{R_0}{r} - 1) + R_0 T_0 (1 - \frac{R_E}{r}]$

$|j(r)| = \lambda \frac{1}{R_0 - R_E} \frac{1}{r^2} [R_E T_E R_0 - R_0 T_0 R_E ]$

At the surface:

$|j(R_0)| = \lambda \frac{1}{R_0 - R_E} \frac{R_E}{R_0} [T_E - T_0 ]$

2D stationary solution – cylinder, pipe [Top]

Cylindrical Laplacian is zero

$\frac{1}{r} \frac{d}{dr} (r \frac{dT}{dr}) = 0$

Same boundary conditions, plugging in

$r \frac{dT}{dr} = A$

$dT = A \frac {dr}{r}$

$\int_{T}^{T_0} dT = A \int_{R_0}^{r} \frac {dr}{r}$

$T(r) = T_0 + A \ln{(\frac{r}{R_0})} = T_0 + A (\ln{r} - \ln{R_0})$

$T(R_E) = T_E = T_0 + A \ln{(\frac{R_E}{R_0})}$

$A = \frac{T_E - T_0}{\ln{(\frac{R_E}{R_0})}}$

Solutions for temperature and energy flow at any r:

$T(r) = T_0 + (T_E - T_0) \frac{\ln{(\frac{r}{R_0})}}{\ln{(\frac{R_E}{R_0})}}$

$|\vec{j(r)}| = |\frac {1}{r} \lambda \frac{T_E - T_0}{\ln{(\frac{R_E}{R_0})}}|$

Expressing r in terms of distance from the surface, $\Delta r = r - R_0$

$|\vec{j(r)}| = |\frac {1}{\Delta r + R_0} \lambda \frac{T_E - T_0}{\ln{(\frac{R_1}{R_0})}}|$

Comparison of overall heat flow: 1D versus 2D [Top]

j is the energy flow per area, and the area traversed by the flow depends on geometry. in the 1D case the area is always the same area, equal to the area of the plane. For a cylinder, the area increases with r.

The integrated energy flow J for a plate with area F is

$J_{Plate} = F \lambda \frac{|T_E - T_0|}{\Delta x}$

If the two temperatures are given, J decreases linearly with increasing thickness of the cylindrical ‘shell’, e.g. a growing layer of ice.

For a cylinder of length l the energy flow J is…

$J_{Cyl} = 2 \pi l r |\frac {1}{r} \lambda \frac{T_E - T_0}{\ln{(\frac{R_E}{R_0})}}|$

$= 2 \pi l \lambda |\frac{T_E - T_0}{\ln{(\frac{R_E}{R_0})}}|$ \par

Factor r has been cancelled, and the for given temperatures J is only decreasing linearly with increasing outer radius $R_E$ . That’s why vendors of plate heat exchangers (in vessels with phase change material) worry more about a growing layer of sold material than user for e.g. ‘ice on coil’ I quoted a related research paper on ‘ice storage powered’ heat pump system in this post – they make exactly this point and provide some data. In addition to conduction also convection at both sides of the heat exchanger should be taken into account, too, in a ‘serial connection’ of heat transferring components.

3 Comments Add yours

Cleon Teunissen says:

September 30, 2018 at 10:40 PM

In general, these days images are likely to be displayed with distortion.

Computer screens are built with ever more pixels per centimeter.

The effect of that: the text of menus and dialog boxes becomes smaller and smaller.

Because of that more and more people are using the Windows option of displaying all texts of menus and dialog boxes at a 125% scale, or a 150% scale. (Also – and probably more importantly – on high pixel density devices, such as laptops, Windows default setting is 125% scaling.)

Our browsers query that display setting and adopt it. That is: when the Windows display is set to 125% then the browser will apply a 125% zoom to the web pages it displays, which includes scaling up every image to 125%.

As you can imagine, scaling up a PNG image of a mathematical formula to 125% gives noticable distortion.

Incidentally, Apple avoided this problem by making a factor 2 jump. The iPad went from 1024×768 to 2048×1536 pixels. If they would have done nothing then images in webpages would have been displayed at half the size of their display on a 1024×768 device. To keep the image size relative to the tablet the same they simply double up the pixels of the image.

Cleon Teunissen says:

September 30, 2018 at 7:00 PM

Hi Elke,

the font size of the equations has come out rather small.

Perhaps you are constrained by some maximum column with of the WordPress theme that you are using, but if not here is an idea that occurred to me some time ago.

I can generate the images displaying the mathematical formula’s in a very much oversized size, so the browser has to downsize them anyway. In the HTML/CSS I declare the target width and height of the image in ’em’ units in such a way that the font in the downsized image comes out about the same size as the ’em’ size of the page. As you know, browsers support that the user sets his/her own preferred font size, overruling any font size in the page’s HTML/CSS. With the image dimensions declared in units of ’em’: whatever the page visitor is doing the size of the formula’s will be in proportion.

I tried that on one of the pages of my website:
– Firefox: looks great
– Google Chrome: looks great
– Internet Explorer: looks bad; bad downsizing
– Microsoft Edge: looks pretty good

Cleon Teunissen

1. elkement says:
  
  September 30, 2018 at 8:39 PM
  
  There is a size parameter I can set when using the LaTex ‘add-in’ provided by wordpress.com – it’s like a query string appended to the actual latex. I experimented with it on other posts, but I only liked the larger size for some of the formulas.
  
  So I’d have to test with different sizes for each formula – some additional work I did not want to do for this post, admittedly. And as you say – one would have to test with different browsers anyway.
  But on principle I can change the size.

elkemental Force