MCS 571 S97 --- Numerical Methods for PDEs --- Hanson

PPDE Computer Problem

Problem program with output due Wednesday 05 March 1997 in class. Write a single, efficient program, in double precision (64 bit), that compares the performance of the three numerical methods:

$u_t = F(x,t,u,u_x,u_{xx}),$

where $ x_0 < x < x_m $ , $ t > 0; $

$F(x,t,u,u_x,u_{xx}) = a(x) u_{xx} + b(x) u_x + c(x) 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 (x-x_0) (x_m-x) (1-0.648 (x-x_0)/(x_m-x_0)),$

$c(x) = 0.387 (x-x0) (x_m-x)/((x_m-x_0)^2 + x^2),$

$d(x,t) = -0.067 \exp(-0.571t) \sin(\pi x/(x_m-x_0)),$

on $ x_0 < x < x_m $ , with initial conditions: $u(x,0) = 100 \exp(-0.17 (x-x_0) (x_m-x)) /\sqrt{\pi/0.571} ,$

on $ x_0 < x < x_m; $ and boundary conditions:

$u(x_0,t) = 3.7t/(1+t^2),$

$u(x_m,t) = 4.2 \exp(-0.5t^2),$

for t > 0 with $ x_0 = 0.3 $ , $ x_m = 6.1 $ , m = 11, n= 101 and mesh requirements,

$h = (x_m-x_0)/m$ ,

$\widehat{r} \equiv \widehat{a} k/h^2 = 0.01$ ,

where $\widehat{a}=\min[a(x)]$ , $ x_i = x_0+i 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(x_0,0) $ and $ u(x_m,0) $ , use the average values of the boundary and initial conditions. Warning: Do not write your code only for numerical values given here, but write it more generally and not just for this test case.

Method Specifications: The Euler's Explicit Method (EEM) can be written as

$U_{i,j+1} = U_{ij} + k F(x_i,t_j,U_{ij},DU_{ij},DDU_{ij}),$

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.

The Crank-Nicolson Implicit Method (CNIM) is

$U_{i,j+1} = U_{ij} + k F(x_i,t_{j+\frac{1}{2}},CU_{ij},DCU_{ij},DDCU_{ij}),$

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.

The 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 (RK1+2 (RK2+RK3)+RK4)_{ij}/6,$

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+\frac{1}{2}},RU2_{ij},DRU2_{ij},DDRU2_{ij}),$

with $RU2 = RU2_{ij} = U_{ij}+\frac{1}{2} k RK1_{ij}$ and corresponding derivatives $DRU2_{ij}$ and $DDRU2_{ij}$ ,

$RK3_{ij} = F(x_i,t_{j+\frac{1}{2}},RU3_{ij},DRU3_{ij},DDRU3_{ij}),$

with $RU3 = RU3_{ij} = U_{ij}+\frac{1}{2} k 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 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.

Output Specifications: Tabulate for j = 0 to 5 step 1 and j=11 to n step 10 (i.e., variable sample) 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.

Some Proper Efficiency Suggestions:

1) Compute once and store the spatially dependent coefficients in an array to use in all time steps.

2) Avoid calculating items like constant expressions in loops, but keep code general.

Added Note: These methods can also be used for nonlinear PPDEs, directly for EEM and RKEM, but with a predictor-corrector modification for CNIM.

Web Source: http://www.math.uic.edu/~hanson/pcp571s97.html

Email Comments or Questions to hanson@uic.edu