(*^
::[ Information =
"This is a Mathematica Notebook file. It contains ASCII text, and can be
transferred by email, ftp, or other text-file transfer utility. It should
be read or edited using a copy of Mathematica or MathReader. If you
received this as email, use your mail application or copy/paste to save
everything from the line containing (*^ down to the line containing ^*)
into a plain text file. On some systems you may have to give the file a
name ending with ".ma" to allow Mathematica to recognize it as a Notebook.
The line below identifies what version of Mathematica created this file,
but it can be opened using any other version as well.";
FrontEndVersion = "NeXT Mathematica Notebook Front End Version 2.2";
NeXTStandardFontEncoding;
fontset = title, inactive, noPageBreakBelow, noPageBreakInGroup, nohscroll, preserveAspect, groupLikeTitle, center, M7, bold, L1, e8, 24, "Times"; ;
fontset = subtitle, inactive, noPageBreakBelow, noPageBreakInGroup, nohscroll, preserveAspect, groupLikeTitle, center, M7, bold, L1, e6, 18, "Times"; ;
fontset = subsubtitle, inactive, noPageBreakBelow, noPageBreakInGroup, nohscroll, preserveAspect, groupLikeTitle, center, M7, italic, L1, e6, 14, "Times"; ;
fontset = section, inactive, noPageBreakBelow, nohscroll, preserveAspect, groupLikeSection, grayBox, M22, bold, L1, a20, 18, "Times"; ;
fontset = subsection, inactive, noPageBreakBelow, nohscroll, preserveAspect, groupLikeSection, blackBox, M19, bold, L1, a15, 14, "Times"; ;
fontset = subsubsection, inactive, noPageBreakBelow, nohscroll, preserveAspect, groupLikeSection, whiteBox, M18, bold, L1, a12, 12, "Times"; ;
fontset = text, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 12;
fontset = smalltext, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 10, "Times"; ;
fontset = input, noPageBreakInGroup, nowordwrap, preserveAspect, groupLikeInput, M42, N23, bold, L1, 12, "Courier"; ;
fontset = output, output, inactive, noPageBreakInGroup, nowordwrap, preserveAspect, groupLikeOutput, M42, N23, L-5, 12, "Courier"; ;
fontset = message, inactive, noPageBreakInGroup, nowordwrap, preserveAspect, groupLikeOutput, M42, N23, L1, 12, "Courier"; ;
fontset = print, inactive, noPageBreakInGroup, nowordwrap, preserveAspect, groupLikeOutput, M42, N23, L1, 12, "Courier"; ;
fontset = info, inactive, noPageBreakInGroup, nowordwrap, preserveAspect, groupLikeOutput, M42, N23, L1, 12, "Courier"; ;
fontset = postscript, PostScript, formatAsPostScript, output, inactive, noPageBreakInGroup, nowordwrap, preserveAspect, groupLikeGraphics, M7, l34, w282, h287, L1, 12, "Courier"; ;
fontset = name, inactive, noPageBreakInGroup, nohscroll, preserveAspect, M7, italic, B65535, L1, 10, "Times"; ;
fontset = header, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, italic, L1, 12, "Times"; ;
fontset = leftheader, 12;
fontset = footer, inactive, nohscroll, noKeepOnOnePage, preserveAspect, center, M7, italic, L1, 12, "Times"; ;
fontset = leftfooter, 12;
fontset = help, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 12;
fontset = clipboard, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 12;
fontset = completions, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 12, "Courier"; ;
fontset = special1, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 12;
fontset = special2, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 12;
fontset = special3, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 12;
fontset = special4, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 12;
fontset = special5, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 12;
paletteColors = 128; automaticGrouping; currentKernel;
]
:[font = title; inactive; Cclosed; preserveAspect; startGroup]
A Heat Flow Problem
:[font = subsection; inactive; Cclosed; preserveAspect; startGroup]
Discussion of the Problem
:[font = text; inactive; preserveAspect]
Our second example has its origins in a physical problem of heat flow.
Suppose we have a rod of material of unit length, say, situated on the
x-axis, for 0 <= x <= 1$. Suppose further that the rod is laterally
insulated, but has a known internal heat source. Also suppose that
both the left and right ends of the rod are held at 0 degrees Celcius.
When sufficient time passes, the temperature of the rod at each point
will settle down to ``steady state'' values, dependent only on position
x. It can be shown from physical principles that if u is the
temperature function of the rod, i.e., u(x) is the temperature of the
rod at x , then u(x) and f(x) are related by the following equations:
-u''(x) = f(x),
u(0) = 0,
u(1) = 0.
Here ' denotes derivatives with respect to x and f(x) is a known
function which depends on the material and the internal heat
distribution of the rod. The type of problem just described is called
a boundary value problem and abbreviated to BVP. Note that if we
knew how to integrate f(x) twice, it would be easy to determine
u(x) : simply use the boundary conditions at 0 and 1 to determine the
two constants of integration. The interesting case is when f(x) is
not easily integrated. In fact, a closed formula for f(x) may not
even be given. How do we solve the problem then?
We can formulate a discrete approximation to the BVP above as follows.
Choose an integer n and divide the interval [0, 1] into n+1 equal
subintervals with endpoints 0 = x_0, x_1, ... ,x_{n+1} = 1 (these
points are called nodes ). Then the width of each subinterval is
h = 1/(n+1) . Next let u_i be our approximation to u(x_i) .
Observe that we can approximate the first derivative by the forward difference formula
(1) u'(x) ~(u(x+h) - u(x)) / h .
Likewise, we can approximate the second derivative by the
backward difference formula
(2) u''(x) ~ (u'(x-h) -u'(x)) / (-h)
Now use the formula (1) with x and x-h replacing x and substitute into Formula (2) to obtain the second derivative formula
(3) u''(x) ~ (u(x+h) - 2u(x) + u(x-h)) / (h^2)
At each unknown node x_i replace the second derivative
in the equation -u''(x_i) = f(x_i) , i = 1,...,n by the approximation
given by Formula (3) ; then multiply the equation by h^2.
Observe that u(x_i + h) = u_{i + 1} and one obtains a system of
n equations in the n unknowns u_1, u_2,...,u_n, where the ith
equation takes the form
-u_{i-1} + 2u_i - u_{i+1} = h^2f(i h)
Remember that u_0 and u_{n+1} are known to be 0, so the terms u_0
and u_{n+1} disappear.
It is reasonable to expect that the smaller h is, the more accurately
u_i will approximate u(x_i) . This is indeed the case. But consider
what we are confronted with when we take n = 5 , i.e., h = 1/6 ,
which is hardly a small value of $h$. The system of five equations in
five unknowns becomes
2x_1 - x_2 = f(1/6)
-x_1 + 2x_2 - x_3 = f(2/6)
-x_2 + 2x_3 - x_4 = f(3/6)
-x_3 + 2x_4 - x_5 = f(4/6)
-x_4 + 2x_5 = f(5/6)
This problem is already about as large as we would want to work by hand.
Now let's try a specific example. Use N[Sin[Pi*x]] for f[x]. Define the required equations yourself and use Mathematica to compute the approximate solution by filling in the following cells with your equations. Use x1, x2, x3, x4, and x5 for your variables. BE SURE TO USE DECIMALS SOMEWHERE IN YOUR EQUATIONS TO FORCE MATHEMATICA TO VIEW EVERYTHING AS NUMBERS. The last cell will show you a graph of your approximate heat function u(x).
:[font = input; preserveAspect]
vars = {x1, x2, x3, x4, x5}
:[font = input; preserveAspect]
equations = { ,
,
,
,
}
:[font = input; preserveAspect]
rules = Solve[equations,vars]
:[font = input; preserveAspect; endGroup]
ListPlot[{{0.,0},{N[1/6],x1},{N[2/6],x2},{N[3/6],x3},
{N[4/6],x4},{N[5/6],x5},{1.,0.}}/.rules [[1]], PlotJoined->True]
:[font = subsection; inactive; preserveAspect; endGroup]
Exercises
^*)