One dimensional steady state diffusion, with and without

source. Effective transfer coefficients

21 mars 2017

For steady state situations (∂t = 0) and if convection is not present or negligible the transport
equation reduces to Laplace’s equation ∆H = 0 or Poisson’s equation ∆H = RH if there is a
source term. We consider here a few cases of one dimensional steady state diffusion and we also
introduce the notions of heat transfer coefficient and thermal resistance.

1 Boundary conditions and transfer coefficients

1.1 Conditions on temperature, concentration and their gradients
At the interface between two materials, we need to specify the boundary conditions which will
apply to the temperature and concentration fields. The temperature field has to be continuous at
an interface, otherwise there would be a diverging temperature gradient, and according to Fourier’s
law, an infinite heat flux. If we use a continuum model of matter, an interface is a two dimensional
object with zero volume and there cannot be a production or storage of heat at the interface. As
a consequence, the heat flux normal to the interface should also be continuous.
Similar conditions apply for mass transport : mass cannot accumulate at an interface of ne-
gligible volume and the mass flux should be continuous across an interface. However, there are
case that are quit specific to mass transport. Semi permeable membranes can restrict the motion
of some chemical species and impose a discontinuity of concentration at an interface. Chemical
reactions can take place at an interface, for example at the surface of a solid material acting as a
catalyst. When a reaction takes place, the mass flux JC will depend on the kinetics : JC = kC m
where m is the order of the reaction and k is a rate constant, usually varying with temperature
k ∝ exp(−E/kB T ), E being an activation energy. We will discuss in more detail these boundary
conditions with chemical reactions when we analyze the design of a biochemical sensor.

1.2 Effective transfer coefficients at an interface

Aside from the continuity of temperature or concentration and their gradients, we can in
many cases specify the value of temperature or concentration (Dirichlet boundary condition) or
specify the value of the gradient of temperature or concentration (Neumann boundary condition).
If transport occurs by diffusion only (i.e. when P e  1), the flux at an interface is simply given
by the gradient of temperature or concentration multiplied by either the thermal conductivity or
the mass diffusivity. However if the transport occurs by a combination of diffusion and convection
(P e ≥ 1), the flux will depend in general of the velocity field and we will need a detailed analysis
of the flow to compute it. Nevertheless, the convection diffusion equation is linear in temperature
and concentration and this linearity implies that fluxes will vary linearly with the temperature
and concentration differences imposed at the boundary of the system. This is why we can write
in general the heat flux JH or mass flux JC at a boundary in the following manner :

JH = hH (T1 − T0 ) or JC = hC (C1 − C0 ) (1)

1
where T1 (resp. C1 ) and T0 (resp. C0 ) are the temperatures (resp. concentrations) on either side
of the boundary and hH (resp. hC ) is an effective heat transfer (resp. mass transfer) coefficient.
The effective transfer coefficients are determined experimentally, but in some cases that we will
discuss later in the course, they can be determined analytically. In the two problems considered
below (heat flux through a composite wall and bioheat equation), we use these transfer coefficients
to specify boundary conditions.

2 Planar geometry
If the geometry of the problem is such that gradients of temperature or concentration occur
only in one directions (say x) in space, the transport equation in steady state is :
∂2T
=0 (2)
∂x2
or
∂2T q̇(x)
+ =0 (3)
∂x2 λ
with a source term q̇(x) giving the amount heat produced par unit volume and per unit time. Here
we consider specifically an heat transfer problem, since there are many examples in applications,
but a steady state 1D mass transfer problem would be formally identical.

2.1 Thermal resistance

Without a source term, the solution of the equation has the general expression T (x) = A + Bx,
A and B depending on the boundary conditions. The flux is given by : J = −λdT /dx = −λA
where λ is the thermal conductivity. If we consider now a material of thickness e in the x direction,
the temperature gradient is (T1 − T0 )/e = ∆T /e where T1 and T0 are the temperatures on the
left and right sides of the material. Hence, the heat flux is J = −λ∆T /e. Using an analogy with
Ohm’s law for electrical conduction I = U/R, the heat flux being analogous to the current and
the temperature difference to the electrical potential, we define a thermal resistance :
e
RT = . (4)
λ

2.2 Heat transfer through a composite planar structure

The notion of thermal resistance can be applied readily to composite planar systems such as
a wall made of three different layers (fig. 1). There is a plaster layer of thickness ep =1cm, with
thermal conductivity λp = 0.2 W/m.K, a layer of glass fiber (ef =10cm, λf = 0.04 W/m.K) and an
outside layer of wood ( (ew =2cm, λf = 0.12 W/m.K). The air inside the house is at temperature
Ti =20◦ C and the air outside at To =-5◦ C. We want to compute the heat flux J through the wall
to estimate the power needed to keep the inside air at 20◦ C.
The heat flux is constant through the wall and is given by :

(T1 − T2 )λp (T2 − T3 )λf (T3 − T4 )λw

J= = =
ep ef ew
T1 , T2 , T3 and T4 being the temperatures at the different interfaces from left to right.
We have for the temperatures :
Jep Jef Jew
T1 = T2 + and T2 = T3 + and T3 = T4 +
λp λf λw
and the total temperature difference is :
 
ep ef ew
T1 − T4 = J + +
λp λf λw

2
 plaster glass fiber wood
T1 T4

air inside air outside

Ti To

ep ef ew

Figure 1 – A composite wall

The system is analogous to three resistors in series, we sum the individual thermal resistances
to get the global resistance. If we examine the different resistances, we see that λf is an order of
magnitude smaller than λp and λw . In addition, ef is an order of magnitude larger than ep and
ew . So we can consider that the thermal resistance of the wall is dominated by the layer of glass
fiber and make the assumption :
ef
T1 − T4 ≈ J .
λf
To close the problem, we need to write the boundary conditions inside and outside the house.
The heat transfer at the inside and outside interfaces is dominated by convection and we use heat
transfer coefficients determined experimentally, namely hi =30 W/m2 .K and ho =60 W/m2 .K. So,
using these heat transfer coefficients, the heat flux is given by :

J = (Ti − T1 )hi = (T4 − To )ho .

Then :
J J
T1 = Ti − and T4 = To + .
hi ho
Finally :  
1 1 ef
T1 − T4 = Ti − To − J + =J .
hi ho λf
and the value of the heat flux is given by :
!
1
J = (Ti − To ) 1 1 ef
hi + ho + λf

Again, if we examine the values of the different terms, we can see that the dominant one is the
thermal resistance of the glass fiber and we have : J ≈ (Ti − To )λf /ef = 10 W/m2 .

3 Cylindrical geometry
In a cylindrical geometry, where the temperature depends only on the radius r, the transport
equation is in steady state :
 
1 d dT
r =0 (5)
r dr dr
Integrating twice, we the general solution : T (r) = A ln r + B where A and B are determined
by the boundary conditions.

3
3.1 Uniform heat source
With a source term, the equation is :
 
1 d dT q̇(r)
r + =0 (6)
r dr dr λ
If the production of heat is uniform in space, the integration with respect to r gives the general
solution :

q̇r2
T (r) = − + A ln r + B (7)
4λ
If we consider the case of a cylindrical metal wire with radius R in which an electric current
generates heat by Joule effect at a rate q̇, we can determine the temperature distribution within
the wire from eqn. 7 by noting that T cannot diverge at the center of wire (r = 0), hence A = 0,
and by specifying the surface temperature at R.
If we prescribe the surface temperature Ts , we have : Ts = −q̇R2 /4λ + B and :

q̇(R2 − r2 )
T (r) − Ts = . (8)
4λ
The temperature distribution within the wire is parabolic with a maximum on the axis Tm =
Ts + q̇R2 /4λ.

3.2 Distributed heat source. The bioheat equation.

Mammals need to regulate precisely their internal temperature. This regulation is achieved by
several mechanisms, one of which being the transport of heat by blood from the core of the body
to the tissues in the limbs.
H. Pennes performed in 1948 a pioneering work in this domain by performing temperature
measurements within forearms of human subjects (fig. 2) and developing a simple heat transfer
model : the bioheat equation. In this model, the forearm is approximated as a cylinder made of a
uniform muscle tissue and the temperature is assumed to depend only on the radial coordinate
r. Heat is generated within the muscles by metabolic activity at a constant rate qm ˙ . In addition
there is a heat transfer from the blood to the surrounding tissue at a rate q˙b proportional to the
difference of temperature between the arterial flow Ta and the local temperature T (r), to the mass
flow rate of blood per unit volume of tissue ρb Vb and to the specific heat per unit mass of blood
Cb , i.e. q˙b = (Ta − T (r))ρb Cb Vb . Putting these two heat sources in the transport equation 6 gives
the bioheat equation in cylindrical coordinates :
 
1 d dT ˙ + (Ta − T (r))ρb Cb Vb
qm
r + =0 (9)
r dr dr λ
The bioheat equation can be rewritten as :
   
1 d dT ρb Cb Vb qm
˙
r − T (r) − Ta − =0 (10)
r dr dr λ ρb Cb Vb
or :
d2 θ dθ
r2 2 + r − Kr2 θ = 0 (11)
dr dr
where θ = T (r) − Ta − qm ˙ /(ρb Cb Vb ) and K = ρb C√b Vb /λ. This last equation is a modified Bessel
equation, with a solution of the form θ = AI0 (r K) where I0 is the modified Bessel function
of the first kind with zeroth √ order. The temperature at the surface of the limb (r = R) is Ts =
Ta + qm
˙ /(ρb Cb Vb ) + AI0 (R K). We assume now that the heat flux on the skin is given by a heat
transfer coefficient hs , we have the following boundary condition :
 
dT
−λ = hs (Ts − T0 ) (12)
dr r=R

4
 37 -1
Vb=0.0005 s
-1
Vb=0.0003 s
-1
Vb=0.0001 s

36

35

)C¡( erutarepmeT
34

33

32

31

30
0 1 2 3 4

Radial position (cm)

Figure 2 – Temperature profiles though a human forearm. Left, experimental results from H.
Pennes, Journal of Applied Physiology, 1, 93 (1948). Different curves correspond to different room
temperatures. Right, solutions of the bioheat equation for a cylinder of radius 5 cm, an outside
temperature of 25◦ C and three different values of the blood flow rate.

where T0 is the temperature of the ambient air. Using the relation for the Bessel functions :
dI0 (r)/dr = I1 (r), the boundary condition gives :
√ √ √
− λA KI1 (R K) = hs [Ta − T0 + qm
˙ /(ρb Cb Vb ) + AI0 (R K)] (13)
yielding the integration constant :

hs [Ta − T0 + qm
˙ /(ρb Cb Vb )]
A=− √ √ √ (14)
λ KI1 (R K) + hs I0 (R K)
The temperature distribution is finally given by :
√
qm
˙ hs I0 (r K)[Ta − T0 + qm
˙ /(ρb Cb Vb )]
T (r) = Ta + − √ √ √ (15)
ρb Cb Vb λ KI1 (R K) + hs I0 (R K)
or
√
hs I0 (r K)θ0
θ(r) = √ √ √ (16)
λ KI1 (R K) + hs I0 (R K)
and the skin temperature is given by :
θ0
θs = θ(R) = √ √ (17)
λ KI1 (R K)
1+ √
hs I0 (R K)
√ √
There are two dimensionless parameters involved√in this result : R K and λ K/hs . From
a dimensional analysis of eqn. 9, we can see that 1/√ K is a length scale such that the diffusive
flux balances the flux due to the√blood flow. If R K is small, this means that the radius R of
the limb is small compared to 1/ K and we expect that losses by diffusion will be large and the
heat convected by the blood flow will not be sufficient to maintain the skin temperature at a high
enough value (say to prevent the skin from freezing).

5
 √
The other dimensionless parameter λ K/hs compare the heat flux generated by the√blood
flow and the heat flux due to the transfer to the surrounding air. If hs gets very large, λ K/hs
goes to zero and the skin temperature goes down to the temperature of the surrounding air.
The typical value of the parameters are the following :
— thermal conductivity of muscles : λ = 0.5 W/m.K
— metabolic heat rate production qm ˙ = 700 W/m3
3
— blood density ρb = 1000 kg/m
— blood specific heat Cb = 3600 J/kg.K
— blood flow rate per unit volume Vb = 5 × 10−4 s−1
2
— heat transfer coefficient h√s =2 W/m .K √
From these values we get 1/ K ≈ 1.7 cm and λ K/hs ≈ 1. The temperature profiles within
a limb of radius 5 cm (corresponding to a forearm), given by eqn. 15 are shown on fig. 2 for
three different values of the blood flow rate. The curve with Vb = 3 × 10−4 s−1 reproduces rather
accurately the experimental data recorded by Pennes.

4 Summary
— boundary conditions at an interface : continuity of temperature and heat fluxes. For mass
transport, chemical reactions specify the mass flux through the kinetics of the reaction
— when convective transport is present, the flux at an interface can be modeled by an effective
transport coefficient
— in a steady state one dimensional or axisymmetric heat transfer problems, we can define
thermal resistances for materials and use analogies with electrical resistance networks to
analyze the heat transfer characteristics of composite systems
— the formalism of heat transfer by diffusion, with internal production of heat, can be used
to derive a bioheat equation modeling heat transfer through living tisues