Graham W. Griffiths

Method of lines

Our physical world is most generally described in scientific and engineering terms with respect to threedimensional space and time which we abbreviate as spacetime. PDEs provide a mathematical description of physical spacetime, and they are therefore among the most widely used forms of mathematics. As a consequence, methods for the solution of PDEs, such as the MOL [Sch-91, Sch-09, Gri-11], are of broad interest in science and engineering.

ODE/PDE analysis of corneal curvature

The starting point for this paper is a nonlinear, two-point boundary value ordinary differential equation (BVODE) that defines corneal curvature according to a static force balance. A numerical solution to the BVODE is computed by first converting the BVODE to a parabolic partial differential equation (PDE) by adding an initial value (t, pseudo-time) derivative to the BVODE. A numerical solution to the PDE is then computed by the method of lines (MOL) with the calculation proceeding to a sufficiently large value of t such that the derivative in t reduces to essentially zero. The PDE solution at this point is also the solution for the BVODE. This procedure is implemented in R (an open source scientific programming system) and the programming is discussed in some detail. A series approximation to the solution is derived from which an estimate for the rate of convergence is obtained. This is compared to a fitted exponential model. Also, two linear approximations are derived, one of which leads to a closed form solution. Both provide solutions very close to that obtained from the full nonlinear model. An estimate for the cornea radius of curvature is also derived. The paper concludes with a discussion of the features of the solution to the ODE/PDE system.

Analysis of cornea curvature using radial basis functions - Part I

We discuss the solution of cornea curvature using a meshless method based on radial basis functions (RBFs). A full two-dimensional nonlinear thin membrane partial differential equation (PDE) model is introduced and solved using the multiquadratic (MQ) and inverse multiquadratic (IMQ) RBFs. This new approach does not rely on radial symmetry or other simplifying assumptions in respect of the cornea shape. It also provides an alternative to corneal topography modeling methods requiring accurate material parameter values, such as Youngs modulus and Poisson ratio, that may not be available. The results show good agreement with published corneal data and allow back calculations for estimating certain physical properties of the cornea, such as tension and elasticity coefficient. All calculations and generation of graphics were performed using the R language programming environment [34] and RStudio, the integrated development environment (IDE) for R [36], both of which are open source and free to download. Part II [48] of this paper demonstrates how the method has been used to provide a very accurate fit to a corneal measured data set.

Linear and nonlinear waves

The study of waves can be traced back to antiquity where philosophers, such as Pythagoras (c.560-480 BC), studied the relation of pitch and length of string in musical instruments. However, it was not until the work of Giovani Benedetti (1530-90), Isaac Beeckman (1588-1637) and Galileo (1564-1642) that the relationship between pitch and frequency was discovered. This started the science of acoustics, a term coined by Joseph Sauveur (1653-1716) who showed that strings can vibrate simultaneously at a fundamental frequency and at integral multiples that he called harmonics. Isaac Newton (1642-1727) was the first to calculate the speed of sound in his Principia. However, he assumed isothermal conditions so his value was too low compared with measured values. This discrepancy was resolved by Laplace (1749-1827) when he included adiabatic heating and cooling effects. The first analytical solution for a vibrating string was given by Brook Taylor (1685-1731). After this, advances were made by Daniel Bernoulli (1700-82), Leonard Euler (1707-83) and Jean d’Alembert (1717-83) who found the first solution to the linear wave equation, see section (3.2). Whilst others had shown that a wave can be represented as a sum of simple harmonic oscillations, it was Joseph Fourier (1768-1830) who conjectured that arbitrary functions can be represented by the superposition of an infinite sum of sines and cosines now known as the Fourier series. However, whilst his conjecture was controversial and not widely accepted at the time, Dirichlet subsequently provided a proof, in 1828, that all functions satisfying Dirichlet’s conditions (i.e. non-pathological piecewise continuous) could be represented by a convergent Fourier series. Finally, the subject of classical acoustics was laid down and presented as a coherent whole by John William Strutt (Lord Rayleigh, 1832-1901) in his treatise Theory of Sound. The science of modern acoustics has now moved into such diverse areas as sonar, auditoria, electronic amplifiers, etc.

Linear Advection Equation

The partial differential equation (PDE) analysis of convective systems is particularly challenging since convective (hyperbolic) PDEs can propagate steep fronts and even discontinuities. To demonstrate this characteristic, we consider in this chapter the numerical and analytical integration of the linear advection equation, possibly the simplest PDE, but ironically, one of the most diffcult to integrate numerically. The propagation of moving fronts is illustrated for several cases, from a smooth Gaussian pulse to a discontinuity; the latter is resolved with flux limiters. Matlab code for various MOL solutions are discussed in detail along with the associated numerical and graphical output. The concluding appendix provides a brief introduction to the theory of flux limiters and a survey of fifteen of these mathematical devices which have been used in a spectrum of convective system applications. All the associated Matlab code is available for download.

Kuramoto–Sivashinsky Equation

The Kuramoto-Sivashinsky equation has a nonlinear convection term and second, third and fourth order spatial (boundary value) derivatives. Thus, this example illustrates the solution of a higher (fourth) order PDE. The four required boundary conditions (BCs) are taken from the analytical solution as two Dirichlet BCs and two Neumann BCs. A numerical solution is computed by the method of lines (MOL), including detailed discussion of the Matlab routines and the numerical and graphical output. In the chapter appendix, diffusion-induced chaos produced by the Kuramoto-Sivashinsky equation is discussed whereby a simple initial condition can lead to chaotic behavior. Additionally, traveling wave solutions are derived using the tanh method, as outlined in the main Appendix. Maple code is presented which performs this procedure automatically to obtain the specific analytical solution used to validate the numerical solution. All the associated computer code is available for download, including additional Maple code that solves the PDE problem using the exp and Riccati methods.

Fitzhugh–Nagumo Equation

The Fitzhugh-Nagumo (F-N) partial differential equation (PDE) is an extension of the diffusion equation of Chapter 3 with a linear and a cubic source term. The BCs include a single pulse and a train of pulses in time. A numerical solution is computed by the method of lines (MOL), including detailed discussion of the Matlab routines and the numerical and graphical output. The appendix to this chapter analyses the spatially distributed, coupled Fitzhugh-Nagumo equations, specifically, how the application of an external stimulus to an axon can result in a traveling wave along an excitable medium (e.g., nerve). In addition, the tanh and exp methods, as outlined in the main Appendix, are used to obtain traveling wave solutions to the single equation form of the F-N equations. Maple code is presented which performs this procedure automatically to obtain the specific solution used to evaluate the numerical solution. All the associated computer code is available for download, including additional Maple code that solves the PDE problem using the Riccati method.

Analysis of cornea curvature using radial basis functions – Part II: Fitting to data-set

In part I we discussed the solution of corneal curvature using a 2D meshless method based on radial basis functions (RBFs). In Part II we use these methods to fit a full nonlinear thin membrane model to a measured data-set in order to generate a topological mathematical description of the cornea. In addition, we show how these results can lead to estimations for corneal radius of curvature and certain physical properties of the cornea; namely, tension and elasticity coefficient. Again all calculations and graphics generation were performed using the R language programming environment. The model describes corneal topology extremely well, and the estimated properties fall well within the expected range of values. The method is straight forward to implement and offers scope for further analysis using more detailed 3D models that include corneal thickness.

Method of lines solutions for the three-wave model of Brillouin equations

Purpose – The purpose of this paper is to present the method of lines (MOL) solution of the stimulated Brillouin scattering (SBS) equations (a system of three first-order hyperbolic partial differential equations (PDEs)), describing the three-wave interaction resulting from a coupling between light and acoustic waves. The system has complex numbers and boundary values. Design/methodology/approach – System of three first-order hyperbolic PDEs are first transformed and then spatially discretized. Superbee flux limiter is proposed to offset numerical damping and dispersion, brought on by the low order approximation of spatial derivatives in the PDEs. In order to increase computational efficiency, the structured structure of the PDE Jacobian matrix is identified and a sparse integration algorithm option of the ordinary differential equation (ODE) solvers is used. The flux limiter based on higher order approximations eliminates numerical oscillation. Examples are presented, and the performance of the Matlab OD...

18 – Boussinesq Equation

Publisher Summary This chapter discusses the Boussinesq equation that has: (1) a second order derivative in the initial value variable, (2) second-order derivatives in the spatial (boundary value) variable with respect to the dependent variable and the square of the dependent variable, and (3) a fourth order derivative in the spatial variable. The four required boundary conditions (BCs) are homogeneous in the dependent variable (Dirichlet BCs) and the second derivative of the dependent variable. The fourth-order spatial derivative is calculated by a finite difference (FD) derived specifically for fourth derivatives and by application of stage-wise differentiation in which a FD for second derivatives is used twice. The numerical solutions are computed by the method of lines (MOL), including detailed discussion of the Matlab routines and the numerical and graphical output. This chapter presents the background to the Boussinesq equation and conditions that must be satisfied for solitary wave solutions to exist. Solutions by application of direct integration and Riccati methods are obtained that match the solution used for verification of the numerical solutions. Maple code for the Riccati based solution is presented.