Douglas Wilhelm Harder MMath LEL Department of Electrical and Computer Engineering University of Waterloo Waterloo Ontario Canada eceuwaterlooca dwharderalumniuwaterlooca 2012 by Douglas Wilhelm Harder Some rights reserved ID: 797357
Download The PPT/PDF document "Heat-conduction/Diffusion Equation" is the property of its rightful owner. Permission is granted to download and print the materials on this web site for personal, non-commercial use only, and to display it on your personal computer provided you do not modify the materials and that you retain all copyright notices contained in the materials. By downloading content from our website, you accept the terms of this agreement.
![Chapter 3: Unsteady State [ Transient ] Heat Conduction](https://thumbs.docslides.com/335260/chapter-3-unsteady-state-transient-heat-conduction.jpg)
Slide1
Heat-conduction/Diffusion Equation
Douglas Wilhelm Harder,
M.Math
. LEL
Department of Electrical and Computer Engineering
University of Waterloo
Waterloo, Ontario, Canada
ece.uwaterloo.ca
dwharder@alumni.uwaterloo.ca
© 2012 by Douglas Wilhelm Harder. Some rights reserved.
Slide2Outline
This topic discusses numerical solutions to the heat-conduction/ diffusion equation:
Background
Discuss the physical problem and propertiesExamine the equationApproximate solutions using a finite-difference equationConsider numerical stabilityExamples
2
Heat-conduction/Diffusion Equation
Slide3Outcomes Based Learning Objectives
By the end of this laboratory, you will:
Understand the heat-conduction/diffusion equation
Understand how to approximate partial differential equations using finite-difference equationsSet up solutions in one spatial dimension3
Heat-conduction/Diffusion Equation
Slide4Background
In the last laboratory, we looked finding solutions to boundary-value problems where we were given:
An ordinary-differential equation,
An interval [a,
b]
, and
Two boundary conditions
4
Heat-conduction/Diffusion Equation
Slide5Background
To approximate the solution:
We broke the interval interval
[a, b
] into n
points
x
i
We approximated
u
i
≈
u
(xi)We used a finite-difference equation to create a system of linear equations in the unknowns u2, u3, ..., un – 1
5
Heat-conduction/Diffusion Equation
Slide6The Heat-conduction/Diffusion Equation
The conduction of heat or diffusion of a material depends both on space and time
We have a partial differential equation in both space and time
In one dimension, this simplifies to
Here
k
> 0
is the
diffusivity coefficient
6
Heat-conduction/Diffusion Equation
Slide7Motivating Examples
Suppose a metal bar is in contact with a body at
80
oCIf the bar is insulated, over time, the entire length of the bar will be at 80
oC
80
o
C
5
o
C
7
Heat-conduction/Diffusion Equation
Slide8Motivating Examples
At time
t =
0, one end of the bar is brought in contact with a heat sink at 5 oC
80
o
C
5
o
C
8
Heat-conduction/Diffusion Equation
Slide9Motivating Examples
At the end closest to the heat sink, the bar will cool rapidly; however, the cooling is not uniform
80
o
C
5
o
C
9
Heat-conduction/Diffusion Equation
Slide10Motivating Examples
Over time, however, the temperature will drop linearly from
80 oC
to 5
o
C
across the bar
80
o
C
5
o
C
10
Heat-conduction/Diffusion Equation
The description of this will be a function u(x
, t)
Slide11Motivating Examples
Other applications—often in higher dimensions—include:
Heat equation
Fick’s laws of diffusionThermal diffusivity in polymers
Steven Byrnes
Oleg
Alexandrov
11
Heat-conduction/Diffusion Equation
Slide12Motivating Examples
It can even be used in modeling human vision
Diffusion is used for enhancing coherence in images
Scale space
J. Weickert and ter Haar Romeny
12
Heat-conduction/Diffusion Equation
Slide13Laplace’s Equation
Last week, we saw a related ordinary-differential equation:
This is the one-dimensional version of Laplace’s equation:
13
Heat-conduction/Diffusion Equation
Slide14Laplace’s Equation
In Laboratory 1, we saw that the only solution to Laplace’s equation in one dimension with boundary values is a straight line
14
Heat-conduction/Diffusion Equation
Slide15Laplace’s Equation
Now, suppose that
u
(x, t)
satisfies Laplace’s equation with respect to x
:
What does this say about heat-conduction or diffusion?
That is, the solution is in a steady state
It does not change with respect to time
15
Heat-conduction/Diffusion Equation
Slide16Rate of Change
Proportional to Concavity
Suppose, however, that our function
u(x,
t) does not satisfy Laplace’s equation…
In one dimension, what does the equation
mean?
Recall we assumed that
k
> 0
16
Heat-conduction/Diffusion Equation
Slide17Rate of Change
Proportional to Concavity
Visually:
Where
u
is concave up, the rate of change over time will be positive
17
Heat-conduction/Diffusion Equation
Slide18Rate of Change
Proportional to Concavity
Visually:
Where
u
is concave down, the rate of change over time is negative
18
Heat-conduction/Diffusion Equation
Slide19Rate of Change
Proportional to Concavity
Visually:
Of course, a function may be concave up and down in different regions
Ultimately, the concavity tends to disappear—that is, as time goes to infinity,
u
(
x
,
t
)
becomes a straight line…
19
Heat-conduction/Diffusion Equation
Slide20Examples
Consider heat conduction problem:
A bar initially at
5 oC with one end touching a heat sink also at 5
oC
At time
t
= 0
, the other end is placed in contact with a heat source at
80
o
C
20
Heat-conduction/Diffusion Equation
Slide21Examples
As time passes, the bar warms up
21
Heat-conduction/Diffusion Equation
Slide22Examples
A plot of the temperature over time across the bar is given by this solution:
The end placed in contact with the
80 oC heats up faster than points closer to the heat sink
As the concavity becomes smaller, the rate of change also gets smaller
22
Heat-conduction/Diffusion Equation
Slide23Examples
Consider another example:
A bar initially at
100 o
C is placed at time t
= 0
in contact with a heat source of
80
o
C
at one end and a heat sink at
5
o
C
23
Heat-conduction/Diffusion Equation
Slide24Examples
The ultimate temperature will be no different:
A uniform change from
80 oC down to 5
oC
24
Heat-conduction/Diffusion Equation
Slide25Examples
A plot of the temperature over time across the bar is given by this solution:
Both ends drop from
100 oC however, the centre cools slower than does the bar at either end point
25
Heat-conduction/Diffusion Equation
Slide26Rate of Change
Proportional to Concavity
Ultimately, the concavity tends to disappear—that is, as
t → ∞,
u(
x
,
t
)
becomes a straight line…
Important note:
This statement is
only
true
in one
spacial dimensionIn higher spacial dimensions, u(x, t) will approach a solution to Laplace’s equation in the given region, but it won’t be so simple26Heat-conduction/Diffusion Equation
Slide27The Problem
For any such problem, we know the initial state:
we would like to find a function
u
(
x
,
t
)
that satisfies the partial differential equation up to some point
t
final
27Heat-conduction/Diffusion Equationtfinal
Slide28The Problem
To do this, we must also know the boundary values:
28
Heat-conduction/Diffusion Equation
t
final
Slide29Approximating the Solution
As with the BVP, we will divide the space interval into discrete points
We are interested, however, also in how u(
x,
t
)
changes over time
Just like the IVPs, we will divide time into discrete steps
29
Heat-conduction/Diffusion Equation
Slide30Approximating the Solution
This creates an
nx
× n
t
grid
We will approximate
u
(
x
i
,
t
k
)
≈ ui,k30Heat-conduction/Diffusion Equationtfinal
Slide31Approximating Partial Derivatives
Questions:
How do we approximate partial derivatives?
What are the required conditions?31Heat-conduction/Diffusion Equation
Slide32Approximating Partial Derivatives
Recall our approximations of the derivatives:
In this case, however, we are dealing with partial derivatives
32
Heat-conduction/Diffusion Equation
Slide33Approximating Partial Derivatives
If you recall the definition of the partial derivative, we focus on one variable and leave the other variables constant
To ensure clarity:
h
is used to indicate small steps in space
D
t
is used to indicate a small step in time
33
Heat-conduction/Diffusion Equation
Slide34Approximating Partial Derivatives
The obvious extension is to define
or
34
Heat-conduction/Diffusion Equation
Forward divided-
difference formulas
Centred divided-
difference formulas
Slide35Approximating Partial Derivatives
Recall, however, that we will be approximating the solution on a grid
u
(xi, t
k)
≈
u
i,k
where
Thus, our formulas turn out to be
35
Heat-conduction/Diffusion Equation
Recall that
x
i
±
h
=
x
i
±
1
and
t
k
±
D
t
=
t
k
± 1
Slide36Approximating Partial Derivatives
Thus, our approximations are:
Similarly,
36
Heat-conduction/Diffusion Equation
Slide37Approximating Partial Derivatives
Problem:
We are only given the state at a single initial time
t0Thus, we cannot use three unknowns in time:
u
i,k
– 1
u
i,k
u
i,k
+ 1
Thus, we will use the formula37Heat-conduction/Diffusion Equation
Slide38Approximating Partial Derivatives
If we substitute the second into the heat-conduction/ diffusion equation, we have
Solving for the event in the future, u
i,
k
+ 1
:
38
Heat-conduction/Diffusion Equation
Slide39Similarity with IVPs
If this looks familiar to you, it should:
Recall Euler’s method: if
y(1)
(t
) =
f
(
t
,
y
(
t
))
and
y
(t0) = y0, it follows that
The slope at
(
t
k
,
y
k
)
39
Heat-conduction/Diffusion Equation
Slide40Approximating the
Second Partial Derivative
The next step is to approximate the concavity
We will use our 2nd-order formula:
to get
40
Heat-conduction/Diffusion Equation
What if this value is very large?
E.g.
,
h
is small
Slide41Approximating the
Second Partial Derivative
This is our finite-difference equation approximating our partial-differential equation:
Graphically, we see
that
u
i
,
k
+ 1
depends on
three previous values:
41
Heat-conduction/Diffusion Equation
Slide42Initial and Boundary Conditions
For a 1
st
-order ODE, we require a single initial condition42
Heat-conduction/Diffusion Equation
Slide43Initial and Boundary Conditions
For a 2
nd
-order ODE, we require either two initial conditions or a boundary condition:
Time variables are usually
associated with initial conditions
Space variables are usually
associated with boundary conditions
43
Heat-conduction/Diffusion Equation
Slide44Initial and Boundary Conditions
For the heat-conduction/diffusion equation,
we have:
A first-order partial derivative with respect to timeA second-order partial derivative with respect to space
44
Heat-conduction/Diffusion Equation
Slide45Initial and Boundary Conditions
For the heat-conduction/diffusion equation,
we require:
For each space coordinate, we require an I.C. for the timeFor each time coordinate, we require a B.C. for the space
45
Heat-conduction/Diffusion Equation
Slide46Initial and Boundary Conditions
In this case, we will require an initial value for each space coordinate
Keeping in the spirit of one dimension:
If we were monitoring the progress of the temperature of a bar, we would know the initial temperature at each point of the bar attime t =
t0
If we were monitoring the diffusion of a gas , we would have to know the initial concentration of the gas
46
Heat-conduction/Diffusion Equation
Slide47Initial and Boundary Conditions
Is a one-dimensional system restricted to bars and thin tubes?
Insulated wires, impermeable tubes,
etc. May systems have symmetries that allow us to simplify the model of the systemConsider the effect of a thin insulator between two
materials that have approximately uniform
temperature
Consider the diffusion of particles across a
membrane separating liquids with different
concentrations
47
Heat-conduction/Diffusion Equation
Slide48Approximating the Solution
Just as we did with BVPs, we will divide the spacial interval
[
a, b] into
nx
points or
n
x
– 1
sub-intervals
48
Heat-conduction/Diffusion Equation
Slide49Approximating the Solution
Just as we did with BVPs, we will divide the spacial interval
[
a,
b] into
n
x
points or
n
x
– 1
sub-intervals
The initial state at time
t
0
would be defined by a function uinit(x) whereuinit:[a, b] → R
49
Heat-conduction/Diffusion Equation
Slide50Approximating the Solution
As with Euler’s method and other IVP solvers, we divide the time interval
[
t0
, t
final
]
into
n
t
points or
n
t
– 1
sub-intervals
Thus, our t-values will be t = (t1, t2, t3, ..., t
n )
50
Heat-conduction/Diffusion Equation
t
Slide51Approximating the Solution
As with Euler’s method, we will attempt to approximate the solution
u
(x, t)
at discrete times in the futureGiven the state at time
t
1
, approximate the solution at time
t
2
51
Heat-conduction/Diffusion Equation
t
1
Slide52Approximating the Solution
Never-the-less, we must still determine the 2
nd
-order component
This requires two boundary values at
a
and
b
52
Heat-conduction/Diffusion Equation
t
1
Slide53Approximating the Solution
Indeed, at each subsequent point, we will require the boundary values
Over time, the boundary values could change
53
Heat-conduction/Diffusion Equation
t
1
Slide54Approximating the Solution
Therefore, we will need two functions
a
bndry(t) and
bbndry
(
t
)
At each step, we would evaluate and determine two the boundary conditions
54
Heat-conduction/Diffusion Equation
t
1
Slide55Approximating the Solution
Therefore, we will need two functions
a
bndry(t) and
bbndry
(
t
)
We would continue from our initial time
t
0
and continue to a final time
t
final
55
Heat-conduction/Diffusion Equationt1
Slide56Approximating the Solution
Therefore, we will need two functions
a
bndry(t) and
bbndry
(
t
)
We would continue from our initial time
t
0
and continue to a final time
t
f
56
Heat-conduction/Diffusion Equationt1
Slide57Approximating the Solution
Now that we have functions that define both the initial conditions and the boundary conditions, we must solve how the heat conducts or the material diffuses
57
Heat-conduction/Diffusion Equation
Slide58Approximating the Solution
Given the state at time
t1, we approximate the state at time t2
58
Heat-conduction/Diffusion Equation
t
1
Slide59Approximating the Solution
Given the state at time
t2, we approximate the state at time t3
59
Heat-conduction/Diffusion Equation
t
1
Slide60Approximating the Solution
Given the state at time
t3, we approximate the state at time t4
60
Heat-conduction/Diffusion Equation
t
1
Slide61Approximating the Solution
We continue this process until we have approximated the solution at each of the desired times
61
Heat-conduction/Diffusion Equation
t
1
Slide62Approximating the Solution
This last image suggests strongly that our solution will be a matrix
Given the heat-conduction/diffusion equation
we will approximate the solution u
(x,
t
)
where
a
<
x
<
b
and
t
0 < t ≤ tfinal by finding an nx × nt matrix U
where U(i
, k) (that is,
u
i
, k
)
will approximate the solution at time
(
x
i
,
t
k
)
with
x
i
=
a
+ (
i
– 1)
h
t
k
=
t
0
+ (
k
– 1)
D
t
62
Heat-conduction/Diffusion Equation
Slide63The diffusion1d
Function
The signature will be
function [
x_out,
t_out
,
U_out
] = ...
diffusion1d( kappa,
x_rng
,
nx
,
t_rng
, nt,
u_init, u_bndry ) where kappa the diffusivity coefficient x_rng the space range
[
a,
b
]
nx
the number of points into which we will divide
[
a
,
b
]
t_rng
the time interval
[
t
0
,
t
final
]
nt
the number of points into which we will divide
[
t
0
,
t
final
]
u_init
a function handle giving the initial state
u
init
:
[
a
,
b
]
→
R
u_bndry
a function handle giving the two boundary conditions
u
bndry
:
[
t
0
,
t
final
]
→
R
2
where
63
Heat-conduction/Diffusion Equation
Slide64Step 1: Error Checking
As with the previous laboratories, you will have to validate the types and the values of the arguments
Recall, however, that we also asked:
What happens if
is too large?A large coefficient could magnify any error in our approximation:
64
Heat-conduction/Diffusion Equation
Slide65Step 1: Error Checking
To analyze this, consider the total sum:
The two approximations of the slope have an
O(h
)
error
To eliminate the contribution from both these errors, we will therefore require that
65
Heat-conduction/Diffusion Equation
Slide66Step 1: Error Checking
Once the parameters are validated, the next step is to ensure
If this condition is not met, we should exit and provide a useful error message to the user:
For example,
The ratio kappa*dt
/h^2 = ??? >= 0.5, consider using
n_t
= ???
where
The first
???
is replaced by the calculated ratio, and
The second
???
is found by calculating the smallest integer for
n
t that would be acceptable to bring this ratio under 0.5You may wish to read up on the floor and ceil commands66Heat-conduction/Diffusion Equation
Slide67Step 1: Error Checking
Essentially, given
k
, t0, tfinal
and h
, find an appropriate value (that is, the smallest integer value) of
n
t
*
to ensure that
Hint: solve for
n
t
*
and choose the next largest integer...
You may wish to review the
floor and ceil commands67Heat-conduction/Diffusion Equation
Slide68Step 2: Initialization
In order to calculate both the initial and boundary values, we will require two vectors:
A column vector of
x values, andA row vector of t values
Both of these can be constructed using
linspace
Both of these will be returned to the user
We will begin by creating the
n
x
×
n
t
matrixThe zeros command is as good as any68Heat-conduction/Diffusion Equation
Slide69Step 2: Initialization
Next, we assign the initial values to the first column
n
x
69
Heat-conduction/Diffusion Equation
Slide70Step 2: Initialization
Next, assign the boundary values in the first and last rows
n
x
n
t
n
t
– 1
70
Heat-conduction/Diffusion Equation
Slide71Step 3: Solving
We still have to calculate the interior values
n
x
n
t
71
Heat-conduction/Diffusion Equation
Slide72Step 3: Solving
For this, we go back to our formula:
72
Heat-conduction/Diffusion Equation
Slide73Step 3: Solving
Using this formula
we can first calculate
u
2,2 through
u
n
– 1,2
73
Heat-conduction/Diffusion Equation
x
Slide74Step 3: Solving
Using this formula
and then
u
2, 3 through u
n
– 1, 3
and so on
x
74
Heat-conduction/Diffusion Equation
Slide75Useful Matlab Commands
Recall that you can both access and assign to sub-matrices of a given matrix:
>> U = [1 2 3 4 5; 6 7 8 9 10; 11 12 13 14 15; 16 17 18 19 20]
U = 1 2 3 4 5 6 7 8 9 10
11 12 13 14 15
16 17 18 19 20
>> U(:,1) = [
0 0 0 0
]';
>> U(1,2:end) = [
42 42 42 42
];
>> U(2:end - 1,3) = [
91 91
]'
U =
0 42 42 42 42 0 7 91
9 10
0
12
91
14 15
0
17 18 19 20
75
Heat-conduction/Diffusion Equation
Slide76Useful Matlab Commands
Matlab actually gives you a lot of power when it comes to assigning entries by row or by column:
>> U = [1 2 3 4 5; 6 7 8 9 10; 11 12 13 14 15; 16 17 18 19 20]
U = 1 2 3 4 5 6 7 8 9 10
11 12 13 14 15
16 17 18 19 20
>> U([1, 3], 2:end) = [
42 42 42 42
;
91 91 91 91
]
U =
1
42 42 42 42
6 7 8 9 10
11 91 91 91 91 16 17 18 19 2076Heat-conduction/Diffusion Equation
Slide77Useful Matlab Commands
The
throw
command in Matlab allows the user to format its message in a manner similar to the sprintf command in Matlab
(and C) function [] =
error_example
( a, b )
throw(
MException
( '
MATLAB:invalid_argument
', ...
'the argument
%f
(
%d
) and %f (%d) were passed', ... a, round(a), b, round(b) ) );
end
When run, the output is:
>>
error_example
( 3, pi )
??? Error using ==> ff at 2
the argument 3.000000 (3) and 3.141593 (3) were passed
77
Heat-conduction/Diffusion Equation
Slide78Useful Matlab Commands
You can, of course, modify this:
function [] =
error_example( a, b ) throw( MException( 'MATLAB:invalid_argument
', ... 'the ratio of
%d
and
%d
is
%f
',
a
,
b
,
a/b
) ); end When run, the output is:>> error_example( 3, 5 )??? Error using ==> error_example at 2the ratio of 3 and 5 is 0.600000
See >> help
sprintf
78
Heat-conduction/Diffusion Equation
Slide79Useful Matlab Commands
A very useful command is the
diff
command:>> v = [1 3 6 7 5 4 2 0]v = 1 3 6 7 5 4 2 0
>> diff( v )
ans
=
2 3 1 -2 -1 -2 -2
>> v(2:end) - v(1:end - 1)
ans
=
2 3 1 -2 -1 -2 -2
>> diff( v, 2 )
ans
=
1 -2 -3 1 -1 0
>> v(3:end) - 2*v(2:end - 1) + v(1:end - 2)ans = 1 -2 -3 1 -1 079Heat-conduction/Diffusion Equation
Slide80Useful Matlab Commands
Why use
diff
instead of a for loop?We could simply use two nested for loops:>> for k = 2:nt for ix = 2:(
nx - 1)
% calculate U(
i
, k)
end
end
however, it would be better to use:
>> for k = 2:nt
% calculate U(2:(
nx
– 1), k) using diff
end
In MatlabMost commands call compiled C codeFor loops and while loops are interpreted80Heat-conduction/Diffusion Equation
Slide81Useful Matlab Commands
How slow is interpreted code?
This function implements matrix-matrix-multiplication
function [M_out] = multiplication( M1, M2 ) [m, n1] = size( M1 );
[n2, p] = size( M2 );
if n1 ~= n; error( 'dimensions must agree' ); end;
M_out
= zeros( m, p );
for
i
=1:m
for j=1:n1
for k = 1:p
M_out
(i, k) = M12(i, k) + M1(i, j)*M2(j, k); end end endend
81
Heat-conduction/Diffusion Equation
Slide82Useful Matlab Commands
Now, consider:
>> M = rand( 1300, 1000 );
>> N = rand( 1000, 1007 ); >> tic; M*N; toc
Elapsed time is
0.196247
seconds.
>> tic; multiplication( M, N );
toc
Elapsed time is
35.530549
seconds.
>> 35.530549/0.196247
ans
=
181.0502
What’s two orders of magnitude between you and your employer?82Heat-conduction/Diffusion Equation
Slide83Useful Matlab Commands
Also, what is faster? Using
diff
or calculating it explicitly?>> v = rand( 1, 1000000 );>> tic; diff( v, 2 ); toc
Elapsed time is 0.009948 seconds.
>> tic; v(3:end) - 2*v(2:end - 1) + v(1:end - 2);
toc
Elapsed time is 0.039661 seconds.
At this point, it’s only a factor of four...
83
Heat-conduction/Diffusion Equation
Slide84Useful Matlab Commands
Suppose we want to implement a function which is only turned on for, say, one period
In
Matlab
, you could use the following:
function [y] = f( x )
y = -
cos
(x) .* ((x >= 0.5*pi) & (x <= 2.5*pi ));
end
84
Heat-conduction/Diffusion Equation
Slide85Useful Matlab Commands
The output of
>> x =
linspace( 0, 10, 101 );
>> plot( x, f( x ), '
r
o' )
In order to understand the output, consider
>>
x =
linspace
( 0, 10, 11 );
>> x >= 0.5*pi
ans =
0 0 1 1 1 1 1 1 1 1 1
>> x <= 2.5*pi
ans = 1 1 1 1 1 1 1 1 0 0 0>> (x >= 0.5*pi) & (x <= 2.5*pi)ans = 0 0 1 1 1 1 1 1 0 0 085
Heat-conduction/Diffusion Equation
Slide86Useful Matlab Commands
Functions may be defined in one of two ways
Functions, and
Anonymous functions Create a new file with File→New→Function:
function [u] = u2a_init( x )
u = (x >= 0.5).*(1 + 0*x);
end
Saving this function to a file
u2a_init.m
:
>> x =
linspace
( 0, 1, 11 )
x =
0 0.1 0.2 0.3 0.4
0.5 0.6 0.7 0.8 0.9 1.0
>> u2a_init ( x )ans = 0 0 0 0 0 1 1 1 1 1 1>> isa( @u2a_init, 'function_handle' );ans =
1
86
Heat-conduction/Diffusion Equation
Slide87Useful Matlab Commands
Functions may be defined in one of two ways
Functions, and
Anonymous functions Anonymous functions may be defined as follows:>> u
2b_init = @(x)((x >= 0.5).*(1 + 0*x))
u2b_init =
@(x)((x>=0.5).*(1+0*x))
>> x =
linspace
( 0, 1, 11 )
x =
0 0.1 0.2 0.3 0.4
0.5 0.6 0.7 0.8 0.9 1.0
>> u2b_init ( x )
ans =
0 0 0 0 0 1 1 1 1 1 1>> isa( u2b_init, 'function_handle' ); ans = 1
87
Heat-conduction/Diffusion Equation
Slide88Example 1
As a first example:
>> [x1, t1, U1] = diffusion1d(
0.1, [0,1],
6,
[1, 3]
,
12
, @u2a_init, @u2a_bndry );
>>
mesh
( t1, x1, U1 )
>> size(x1)
ans =
6 1
>> size(t1)
ans = 1 12>> size(U1)ans = 6 12
x
t
b
= 1
a
= 0
t
0
= 1
t
final
= 3
u
init
(
x
) = 0.9
b
bndry
(
t
) = 4.1
a
bndry
(
t
) = 1.7
n
t
= 12
n
x
= 6
k
= 0.1
88
Heat-conduction/Diffusion Equation
Slide89Example 1
The boundary conditions are specified by functions:
>> [x1, t1, U1] = diffusion1d(
0.1, [0,1],
6,
[1, 3]
,
12
, @u2a_init, @u2a_bndry );
function
[u] = u2a_init ( x )
u = x*0 + 0.9;
end
function
[u] = u2a_bndry ( t )
u = [t*0 + 1.7; t*0 + 4.1];end
x
t
b
= 1
a
= 0
t
0
= 1
t
final
= 3
u
init
(
x
) = 0.9
b
bndry
(
t
) = 4.1
a
bndry
(
t
) = 1.7
n
t
= 12
n
x
= 6
k
= 0.1
89
Heat-conduction/Diffusion Equation
Slide90Example 1
The boundary conditions are specified by functions:
>> [x1, t1, U1] = diffusion1d(
0.1, [0,1],
6,
[1, 3]
,
12
, @u2a_init, @u2a_bndry );
>> x =
linspace
(
0
,
1
,
6 )';>> u2a_init ( x )ans = 0.9 0.9 0.9 0.9 0.9 0.9
>> t =
linspace(
1
,
3
,
12
);
>> u2a_bndry( t(2:end) )
ans =
1.7 1.7 1.7 1.7 1.7 1.7 1.7 1.7 1.7 1.7 1.7
4.1 4.1 4.1 4.1 4.1 4.1 4.1 4.1 4.1 4.1 4.1
90
Heat-conduction/Diffusion Equation
Slide91Example 1
Because the large size of both
h
and Dt
, there is still some numeric fluctuations
x
t
91
Heat-conduction/Diffusion Equation
Slide92Example 2
Decreasing
h
and Dt produces a better result>> [x2, t2, U2] = diffusion1d( 0.1, [0,1], 100, [1, 3], 4000, @u2a_init, @u2a_bndry );
>> size(x2)
ans =
100 1
>> size(t2)
ans =
1 4000
>> size(U2)
ans =
100 4000
>>
mesh
( t2, x2, U2 )
x
t
h
= 0.01
D
t
= 0.0005
92
Heat-conduction/Diffusion Equation
Slide93Example 3
If we exceed
0.5
, the approximation begins to fail>> [x2, t2, U2] = diffusion1d( 0.1, [0,1], 100, [1, 3], 3916, @u2a_init, @u2a_bndry );
>>
mesh
( t2, x2, U2 )
x
t
Stable when
n
t
= 4000
93
Heat-conduction/Diffusion Equation
Slide94Example 4
As the ratio exceeds
0.5
, it quickly breaks down:>> [x2, t2, U2] = diffusion1d( 0.1, [0,1], 100, [1, 3], 3900, @u2a_init, @u2a_bndry );
>>
mesh
( t2, x2, U2 )
x
t
Stable when
n
t
= 4000
2
×
10
14
94
Heat-conduction/Diffusion Equation
Slide95Animations
Now, these images show us
u
(x, t)
, but what about visualizing how the temperature changes over time?
In the source directory of the laboratory, you will find:
animate.m
frames2gif.m
isosurf.m
rgb2gif8bit.m
>> [x3, t3, U3] = diffusion1d( 0.1, [0,1], 20, [1, 3], 150, ...
@u2a_init, @u2a_bndry );
>> frames =
animate
( U3 ); %
Animate
the output file>> frames2gif( frames, 'U3.gif' ); % Save it as an animated gif
95Heat-conduction/Diffusion Equation
Slide96Animations
Now we can visualize the state as it changes over time:
96
Heat-conduction/Diffusion Equation
mesh
( t3, x3, U3 );
frames =
animate
( U3 );
frames2gif( frames, 'U3.gif' );
Slide97Summary
We have looked at the heat-conduction/diffusion equation
Considered the physical problem
Considered the effects in one dimensionFound a finite-difference equation approximating the PDESaw that
ExamplesAnimations
97
Heat-conduction/Diffusion Equation
Slide98References
[1] Glyn James, Modern Engineering Mathematics, 4
th
Ed., Prentice Hall, 2007, p.778.[2] Glyn James, Advanced Modern Engineering Mathematics, 4th Ed., Prentice Hall, 2011, p.164.
98
Heat-conduction/Diffusion Equation