0% found this document useful (0 votes)
51 views9 pages

1 Heat Equation On The Plate: Weston Barger

You are on page 1/ 9

Lecture 15

Weston Barger

1 Heat Equation on the Plate


Suppose we are given a plate with thermal conductivity σ > 0 and an initial temperature
distribution f (x , y), and we would like to find the temperature distribution u(x , y, t ) on
the plate at a later time. Suppose further that the edges of the plate are held constant.
Then we would like to solve the problem




 (x , y) ∈ (0, L) × (0, M ), 0 < t

y ∈ [0, M ], 0 < t




 u(0, y, t ) = 0,


u(L, y, t ) = 0, y ∈ [0, M ], 0 < t
  

ut = σ uxx + uyy , (1.1)



 u(x , 0, t ) = 0, x ∈ [0, L], 0 < t

x ∈ [0, L], 0 < t




 u(x , M , t ) = 0,



 u(x , y, 0) = f (x , y), (x , y) ∈ (0, L) × (0, M ).

Well, let us try separation of variables. This amounts to assuming u(x , y, t ) = X (x )Y (y)T (t ).
Inserting this form of u(x , y, t ) into the PDE in (1.1) yields

X (x )Y (y)T 0 (t ) = σX 00 (x )Y (y)T (t ) + σX (x )Y 00 (y)T (t ).

Dividing through by X (x )Y (y)T (t ) yields


X (x )Y (y)T 0 (t ) X 00 (x )Y (y)T (t ) X (x )Y 00 (y)T (t )
=σ +σ
X (x )Y (y)T (t ) X (x )Y (y)T (t ) X (x )Y (y)T (t )
00
X (x ) 00
Y (y)
=σ +σ .
X (x ) Y (y)
Now, the left and right hand side of the equation depend on distinct variables. Therefore,
for some constant λ, we have
T 0 (t ) Y 00 (y) X 00 (x )
– = = λ.
σT (t ) Y (y) X (x )

1
We also see that
T 0 (t ) Y 00 (y)
– =λ
σT (t ) Y (y)

leads to
Y 00 (y) T 0 (t )
= – λ.
Y (y) σT (t )

By the same argument, for some constant κ we have

Y 00 (y) T 0 (t )
= – λ = κ.
Y (y) σT (t )

In summary, we have the following ODEs

X 00 (x ) = λX (x ), X (0) = X (L) = 0
Y 00 (y) = κY (y), Y (0) = Y (M ) = 0
T 0 (t ) = σ(λ + κ)T (t ).

Let us start with the ODE in X . We that the boundary conditions X (0) = X (L) = 0
can only be satisfied when λ < 0. Therefore,
√  √ 
X (x ) = c1 cos –λx + c2 sin –λx .

Applying the boundary conditions gives

X (0) = c1 = 0 ⇒ c1 = 0,
√  m 2π2
X (L) = c2 sin –λL = 0 ⇒ λ=– .
L2
Thus,

mπx m 2π2
 
Xm (x ) = cm sin , λn = – , m = 1, 2, 3, · · · .
L L2
Now let us solve the ODE in Y . This is actually the same problem we solved in X
(except with boundary M instead of L). Therefore,

nπy n 2π2
 
Yn (y) = dn sin , κn = – , n = 1, 2, 3, · · · .
M M2

2
Finally, we are ready to solve the ODE in T . This has the solution
2 2 2 2 2
Tmn = emn e σ(λm +κn )t = emn e –σπ (m /L +n /M )t , m, n = 1, 2, 3, · · · .

Therefore, letting Amn = cm dn emn , for each n and m we have


mπx nπy
   
2 2 2 2 2
umn (x , y, t ) = Amn e –σπ (m /L +n /M )t sin sin .
L M
By superposition,
∞ X
∞ 2 2 2 2 2 mπx
  
nπy

Amn e –σπ (m /L +n /M )t sin
X
u(x , y, t ) = sin .
m=1 n=1 L M

We would now like to compute the coefficients Amn . Well, we do so by exploiting the or-
thogonality relationships of the functions sin(mπx /L) for m = 1, 2, 3, · · · and sin(nπy/M )
for n = 1, 2, 3, · · ·. Recall that for integers i , j = 1, 2, 3, · · ·,
Z L ! !
i πx j πx L
sin sin dx = δij
0 L L 2
Z M ! !
i πy j πy L
sin sin dx = δij .
0 M M 2

Applying the initial condition yields


∞ X
X ∞ mπx
 
nπy
 
u(x , y, 0) = f (x , y) = Amn sin sin .
m=1 n=1 L M

We choose integers i , j ≥ 1 and multiply both sides by sin(i πx /L) sin(j πy/M ) to get

i πx
!
j πy
! ∞ X
X ∞ i πx
!
j πy
!
mπx
  
nπy

f (x , y) sin sin = Amn sin sin sin sin .
L M m=1 n=1 L M L M

Integrating with respect to x over [0, L] gives


Z L
i πx
!
j πy
! ∞
X L j πy
! 
nπy

f (x , y) sin sin dx = Ain sin sin .
0 L M n=1 2 M M

Integrating with respect to y over [0, M ] gives


Z MZ L ! ! !2
i πx j πy L
f (x , y) sin sin dx dy = Aij .
0 0 L M 2

3
Therefore,

2 2Z M Z L
! !
i πx j πy
 
Aij = f (x , y) sin sin dx dy,
L 0 0 L M

which completely determines the solution (1.1). This technique generalizes to the n-
dimensional “plate” with homogeneous Dirichlet boundary conditions.
The Fourier series is defined for a function g(x , y) on [–π, π] × [–π, π] as
∞ ∞ 
1 2Z π Z π

cjk e –ijx e –iky , g(x , y)e ijx e iky dx dy.
X X
g(x , y) = cjk =
j =–∞ k =–∞
2π –π –π

It is easier to define the complex Fourier series in higher dimensions because it requires
computing less coefficients.

2 Return to shockwaves
Example 2.1. Suppose we would like to model traffic density at a red light. As before,
we denote by u(x , t ) the density of cars at a location x and time t . Let the flow of traffic
be rightwards, and suppose that at time t = 0 the end of the line of the cars backed up
our our light is at x = 0. Our initial profile is thus

 u0 , x < 0
u0 (x ) =
 u1 , x ≥ 0,

where u1 is the maximum automobile density and 0 < u0 < u1 . Denoting the speed limit
as v1 , recall that our constitutive relation was

u2
!
φ = v1 u – , (2.1)
u1

and, with no entrances or exits, our conservation law was


!
2u
ut + v 1 1 – ux = 0.
u1

Thus, our IVP is


! 
2u  –∞ < x < ∞, 0 < t
ut + v1 1 – ux = 0,
u1  u(x , 0) = u0 (x ),

4
where u1 , v1 and u0 (x ) are defined above.
Since our conservation law is in the form ut + c(u)ux = 0, our characteristic lines take
the form
!
2u (x )
x = c(u0 (x0 ))t + x0 = v1 1 – 0 0 t + x0
u1

 –v1 t + x0 , x0 ≥ 0,
=
 v (1 – 2u /u ) t + x , x < 0.
1 0 1 0 0

Now, for x0 ≥ 0, the slope of the characteristic lines is –1/v1 . Since 0 < u0 < u1 , the
slope of the characteristic lines for x0 < 0 is strictly between –1/v1 and 1/v1 . For x0 < 0,
if the slope is positive, then it is clear that tb = 0. Now, if the slope is negative, it is
strictly greater that –1/v1 , so we still have that tb = 0. Thus, tb = 0 and the first crossing
occurs at the origin x = 0. Therefore, our curve (xs (t ), t ) will have the initial condition
xs (0) = 0.
Using the RH condition, we see that

dxs φ(u(xs (t )– , t )) – φ(u(xs (t )+ , t ))


=
dt u(xs (t )– , t ) – u(xs (t )+ , t )
φ(u0 ) – φ(u1 )
=
u0 – u1
φ(u1 ) – φ(u0 )
= .
u1 – u0

By (2.1), φ(u1 ) = 0 and φ(u0 ) = v1 (u0 – u02 /u1 ). Therefore, the IVP for xs (t ) is
 
dxs 0 – v1 u0 – u02 /u1 u
= = –v1 0 , xs (0) = 0.
dt u1 – u0 u1
We obtain that
u
xs (t ) = –v1 0 t .
u1
Therefore,

 u0 , x < –v1 (u0 /u1 )t
u(x , t ) =
 u1 , x ≥ –v1 (u0 /u1 )t .

5
3 Viscosity solutions
Another approach for locating shock waves is through the method of viscosity solutions.
We make a model modification to the traffic model in which we introduce a “smoothing
parameter.” As we let this parameter tend to 0 we recover our shock waves from the
previous section.

3.1 Another traffic flow model


Suppose now that we would like to achieve “smoother” behavior. One way to do this is
it account for a driver’s ability to see higher or lower density traffic ahead and adjust
accordingly. If at a time t > 0, a driver at position x sees that the density of traffic
is higher then ux (t , x ) > 0. Similarly, if traffic is clearing out ahead of the driver then
ux (x , t ) < 0. Then we could modify our velocity function v (u) function to be
!
u ux
v (u) = v1 1 – –r ,
u1 u
where r ≥ 0 is a parameter. Recalling that φ = vu, we have
u2
!
φx (u) = v1 u – – rux .
u1
Therefore, our conservation law ut + φx = 0 is now a the second order equation
!
2u
ut + v1 1 – ux = ruxx .
u1
When r = 0, our previous model is recovered. We place the assumptions on the model
that far in front of the end of the traffic line, the car density is nearly constant u1 . Far
behind the pile of cars waiting at the light the traffic has almost constant density u0 . This
is akin to assuming

lim u(x , 0) = u1 , lim u (x , 0) = 0, (3.1)


x →∞ x →∞ x
lim u(x , 0) = u0 , lim u (x , 0) = 0. (3.2)
x →–∞ x →–∞ x

3.2 Traveling wave solutions


We look for traveling wave solutions to our new traffic flow model
!
2u
ut + v1 1 – ux = ruxx . (3.3)
u1

6
We look for traveling wave solutions to (3.3). Assuming that for some function f (z ),
u(x , t ) = f (x – ct ) and inserting into (3.3), we get
2v1 0
–cf 0 + v1 f 0 – ff = rf 00 .
u1
Integrating gives
v
–cf + v1 f – 1 f 2 = rf 0 + k , (3.4)
u1
where k is our integration constant.
We can solve for c and k by imposing our limit conditions (3.1) and (3.2). Since
u(x , 0) = f (x ), we get

lim f (z ) = u1 , lim f 0 (z ) = 0,
x →∞ x →∞
lim f (z ) = u0 , lim f 0 (z ) = 0.
x →–∞ x →–∞

Taking z → ∞ in (3.4) yields

–cu1 + v1 u1 – v1 u1 = k ,

and taking z → –∞ in (3.4) gives

v u2
–cu0 + v1 u0 – 1 0 = k .
u1
Solving these two algebraic equations for c and k gives
u
c = –v1 0 , k = v1 u0 .
u1
With these values of c and k , the differential equation (3.4) becomes
u v
v1 0 f + v1 f – 1 f 2 = rf 0 + v1 u0 .
u1 u1
We rearrange the above equation to get
ru1 0
(f – u0 )(f – u1 ) = – f .
v1
We rearrange the above equation to get
f0 v
=– 1 .
(f – u0 )(f – u1 ) ru1

7
Figure 1: A wave profile of the viscosity solution (3.5).

Using partial fractions, we integrate both sides of the equation to get



1 f – u
1 v1
log =– z + k1 ,

u1 – u0 f – u0 ru1

where k1 is the constant of integration. We make the assumption that u0 < f < u1 . This
makes sense with our model, since f represents the density of cars. With this additional
assumption, we have
" #
1 u –f v
log 1 = – 1 z + k1 .
u1 – u0 f – u0 ru1

We take k1 = 0 and solve for f to get


 
(u –u )v
u1 + u0 exp – 1 ru10 1 z u0 – u1
f (z ) =   = u1 +  ,
(u –u )v v1 (u1 –u0 )
1 + exp – 1 ru10 1 z 1 + exp ru1 z

where the second equality comes from adding and subtracting u1 exp[∗] in the numerator
of the first fraction. Therefore, recalling c = –v1 (u0 /u1 ), we get
u0 – u1
u(x , t ) = u1 +   . (3.5)
v1 (u1 –u0 ) v1 u0

1 + exp ru1 x+ u1 t

A wave profile for (3.5) is plotted in Figure 1.


The term ruxx added to form the traffic flow equation is called the viscosity term
with viscosity parameter r . The result (3.5) is called the viscosity solution.
Note that both the solution we found in previous section and our viscosity solution are
wave fronts moving against the flow of incoming traffic with speed v0 u0 /u1 . Their profile

8
shapes are similar too, except that the viscosity solution does not have a discontinuity.
The viscosity term had the effect of “smoothing out” the discontinuous shock present in
solutions of the first traffic equation.
We note that
 
 u0 – u1 
lim u(x , t ) = lim u1 +   
r →0+ r →0+ v1 (u1 –u0 )
x + vu1 u1 0 t
 
1 + exp ru1

 u0 , x < –v1 (u0 /u1 )t
=
u1 , x ≥ –v1 (u0 /u1 )t .

Therefore, as we let r → 0+ , we recover our shock wave solution.

You might also like