\documentstyle[12pt,qss12]{article}
\hoffset=-1.025in %0.5in at left if 12pt
\voffset=-1.150in %0.75in at top if 12pt
%\hoffset=-0.525in %1.0in at left if 12pt
%\voffset=-0.900in %1.0in at top if 12pt
%\hoffset=-0.025in %1.5in at left if 12pt
%\voffset=-0.65in %1.25in at top if 12pt
\setlength{\textheight}{9.65in}
\setlength{\textwidth}{7.5in}
\newcommand{\PI}{\mbox{\large$\pi$}}
\begin{document}
\begin{center}
{\large\bf MCS 571 S97 --- Numerical PDE --- PPDE Computer Problem --- Hanson}
\\*[1ex]
NUMERICAL METHODS FOR PPDEs\\*[1ex]
Problem program with output due Wednesday 05 March 1997 in class.
\end{center}
     Write a single, efficient program, in double precision (64 bit),
that compares the performance of the three numerical methods:
1) the Euler's Explicit Method (EEM),
2) the Crank-Nicolson Implicit Method (CNIM), and  3) the
Runge-Kutta Explicit Method (RKEM) modified for PPDEs.  Apply
these methods to the inhomogeneous, linear Parabolic Partial
Differential Equation:
{\large
$$
u_t = F(x,t,u,u_x,u_{xx}),~~ x_0 < x < x_m,~~ t > 0;
$$
\hfill(PPDE)
$$
F(x,t,u,u_x,u_{xx}) = a(x)\cdot u_{xx} + b(x)\cdot u_x + c(x)\cdot u + d(x,t).
$$}
Consider the special test example with spatially variable coefficients,
$$
a(x) = 5.712*x_m^2/(1 + x^2),~~~
b(x) = 0.571\cdot(x-x_0)\cdot(x_m-x)\cdot(1-0.648\cdot(x-x_0)/(x_m-x_0)),
$$
\hfill (COEFF)
$$
c(x) = 0.387\cdot(x-x0)\cdot(x_m-x)/((x_m-x_0)^2 + x^2), ~
d(x,t) = -0.067\cdot \exp(-0.571t)\cdot\sin(\PI\cdot x/(x_m-x_0)),
$$
on $x_0 < x < x_m$, with conditions,
$$
u(x,0) = 100\cdot \exp(-0.17\cdot (x-x_0)\cdot(x_m-x))
/\sqrt{\PI/0.571},~~~
x_0 < x < x_m;~~~~~~~~~~~~  (IC)
$$
$$
u(x_0,t) = 3.7t/(1+t^2),~~~
u(x_m,t) = 4.2\cdot\exp(-0.5t^2),~~~ t > 0; ~~~~~~~~~~~~~~~~~~  (BC)
$$
with $x_0 = 0.3$, $x_m = 6.1$, $m = 11$, $n= 101$ and mesh requirements,
$h = (x_m-x_0)/m$, $\widehat{r} = \widehat{a}\cdot k/h^2 = 0.01$, where 
$\widehat{a}=\min[a(x)]$,
$x_i = x_0+i\cdot h$ for $i=0$ to $m$, $t_j = j\cdot k$ for $j=0$ to $n$.
Use an $O(h^2)$ high order Forward Finite Difference approximation
for $u_x(x_0,t)$.
If your method needs corner values, i.e., u(x0,0) and u(xm,0), use the
average values of the boundary and initial conditions.
 
{\it Warning:  Do not write your code only for numerical values
given here, but write it more generally and not just for this test case.}
 
    {\bf Method Specifications:}
 
The {\it Euler's Explicit Method (EEM)} can be written as
$$
U_{i,j+1} = U_{ij} + k\cdot F(x_i,t_j,U_{ij},DU_{ij},DDU_{ij}), ~~~~~~ (EEM)
$$
where here $DU_{ij}$ and $DDU_{ij}$ are the
Central Finite Difference (CFD) approximations using $U_{ij}$ corresponding
to the derivatives $u_x$ and $u_{xx}$, respectively.
\par
The {\it Crank-Nicolson Implicit Method (CNIM)} is
$$
U_{i,j+1} = U_{ij} + k\cdot F(x_i,t_{j+\half},CU_{ij},DCU_{ij},DDCU_{ij}),
~~~~~~ (CNIM)
$$
where $CU = CU_{i,j} = U_{i,j+1/2} = (U_{i,j+1}+U_{ij})/2$ and the corresponding
CFD derivatives are $DCU_{ij}$ and $DDCU_{ij}$, all using $CU_{ij}$ as discrete
values.  Use the Thomas Tridiagonal Solver Algorithm to solve.
\par
The {\it 4th order Runge-Kutta Explicit method (RKEM)}
(modified PPDE for this Computer Problem from the 4th order RKEM ODE method)
is written,
$$
U_{i,j+1} = U_{ij} + k\cdot (RK1+2\cdot (RK2+RK3)+RK4)_{ij}/6, ~~~~~~ (RKEM)
$$
for interior points, where
$$
RK1_{ij} = F(x_i,t_j,RU1_{ij},DRU1_{ij},DDRU1_{ij}),
$$
with $RU1 = RU1_{ij} = U_{ij}$ with corresponding
derivatives $DRU1_{ij}$ and $DDRU1_{ij}$ (all indexed at $(i,j)$),
$$
RK2_{ij} = F(x_i,t_{j+\half},RU2_{ij},DRU2_{ij},DDRU2_{ij}),
$$
with $RU2 = RU2_{ij} = U_{ij}+\half\cdot k\cdot RK1_{ij}$ and corresponding
derivatives $DRU2_{ij}$ and $DDRU2_{ij}$,
$$
RK3_{ij} = F(x_i,t_{j+\half},RU3_{ij},DRU3_{ij},DDRU3_{ij}),
$$
with $RU3 = RU3_{ij} = U_{ij}+\half\cdot k\cdot RK2_{ij}$ and corresponding
derivatives $DRU3_{ij}$ and $DDRU3_{ij}$,
and
$$
RK4_{ij} = F(x_i,t_{j+1},RU4_{ij},DRU4_{ij},DDRU4_{ij}),
$$
with $RU4 = RU4_{ij} = U_{ij}+k\cdot RK3_{ij}$ and corresponding
derivatives $DRU4_{ij}$ and $DDRU4_{ij}$, using $RU4_{ij}$ as discrete values.
Each $DRU[p]_{ij}$ should satisfy derivative boundary conditions where they are
specified, for $[p] = 1$, 2, 3, and 4.
 
    {\bf Output Specifications:}
 
Tabulate for $j = 0$ to 5 step 1 and $j=11$ to $n$ step 10
(i.e., variable sample)
$$j ~~~~~~~~~  t_j$$
while for $i = 0$,2,5,9,m,
$$
i ~~~~~ EEM[U] ~~ CNIM[U] ~~ RKEM[U] ~~ RELDIF[EEM] ~~ RELDIF[CNIM]
$$
in 4 significant digit exponential format for the above floating point
headings, where $EEM[U]$, $CNIM[U]$ \& $RKEM[U]$ are the respective Euler's
Explicit, Crank-Nicolson Implicit \& Runge-Kutta Explicit Methods results for
$U_{ij}$, while $RELDIF[EEM]$ \& $RELDIF[CNIM]$ are the relative difference of
the first two methods compared to $RKEM$, i.e.,
$$RELDIF[EEM] = EEM[U]/RKEM[U]-1,$$
$$RELDIF[CNIM] = CNIM[U]/RKEM[U]-1,$$
provided $RKEM[U] \neq 0$.  Also print out the outputted locations,
$x_i$, at least once at the beginning of the tabulated numerical output.  In
addition, report the timings in seconds for each method, counting
common initializations of terms like coefficient arrays
for $\{a, b, c, d\}$ and output.
\par
     {\bf Some Proper Efficiency Suggestions:}
\par
1)  Compute once and store the spatially dependent coefficients
in an array to use in all time steps.
\par
2)  Avoid calculating items like constant expressions in loops, but keep
code general.
\par
Added Note:  These methods can also be used for nonlinear PPDEs,
directly for EEM and RKEM, but with a predictor-corrector modification
for CNIM.
\end{document}