Thermal Stacks, Transmission Lines, high power LED’s and Temperature Sensitive Paints (TSP)

I’m writing this because even though it’s already 2020, this kind of stuff still cannot be found anywhere in the internet!

Stacking layers of materials is something we all do in so many engineering applications! Electronic components, batteries, constructions panels, ovens, refrigerators and so many other layered materials that are inevitably under heat transfer. We want to predict the temperatures on the stack, at all layers, to ensure the materials are going to withstand the transfer of heat. In order to do that, any engineer worth their salt will already have used the Thermal/Electrical Circuit Analogy. And it works wonders!

But what if we want to do transient analyses? What if the materials are going to be stressed by some quick burst of power; or some sequence of bursts at a certain frequency? Of course, we can always just do some FEM and solve our issue – but one might find that this can become a multi-scale problem if there’s low frequencies AND high frequencies involved.

This is commonly the case in semiconductors. Applications I’m familiar with are related to LED pulsing, but I’m sure one will find many others. Another more obscure application I found in my research has to do with unsteady measurements of surface temperature when fluid is flowing on a surface. Some frequencies will be too fast and the stack of materials will not respond, so how do we proceed in that scenario? How to predict the temperature rises and – more importantly – understand what we need to change in the design in order to improve its frequency response?

Well, it turns out we can extend the thermal circuit analogy to transients by including the component’s thermal capacitance. Most layers in a stack would then be easily lumped into RC low-pass filters. This is the “basic approach” and it works – to a point.

So let’s examine it, to understand its pros and also its limitations:

The thermal capacitor

Capacitance in electrical circuits is the ability of a given electrical component to store charge. In the thermal analogy, this “charge” is heat – and the thermal capacitance is just the mass of the component times its heat capacity (C=m \times c_p). If we look at the energy conservation equation for a lumped element:

\displaystyle q = m c_p \frac{\partial T}{\partial t}

Which means that, for a given heat flow q, a small capacitance C leads to a large rate of temperature rise (large dT/dt). I will make this problem more 1-dimensional to consider the case where the stack of layers is much thinner than wide, because that’ll lead to interesting observations when we drop the “lumped element” part of this model. For now, let’s simply divide both left-hand and right-hand sides by the sheet area A, to get everything in terms of heat fluxes:

\displaystyle \frac{q}{A} = q" = \frac{m}{A} c_p \frac{\partial T}{\partial t}

\displaystyle q" = \frac{\rho A L}{A} c_p \frac{\partial T}{\partial t}

\displaystyle q" = \rho c_p L \frac{\partial T}{\partial t}

The terms \rho c_p L then are what we will call the specific component’s capacitance (C^*), in the SI units of J/m2.K. One can then represent a layer of a stack like the following RC circuit:

Obviously, this single-layer model is rather simple and does not constitute of much of an issue. The RC circuit has a time constant, \tau=R^*C^*, where R^*=L/k. You can check the units: R^* has units of m2.K.s/J, which times the units of C^* gives us a time constant in seconds. The cut-off frequency is f_{cutoff}=1/\tau, and we get a -20dB/decade low-pass filtering behavior.

It is really worth observing, though, that the capacitor is not connected to T_2, but to an elusive “ground” temperature. This is because, as long as you have a non-zero \partial T/\partial t, you will have heat transfer to the thermal capacitance. Therefore, the ground level can be considered absolute zero temperature – for didactic purposes. The temperature of this ground level is not important, though: since capacitors always block DC, its absolute value will never play a role.

Stacking them up (plus all the ways you will mess this up)

For thermal stacks, it is not straightforward how to assemble the circuits together. Do we put the RC circuits in series with each other? As it turns out – not really. To understand why, let’s think about what happens in a two-layer stack:

Looking at the picture above, we have a stack of two layers and 3 different circuits that could (potentially) represent it. The circuit (a) represents how we would construct the thermal RC network using the same strategy as for the single layer explained before. (a.1) is simply another representation of the same circuit, perhaps less “bulky”. The circuits (b) and (c) are wrong because they do not represent the thermal equations that describe the stack. Let’s see why:

The first layer of the stack has an input heat flux of q_1", which according to the energy conservation equation applied to this lumped element, can be split into the heat flux that stays in the element (q_{c1}") and the heat flux that goes through the element (q_i"). This heat flux then goes to the second layer and splits again into q_{c2}" and q_2". The energy conservation equations here are equivalent to the Kirchoff current laws applied to the nodes at T_1 and T_i in the circuit representation (a):

\displaystyle q_1"=q_i" + q_{c1}"

\displaystyle q_i"=q_2" + q_{c2}"

This sounds easy because you already have the answer spoon-fed to you. But if you were left to your devices, I would bet many of you (including past me) would come up with the circuit representations (b) and (c). (b) doesn’t make any sense because capacitors block DC, therefore there would be no DC heat flux through the stack – total non-sense!

The circuit in (c) also does not make sense, but it is not clear why at a first glance. If you follow Kirchoff’s laws, however, you will see that (c) cannot be right, because in the node at T_i the heat flux equation would have to be q_1" + q_{c1}=q_i" + q_{c2}" – which is simply not the correct energy balance.

It is not really very intuitive that the capacitor of the first layer of the stack should be connected to ground, instead of T_i. In fact, I made quite a few times this mistake. But looking at the energy balance, we can’t really deny there is a correct equivalent circuit.

We can also solve the full heat pde with circuits!

Ok, so we tackled the single wall and the stacked double-wall problem. I think having gone through how this can go wrong, you can now build this to multi-walled systems easily. As you probably already realized, each new layer will behave as a new first-order system, and some very interesting behavior can appear as extra poles are introduced by the new layers. One thing, however, that is worth noting, is that this approach completely fails at high enough frequencies.

To understand why, let’s think about it a little further. If the temperature oscillations in T_1 are high frequency (i.e., their period is below the R^* C^* time constant), the lumped element model would predict a decay rate of -20dB/decade, which is what is expected from a standard first-order system.

But the problem in the (zero-dimensional) lumped element model is that the whole element’s temperature is assumed to be a constant value. At high enough frequencies, this won’t be the case anymore. We can see why the lumped element model becomes inaccurate by modeling the standard 1D heat equation over a slab:

\displaystyle -\frac{\partial T}{\partial x} = R' q"

\displaystyle -\frac{\partial q"}{\partial x} = C' \frac{\partial T}{\partial t}

The negative sign in the left-hand side is related to the direction of the coordinate system. One can recover the same equations by distributing the lumped element model related to each layer into the following circuit:

Where R'=1/k (thermal resistance per unit length) and C'=\rho c_p (thermal capacitance per unit length). This circuit immediately follows from the circuit presented in (a.1) in the previous figure; since one can consider a single layer can be decomposed into an infinite stack of infinitely thin layers. Not all mechanical engineers will recognize, however, that this corresponds to a (simple) version of the 1-dimensional transmission line model. Robertson and Gross (1958, “An Electrical-Analog Method for Transient Heat-Flow Analysis”, Journal of Research of the National Bureau of Standards) did recognize this in their paper, for example. Funnily enough, their figure shows they also got the RC network wrong, since the first resistor is shown between the first capacitor and T1 – an indicator of how tricky this business can be!

A thermal layer is actually a transmission line

I invite you to read the article on Wikipedia about transmission lines if you’re not familiar with this concept, or read the book “Time-Harmonic Electromagnetic Fields”, by Harrington, 1961. Transmission line is an electrical model for pairs of conductors at high frequencies. At high enough frequencies, the pair of conductor’s distributed capacitance and inductance creates waves inside the conductor, which can be modeled by the wave equation. The element of the electrical transmission line is displayed below (figure from the Wikipedia article):

The electrical transmission line follows the 1-dimensional damped wave equation, which is derived from the Telegrapher’s equations:

\displaystyle -\frac{\partial V}{\partial x} = L' \frac{\partial I}{\partial t} + R' I

\displaystyle -\frac{\partial I}{\partial x} = C' \frac{\partial V}{\partial t} + G' V

We can easily see this is the wave equation by taking the Fourier Transform of the two equations:

\displaystyle -\frac{\partial V}{\partial x} = (j\omega L'+ R') I

\displaystyle -\frac{\partial I}{\partial x} = (j\omega C'+ G') V

And then replacing, say, the x derivative of the first equation in the second:

\displaystyle -\frac{\partial^2 V}{\partial x^2} = (j\omega L'+ R') \frac{\partial I}{\partial x}

\displaystyle \frac{\partial I}{\partial x} = -\frac{1}{j\omega L'+ R'} \frac{\partial^2 V}{\partial x^2} = (j\omega C'+ G') V

\displaystyle \frac{\partial^2 V}{\partial x^2} - (j\omega C'+ G')(j\omega L'+ R') V = 0

\displaystyle \frac{\partial^2 V}{\partial x^2} + \omega^2 C' L' V - j\omega (L'G'+ C'R') V - R' G' V = 0

Taking the IFT:

\displaystyle \frac{\partial^2 V}{\partial x^2} - C' L' \frac{\partial^2 V}{\partial t^2} - (L'G'+ C'R') \frac{\partial V}{\partial t} - R' G' V = 0

Which, obviously, is the wave equation. It’s damped by the distributed resistance R' and conductance G'. If we come back to the thermal stack and consider there’s no conductance (G'=0) and no inductance (L'=0), we recover the heat equation:

\displaystyle \frac{\partial^2 V}{\partial x^2} - C'R' \frac{\partial V}{\partial t} = 0

Where the thermal diffusivity is \alpha=1/C'R' = k/\rho c_p ; as expected.

how this interpretation can be useful

Ok, by now you might be finding this all quite trivial and, perhaps, an over-complication of the heat transfer problem. Why would one need to invoke the far more complex “transmission line” to describe the behavior of such a simple system? Well, it basically allows us to easily determine the temperature transfer functions at different points in the stack (perhaps, easily is not the best word, but with a couple easily-programmable equations), without having to solve the differential equations numerically. Plus, it allows us a very cool insight about the behavior of the themal stack at higher frequencies. As it turns out, at high enough frequencies the first layer’s characteristic impedance – a material property – dominates the transfer function.

Thermal transmission lines in practice: the pulsed led

Let’s first analyze the following case that I find is the most practical for high frequency thermal stacks: An electronic power dissipation device that sources power in short bursts. I mentioned in a previous post that you can massively overdrive LEDs to produce a lot more light than “designed” for short bursts, which enables one to expose high speed cameras – like a camera flash. We’ll look at the Luminus Devices Phlatlight CBT-140 for this analysis. The reason why AC thermal analysis is important in pulsed LEDs is because, even though the mean temperatures of each component might be within spec due to the high thermal capacitance of the heat sink, the power burst will temporarily increase the temperature of the LED die, and we obviously also don’t want that to exceed the specs of the LED. Concerns about thermal expansion will not be addressed in this discussion – but you should also do your research on that!

The CBT-140 die looks like shown below. The named squares in the electrical analogy (bottom right) represent the “transmission lines” for each component.

We will consider the die is made of sapphire (many high power LED dies are made of it). This might be incorrect for this model – but I could not find any specs on that. The die has 0.31mm thickness; whereas the copper has 3.44mm thickness. The aluminium heat sink is pretty large. For the sake of this analysis we’ll consider 10mm. As we’ll detail later – it won’t matter much for the AC analysis.

So what we want to find is the “input impedance” as seen by the current source just above the die in the figure. The current source is a heat flux source, which in this case would be a square wave with some (low) duty cycle, since power is generated by the (micrometer thin) LED just above the sapphire die. Because there’s a cascade of 3 “coupled transmission lines”, we have to use the transmission line input impedance formula three times. From Wikipedia:

\displaystyle Z_a = Z_{0,a}\frac{1+\Gamma_L e^{-2\gamma l_a}}{1-\Gamma_L e^{-2\gamma l_a}}

Where Z_a is the input impedance as seen by the node at T_a and \Gamma_L=(Z_L-Z_0)/(Z_L+Z_0) is the reflection coefficient due to the impedance mismatch of the aluminium slab and the convective boundary (also the load impedance) R_{conv}=Z_L. Z_0 in this case is the aluminium’s characteristic impedance, which in general has the form:

\displaystyle Z_0 = \sqrt{\frac{R'}{j\omega C'}} = \sqrt{\frac{1}{j\omega (k \rho c_p)}}

Note the characteristic impedance is inversely proportional to the square root of the frequency. This will explain the “half-order” decay (-10dB/decade) behavior in the frequency domain we’ll see later (in contrast with a first-order expected from these RC networks). The square root of 1/j predicts a phase lag of 45º at high frequencies (not 90º!). Also referring to the Wikipedia article, the propagation constant \gamma is a material property:

\displaystyle \gamma = \sqrt{j\omega R' C'} = \sqrt{j\omega \alpha^{-1}}

As we discussed, the math is not the “prettiest”, but this can all easily be hidden behind an Excel or Matlab formula. In this case, using the properties of the three materials shown in the table below and a convection coefficient of h=100 W/m^2 K \rightarrow R_{conv}=0.01 m^2 K/W, we get:

Aluminium 6061CopperSapphire
k180 W/m K287 W/m K34.6 W/m K
c_p1256 J/Kg K376 J/Kg K648 J/Kg K
\rho2710 Kg/m38800 Kg/m33930 Kg/m3
Z_04.04\times 10^{-5} \sqrt{-j\omega^{-1}}3.24\times 10^{-5} \sqrt{-j\omega^{-1}}1.06\times 10^{-4} \sqrt{-j\omega^{-1}}
\gamma137.5 \sqrt{j\omega}107.4 \sqrt{j\omega}271.3 \sqrt{j\omega}

In order to find the “input impedance” at the node T_c, we use the Z_a as the input load (Z_L), and then iterate. Since we’re looking for the maximum temperature in the stack, T_d, let’s calculate the overall input impedance of the stack Z_d:

In this case we are driving the LED with a 500 nanosecond square pulse that repeats at 100 kHz. The drive current of the LED is 200A, and given the die area of 14 mm2 and a voltage drop of 20V in the LED, we have an instantaneous power of 285 MW/m2. That’s a lot of power!! Well, as it turns out, the mean flux of power (in this case, at a duty cycle of 5%, 14.2 MW/m2) appears in the spectrum at the DC component, but the sharp bursts of power also generate incredibly high amounts of high frequency power flux that really need to be considered.

Looking at the figure above, the solid, thick black line is the actual frequency transfer function of the whole thermal stack. The dashed blue line represents the lumped element assumption and it is a good model until about a frequency of 0.05 Hz. It then starts to fail when the aluminium starts to kill the waves very close to its surface due to its “large characteristc impedance”, (Z_0). You can think about it like; its large capacitance combined with its (relatively) low conductivity prevents the thermal oscillations on one side of the material to immediately propagate to the opposite side – invalidating the lumped element assumption.

This issue happens once again for the copper (not clear exactly where) and then finally for the die itself (the sapphire) above the frequency of 40 Hz. Then the characteristic impedance of the sapphire dominates the transfer function – i.e., the sapphire is thick enough to be effectively “infinite”, as far as high frequency heat transfer is concerned (the sapphire is 0.31mm thick!!). The temperature waves then decay inside the sapphire, and very little oscillations can be detected at its opposite side. At those frequencies, it doesn’t matter anymore which material you put at the interface. It will not see any temperature oscillations.

All of this insight we get without needing to numerically solve the heat equation is quite fascinating. We now know that at high frequencies, the dominant part of the frequency response is a material property – meaning that a material choice here is crucial for high frequency heat transfer!

How important? Well, let’s run the numbers and see what the wave looks like. I removed the DC part of the wave, because this LED really needs some thermal mass attached to it (the aluminium block) as no amount of convection or heat-sinking will get rid of this insane amount of power. Therefore, let’s only focus on the oscillatory part of the wave:

We can see that we predict an overtemperature of almost 25 K for the die surface for this case. This is on top of the die heating up from the constant heat flux at DC. It is quite a lot and could (in fact, DOES) easily lead to die failure if one only reads the copper temperature. Therefore, one needs to account for this overtemperature when designing an emergency shutoff system – otherwise, it’ll only think about shutting off 25 degrees after the LED failed!

Temperature sensitive paint – high frequency measurements

So this one is quite an interesting application from the research standpoint! We use temperature sensitive paints to perform measurements in aerodynamic models. We paint an aerodynamic surface with this special paint, and it fluoresces under UV illumination, emitting a strong (usually red) light signal at colder temperatures and a weaker signal as the paint surface becomes hotter. An example of supplier of these paints is given here. Though this is quite a specialized application used mostly in high end aerodynamics research, there’s still a lot of interest into understanding its behavior.

So let’s say we want to perform unsteady measurements of heat transfer coefficient at the surface. We put the surface under a flow of nominally known temperature (for this example, a jet impinging on a surface) and then we use the TSP technique to image the temperature map over the impinging surface; getting a full temperature map. As the flow oscillates, we have fluctuations of this convection coefficient – which translate to measured fluctuations in temperature. Obviously, we might not be able to see them if the temperature does not fluctuate much because the substrate of the surface has too much thermal capacitance (as discussed throughout this article). In fact, now we know we have waves of temperature going through the material, which will be rapidly decaying. It might perhaps be the case that the waves are so weak that we cannot detect them. Therefore, we need to assess whether we can measure these waves by comparing the expected signal generated by the oscillations with the dynamic range of our camera.

Let’s start with a rather simplistic setup – to illustrate the issue. Say we have an impinging jet as in the figure below. The impinging jet is the flowing arrows at T_{\infty}.

A “thin” layer of temperature sensitive paint is painted on top of a 10mm thick aluminium plate. The layer of paint, although thin, will have an effect on the transfer function. Here’s what happens as we change the paint thickness. The solid black line represents a reasonable paint thickness (5 microns).

As we can see, the transfer function is highly dependent on the thickness of the paint: If the paint is too thin, we get dominance of the aluminium block on the transfer function, which reduces the signal gain very significantly at higher frequencies. As the paint layer thickens, it starts dominating the transfer function – but knowing the paint thickness might be a challenge.

A solution to this issue is to use an insulating film of similar properties to the paint binder (in this case, the binder used for the calculations was a Polyacrylic Acid polymer). Most importantly, if the two materials have similar characteristic impedances (Z_0), then the transfer function obtained will not be as dependent on the paint thickness.

Polyacrylic acid has the following thermal properties (sources here and here). The film considered here is a PET film (this one).

Polyacrylic acidPET Film
k0.543 W/m K0.21 W/m K
c_p1248 J/Kg K1465 J/Kg K
\rho1440 Kg/m31170 Kg/m3
Z_01.01\times 10^{-3} \sqrt{-j\omega^{-1}}1.67\times 10^{-3} \sqrt{-j\omega^{-1}}
\gamma1819 \sqrt{j\omega}2856 \sqrt{j\omega}

As we can see, the thermal impedance of the two materials is somewhat similar. This means less reflection will occur at the interface; which will help with a reduced dependence of the transfer function on the paint coat thickness. The circuit analyzed becomes the following:

This graph shows a great improvement; because the transfer function at higher frequencies is less dependent in the paint thickness. Not perfect – this would only be possible with a perfect match of thermal impedances. Nevertheless, this can definitely be used to make some measurements.

As a remark – the frequencies of interest in this work are around 20 kHz. The gain is then about -38dB (5 microns paint thickness); which is defined as 20 \log_{10}[(T_p-T_h)/(T_\infty -T_h)]. This corresponds to a signal gain of 0.011 – or; for every 1 degree of temperature fluctuation in T_\infty, we get one hundreth of a degree of fluctuation in the paint surface temperature; which is quite a small signal. Since the TSP gives a signal of about 1.4%/ºC, we are talking about an intensity signal of I/I_0 = 0.014 \times 0.011 = 1.54 \times 10^ {-4}, or 0.6 counts/ºC of a 12-bit camera where I_0 exactly saturates the camera (most likely not). I expect this paint measurement system to (hardly) detect temperature fluctuations of 10ºC pk-pk; considering good post-processing and noise reduction to reject camera noise and a very tame signal (like a pure frequency).

I’ll post updates as we perform experiments!

Tying it all back to the fourier number

So we see this interesting transfer function that switches between the lumped RC circuit assumption (with a -20dB/decade decay) to the continuous transmission line (with a -10dB/decade decay). Why does that happen and where is the “kink” in the graph?

Well, we’ll see this ties back to the thermal theory; once we find the kink frequency: The kink is where the -20dB/decade and the -10dB/decade cross. At the -20dB/decade regime, a slab of material has behaves as a capacitor and its thermal impedance is:

\displaystyle Z_C=\frac{1}{j \omega C} = \frac{1}{j \omega C' L}

We have the definition of C' in the previous sections – so let’s focus on the math. The thermal impedance of the slab in the -10dB/decade regime is simply its characteristic impedance (another interpretation is that in the -10dB/decade regime the slab is effectively infinite).

\displaystyle Z_0=\sqrt{\frac{R'}{j \omega C'} }

The kink is located about where |Z_C|=|Z_0|:

\displaystyle \bigg|\sqrt{\frac{R'}{j \omega C'} }\bigg| = \bigg|\frac{1}{j \omega C' L}\bigg|

\displaystyle \frac{R'}{\omega C'} = \frac{1}{\omega^2 C'^2 L^2}

This gives us the radian frequency of the kink:

\displaystyle \omega_{kink}=\frac{1}{L^2} \frac{1}{R' C'}

Now you should already be recognizing that we have an interesting non-dimensional coming up, as \alpha = 1/R' C':

\displaystyle \omega_{kink}=\frac{\alpha}{L^2} \rightarrow \boxed{f_{kink}=\frac{\alpha}{2 \pi L^2}}

The non-dimensional is; obviously; the Fourier number:

\displaystyle Fo=\frac{\alpha}{\omega_{kink} L^2}=1

We usually see the Fourier number with a time in the numerator (Fo=\alpha t/L^2), but this is the same kind of transient phenomenon. If the Fourier number is larger than 1, the “time scales” are large (low frequencies) and the heat transfer of the slab can be modeled by a simple capacitor. If the Fourier number is smaller than 1, then the “time scales” are small (high frequencies) and the heat transfer needs the heat equation to be modeled, which we are simply doing by using the transmission line model. The analogy is quite powerful and enables simple, quick analysis of these systems.

I boxed the f_{kink} equation because it is quite useful to determine whether our AC analysis needs to consider the transmission line model; knowing the frequencies of interest. Hope this provides the insight you needed for your analyses!

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Google photo

You are commenting using your Google account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s