Applications - A Guide to MATLAB for Beginners and Experienced Users - page 188

Here we used surf to show the entire solution u ( x , t ). The output is clearly
unrealistic; notice the scale on the u axis! The numerical solution of partial
differential equations is fraught with dangers, and instability like that seen
above is a common problem withfinite difference schemes.For many partial
differential equations a finite difference scheme will not work at all, but for
the heat equation and similar equations it will work well with proper choice
of t and x . One might be inclined to think that since our choice of x was
larger, it should be reduced, but in fact this would only make matters worse.
Ultimately the only parameter in the iteration we're using is the constant s ,
and one drawback of doing all the computations in an M-file as we did above
is that we do not automatically see the intermediate quantities it computes.
In this case we can easily calculate that s = 2(0 . 1) / (0 . 5) 2 = 0 . 8. Notice that

this implies that the coefficient 1 − 2 s of u j in the iteration above is negative.
Thus the “weighted average” we described before in our interpretation of the
iterative step is not a true average; eachsection of wire is transferring more
energy than it has at each time step!
The solution to the problem above is thus to reduce the time step t ; for
instance, if we cut it in half, then s = 0 . 4, and all coefficients in the iteration
are positive.

tvals = linspace(0, 4, 81);
uvals = heat(2, xvals, tvals, init, [15 25]);
surf(xvals, tvals, uvals)
xlabel x; ylabel t; zlabel u

25

20

15

4

3 5

2

0

1

0 − 5

t

x