Chapter 1
Introduction

The principal objective of this book is to present the method of lines (MOL) numerical integration of partial differential equations (PDEs), with spline collocation (SC) approximation of the PDE boundary-value derivatives. This approach is therefore termed spline collocation method of lines (SCMOL).

The details of SCMOL computer implementation are presented in terms of a series of applications, first for one-dimensional (1D) PDEs, then for two-dimensional (2D) PDEs, and finally for a series of legacy PDEs to illustrate the broad applicability of SCMOL. The approach is not with formal mathematics, for example, theorems and proofs, but rather, by examples of SCMOL discussed in detail, including routines in R.1

In this introduction, some basic properties of splines are reviewed, including example applications of the utilities (functions) for splines in R. The R spline utilities are then applied to PDEs in the following chapters within the SCMOL setting.

Splines are polynomials that can be used for the functional approximation of a set of numerical pairs

Table 1.1 Data pairs for spline interpolation

A cubic spline is of the form

The coefficients are evaluated to (i) return the original data of Table 1.1 and (ii) provide continuity in the computed first and second derivatives of at . The sequence in does not have to be uniformly spaced, which is a particularly useful feature in the SCMOL solution of PDEs, that is, the points can be placed as required to achieve good resolution in .

can be differentiated to give the first and second derivatives

The derivatives of Eqs. (1.1a)-(1.1d) can then be used to approximate first, second, and third derivatives in PDEs.

All of the operations reflected in Eqs. (1.1) are implemented in the R function splinefun,2 which returns the coefficients for a set of data as listed in Table 1.1. are computed by the solution of a tridiagonal algebraic system of equations. This set of coefficients therefore requires the specification of two additional conditions.3

A variety of additional conditions can be specified. For example, if the two second derivatives at the end points of the data set are set to zero, so-called natural cubic splines result. The details of the splines resulting from various sets of two additional conditions as implemented in splinefun are given in Appendix A1 and in [1].

The use of splinefun is illustrated with a series of examples that follows.

1.1 Uniform Grids

Application of splinefun to the function is illustrated with the following code:

Listing 1.1 Spline approximation of .

We can note the following details about this listing:

• Previous workspaces are cleared to avoid the unintended use of out-of-date files. # # Previous workspaces are cleared rm(list=ls(all=TRUE))
• A uniform grid is defined with 11 points, so that . # # Define uniform grid xl=0;xu=1;n=11; x=seq(from=xl,to=xu,by=(xu-xl)/(n-1)); The seq utility is used for this purpose.
• The function is defined on the grid in . # # Define function to be approximated u=sin(pi*x); This definition of u4 illustrates the vectorization available in R, that is, u is an -vector since is an -vector. Also, the function sin can operate on a vector to produce a vector.
• splinefun is used to define a table of spline coefficients, for example, in Eq. (1.1a). # # Set up spline table utable=splinefun(x,u); splinefun is part of the basic R and does not have to be accessed from an external library. Here the default spline mmf is used in which an exact cubic polynomial is fitted to the four points at each end of the data of Table 1.1 [1], p. 73.
• The table utable is used to compute spline approximations to at the values of in x. # # Compute spline approximation us=utable(x); Other values of within the interval xl to xu could be used for interpolation between the grid points.
• The values of u, their spline approximations us, and the difference between the two are displayed. # # Display comparison of function and its spline # approximation cat(sprintf("\n x u us diff")); for(i in 1:n){ cat(sprintf("\n%5.2f%10.5f%10.5f%12.7f", x[i],u[i],us[i],us[i]-u[i])); } The use of the combination cat(sprintf()) provides detailed formatting of the output.

Execution of the code in Listing 1.1 gives the following output:

Table 1.2 Output from Listing 1.1

This output demonstrates that the spline in splinefun returns the original data (an important feature of splines in general).

Equations (1.1a)-(1.1d) can also give numerical approximations to the first to third derivatives, as demonstrated by the following routine:

Listing 1.2 Spline approximation of and the first to third derivatives.

We can note the following details about Listing 1.2:

• Previous workspaces are cleared and the same uniform grid in as in Listing 1.1 is defined.
• The function and its first three derivatives are computed. Again, u, ux, uxx, uxxx are -vectors. # # Define function to be approximated, and its # derivatives u=sin(pi*x); ux=pi*cos(pi*x); uxx=-pi^2*sin(pi*x); uxxx=-pi^3*cos(pi*x); The derivatives are used to evaluate the numerical approximations from splinefun.
• The table of spline coefficients is defined as in Listing 1.1. # # Set up spline table for function, derivatives utable=splinefun(x,u);
• The function us and its first three derivatives are computed from utable. # # Compute spline approximation of function # derivatives us=utable(x); usx=utable(x,deriv=1); usxx=utable(x,deriv=2); usxxx=utable(x,deriv=3); us, usx, usxx, usxxx are -vectors. The argument deriv specifies the order of the derivative to be computed.
• The exact values of and the numerical approximations are compared and displayed. # # Display comparison of function and its spline # approximation, and its derivatives # # u cat(sprintf("\n x u us diff")); for(i in 1:n){ cat(sprintf("\n%5.2f%10.5f%10.5f%12.7f", x[i],u[i],us[i],us[i]-u[i])); } # # ux cat(sprintf("\n x ux usx diff")); for(i in 1:n){ ...