Preface
Vibration ODEs
Finite difference discretization
A basic model for vibrations
A centered finite difference scheme
Implementation
Making a solver function
Verification
Scaled model
Visualization of long time simulations
Using a moving plot window
Making animations
Using Bokeh to compare graphs
Using a line-by-line ascii plotter
Empirical analysis of the solution
Analysis of the numerical scheme
Deriving a solution of the numerical scheme
The error in the numerical frequency
Empirical convergence rates and adjusted \( \omega \)
Exact discrete solution
Convergence
The global error
Stability
About the accuracy at the stability limit
Alternative schemes based on 1st-order equations
The Forward Euler scheme
The Backward Euler scheme
The Crank-Nicolson scheme
Comparison of schemes
Runge-Kutta methods
Analysis of the Forward Euler scheme
Energy considerations
Derivation of the energy expression
An error measure based on energy
The Euler-Cromer method
Forward-backward discretization
Equivalence with the scheme for the second-order ODE
Implementation
The Stoermer-Verlet algorithm
Staggered mesh
The Euler-Cromer scheme on a staggered mesh
Implementation of the scheme on a staggered mesh
Exercises and Problems
Problem 1.1: Use linear/quadratic functions for verification
Exercise 1.2: Show linear growth of the phase with time
Exercise 1.3: Improve the accuracy by adjusting the frequency
Exercise 1.4: See if adaptive methods improve the phase error
Exercise 1.5: Use a Taylor polynomial to compute \( u^1 \)
Problem 1.6: Derive and investigate the velocity Verlet method
Problem 1.7: Find the minimal resolution of an oscillatory function
Exercise 1.8: Visualize the accuracy of finite differences for a cosine function
Exercise 1.9: Verify convergence rates of the error in energy
Exercise 1.10: Use linear/quadratic functions for verification
Exercise 1.11: Use an exact discrete solution for verification
Exercise 1.12: Use analytical solution for convergence rate tests
Exercise 1.13: Investigate the amplitude errors of many solvers
Problem 1.14: Minimize memory usage of a simple vibration solver
Problem 1.15: Minimize memory usage of a general vibration solver
Exercise 1.16: Implement the Euler-Cromer scheme for the generalized model
Problem 1.17: Interpret \( [D_tD_t u]^n \) as a forward-backward difference
Exercise 1.18: Analysis of the Euler-Cromer scheme
Generalization: damping, nonlinearities, and excitation
A centered scheme for linear damping
A centered scheme for quadratic damping
A forward-backward discretization of the quadratic damping term
Implementation
Verification
Visualization
User interface
The Euler-Cromer scheme for the generalized model
The Stoermer-Verlet algorithm for the generalized model
A staggered Euler-Cromer scheme for a generalized model
The PEFRL 4th-order accurate algorithm
Exercises and Problems
Exercise 1.19: Implement the solver via classes
Problem 1.20: Use a backward difference for the damping term
Exercise 1.21: Use the forward-backward scheme with quadratic damping
Applications of vibration models
Oscillating mass attached to a spring
General mechanical vibrating system
A sliding mass attached to a spring
A jumping washing machine
Motion of a pendulum
Dynamic free body diagram during pendulum motion
Motion of an elastic pendulum
Vehicle on a bumpy road
Bouncing ball
Two-body gravitational problem
Electric circuits
Exercises
Exercise 1.22: Simulate resonance
Exercise 1.23: Simulate oscillations of a sliding box
Exercise 1.24: Simulate a bouncing ball
Exercise 1.25: Simulate a simple pendulum
Exercise 1.26: Simulate an elastic pendulum
Exercise 1.27: Simulate an elastic pendulum with air resistance
Exercise 1.28: Implement the PEFRL algorithm
Wave equations
Simulation of waves on a string
Discretizing the domain
The discrete solution
Fulfilling the equation at the mesh points
Replacing derivatives by finite differences
Formulating a recursive algorithm
Sketch of an implementation
Verification
A slightly generalized model problem
Using an analytical solution of physical significance
Manufactured solution and estimation of convergence rates
Constructing an exact solution of the discrete equations
Implementation
Callback function for user-specific actions
The solver function
Verification: exact quadratic solution
Verification: convergence rates
Visualization: animating the solution
Running a case
Working with a scaled PDE model
Vectorization
Operations on slices of arrays
Finite difference schemes expressed as slices
Verification
Efficiency measurements
Remark on the updating of arrays
Exercises
Exercise 2.1: Simulate a standing wave
Exercise 2.2: Add storage of solution in a user action function
Exercise 2.3: Use a class for the user action function
Exercise 2.4: Compare several Courant numbers in one movie
Exercise 2.5: Implementing the solver function as a generator
Project 2.6: Calculus with 1D mesh functions
Generalization: reflecting boundaries
Neumann boundary condition
Discretization of derivatives at the boundary
Implementation of Neumann conditions
Index set notation
Verifying the implementation of Neumann conditions
Alternative implementation via ghost cells
Generalization: variable wave velocity
The model PDE with a variable coefficient
Discretizing the variable coefficient
Computing the coefficient between mesh points
How a variable coefficient affects the stability
Neumann condition and a variable coefficient
Implementation of variable coefficients
A more general PDE model with variable coefficients
Generalization: damping
Building a general 1D wave equation solver
User action function as a class
Pulse propagation in two media
Exercises
Exercise 2.7: Find the analytical solution to a damped wave equation
Problem 2.8: Explore symmetry boundary conditions
Exercise 2.9: Send pulse waves through a layered medium
Exercise 2.10: Explain why numerical noise occurs
Exercise 2.11: Investigate harmonic averaging in a 1D model
Problem 2.12: Implement open boundary conditions
Exercise 2.13: Implement periodic boundary conditions
Exercise 2.14: Compare discretizations of a Neumann condition
Exercise 2.15: Verification by a cubic polynomial in space
Analysis of the difference equations
Properties of the solution of the wave equation
More precise definition of Fourier representations
Stability
Numerical dispersion relation
Extending the analysis to 2D and 3D
Finite difference methods for 2D and 3D wave equations
Multi-dimensional wave equations
Mesh
Discretization
Implementation
Scalar computations
Vectorized computations
Verification
Visualization
Exercises
Exercise 2.16: Check that a solution fulfills the discrete model
Project 2.17: Calculus with 2D mesh functions
Exercise 2.18: Implement Neumann conditions in 2D
Exercise 2.19: Test the efficiency of compiled loops in 3D
Applications of wave equations
Waves on a string
Elastic waves in a rod
Waves on a membrane
The acoustic model for seismic waves
Sound waves in liquids and gases
Spherical waves
The linear shallow water equations
Waves in blood vessels
Electromagnetic waves
Exercises
Exercise 2.20: Simulate waves on a non-homogeneous string
Exercise 2.21: Simulate damped waves on a string
Exercise 2.22: Simulate elastic waves in a rod
Exercise 2.23: Simulate spherical waves
Problem 2.24: Earthquake-generated tsunami over a subsea hill
Problem 2.25: Earthquake-generated tsunami over a 3D hill
Problem 2.26: Investigate Mayavi for visualization
Problem 2.27: Investigate visualization packages
Problem 2.28: Implement loops in compiled languages
Exercise 2.29: Simulate seismic waves in 2D
Project 2.30: Model 3D acoustic waves in a room
Project 2.31: Solve a 1D transport equation
Problem 2.32: General analytical solution of a 1D damped wave equation
Problem 2.33: General analytical solution of a 2D damped wave equation
Diffusion equations
An explicit method for the 1D diffusion equation
The initial-boundary value problem for 1D diffusion
Forward Euler scheme
Implementation
Verification
Numerical experiments
Implicit methods for the 1D diffusion equation
Backward Euler scheme
Sparse matrix implementation
Crank-Nicolson scheme
The unifying \( \theta \) rule
Experiments
The Laplace and Poisson equation
Analysis of schemes for the diffusion equation
Properties of the solution
Analysis of discrete equations
Analysis of the finite difference schemes
Analysis of the Forward Euler scheme
Analysis of the Backward Euler scheme
Analysis of the Crank-Nicolson scheme
Analysis of the Leapfrog scheme
Summary of accuracy of amplification factors
Analysis of the 2D diffusion equation
Explanation of numerical artifacts
Exercises
Exercise 3.1: Explore symmetry in a 1D problem
Exercise 3.2: Investigate approximation errors from a \( u_x=0 \) boundary condition
Exercise 3.3: Experiment with open boundary conditions in 1D
Exercise 3.4: Simulate a diffused Gaussian peak in 2D/3D
Exercise 3.5: Examine stability of a diffusion model with a source term
Diffusion in heterogeneous media
Discretization
Implementation
Stationary solution
Piecewise constant medium
Implementation of diffusion in a piecewise constant medium
Axi-symmetric diffusion
Spherically-symmetric diffusion
Diffusion in 2D
Discretization
Numbering of mesh points versus equations and unknowns
Algorithm for setting up the coefficient matrix
Implementation with a dense coefficient matrix
Verification: exact numerical solution
Verification: convergence rates
Implementation with a sparse coefficient matrix
The Jacobi iterative method
Implementation of the Jacobi method
Test problem: diffusion of a sine hill
The relaxed Jacobi method and its relation to the Forward Euler method
The Gauss-Seidel and SOR methods
Scalar implementation of the SOR method
Vectorized implementation of the SOR method
Direct versus iterative methods
The Conjugate gradient method
What is the recommended method for solving linear systems?
Random walk
Random walk in 1D
Statistical considerations
Playing around with some code
Equivalence with diffusion
Implementation of multiple walks
Demonstration of multiple walks
Ascii visualization of 1D random walk
Random walk as a stochastic equation
Random walk in 2D
Random walk in any number of space dimensions
Multiple random walks in any number of space dimensions
Applications
Diffusion of a substance
Heat conduction
Porous media flow
Potential fluid flow
Streamlines for 2D fluid flow
The potential of an electric field
Development of flow between two flat plates
Flow in a straight tube
Tribology: thin film fluid flow
Propagation of electrical signals in the brain
Exercises
Exercise 3.6: Stabilizing the Crank-Nicolson method by Rannacher time stepping
Project 3.7: Energy estimates for diffusion problems
Exercise 3.8: Splitting methods and preconditioning
Problem 3.9: Oscillating surface temperature of the earth
Problem 3.10: Oscillating and pulsating flow in tubes
Problem 3.11: Scaling a welding problem
Exercise 3.12: Implement a Forward Euler scheme for axi-symmetric diffusion
Advection-dominated equations
One-dimensional time-dependent advection equations
Simplest scheme: forward in time, centered in space
Analysis of the scheme
Leapfrog in time, centered differences in space
Upwind differences in space
Periodic boundary conditions
Implementation
A Crank-Nicolson discretization in time and centered differences in space
The Lax-Wendroff method
Analysis of dispersion relations
One-dimensional stationary advection-diffusion equation
A simple model problem
A centered finite difference scheme
Remedy: upwind finite difference scheme
Time-dependent convection-diffusion equations
Forward in time, centered in space scheme
Forward in time, upwind in space scheme
Two-dimensional advection-diffusion equations
Applications of advection equations
Transport of a substance
Transport of a heat
Exercises
Exercise 4.1: Analyze 1D stationary convection-diffusion problem
Exercise 4.2: Interpret upwind difference as artificial diffusion
Nonlinear problems
Introduction of basic concepts
Linear versus nonlinear equations
A simple model problem
Linearization by explicit time discretization
Exact solution of nonlinear algebraic equations
Linearization
Picard iteration
Linearization by a geometric mean
Newton's method
Relaxation
Implementation and experiments
Generalization to a general nonlinear ODE
Systems of ODEs
Systems of nonlinear algebraic equations
Picard iteration
Newton's method
Stopping criteria
Example: A nonlinear ODE model from epidemiology
Linearization at the differential equation level
Explicit time integration
Backward Euler scheme and Picard iteration
Backward Euler scheme and Newton's method
Crank-Nicolson discretization
1D stationary nonlinear differential equations
Finite difference discretization
Solution of algebraic equations
Multi-dimensional nonlinear PDE problems
Finite difference discretization
Continuation methods
Operator splitting methods
Ordinary operator splitting for ODEs
Strang splitting for ODEs
Example: Logistic growth
Reaction-diffusion equation
Example: Reaction-Diffusion with linear reaction term
Analysis of the splitting method
Exercises
Problem 5.1: Determine if equations are nonlinear or not
Problem 5.2: Derive and investigate a generalized logistic model
Problem 5.3: Experience the behavior of Newton's method
Exercise 5.4: Compute the Jacobian of a \( 2\times 2 \) system
Problem 5.5: Solve nonlinear equations arising from a vibration ODE
Exercise 5.6: Find the truncation error of arithmetic mean of products
Problem 5.7: Newton's method for linear problems
Problem 5.8: Discretize a 1D problem with a nonlinear coefficient
Problem 5.9: Linearize a 1D problem with a nonlinear coefficient
Problem 5.10: Finite differences for the 1D Bratu problem
Problem 5.11: Discretize a nonlinear 1D heat conduction PDE by finite differences
Problem 5.12: Differentiate a highly nonlinear term
Exercise 5.13: Crank-Nicolson for a nonlinear 3D diffusion equation
Problem 5.14: Find the sparsity of the Jacobian
Problem 5.15: Investigate a 1D problem with a continuation method
Appendix: Useful formulas
Finite difference operator notation
Truncation errors of finite difference approximations
Finite differences of exponential functions
Finite differences of \( t^n \)
Software
Appendix: Truncation error analysis
Overview of truncation error analysis
Abstract problem setting
Error measures
Truncation errors in finite difference formulas
Example: The backward difference for \( u'(t) \)
Example: The forward difference for \( u'(t) \)
Example: The central difference for \( u'(t) \)
Overview of leading-order error terms in finite difference formulas
Software for computing truncation errors
Exponential decay ODEs
Forward Euler scheme
Crank-Nicolson scheme
The \( \theta \)-rule
Using symbolic software
Empirical verification of the truncation error
Increasing the accuracy by adding correction terms
Extension to variable coefficients
Exact solutions of the finite difference equations
Computing truncation errors in nonlinear problems
Vibration ODEs
Linear model without damping
Model with damping and nonlinearity
Extension to quadratic damping
The general model formulated as first-order ODEs
Wave equations
Linear wave equation in 1D
Finding correction terms
Extension to variable coefficients
1D wave equation on a staggered mesh
Linear wave equation in 2D/3D
Diffusion equations
Linear diffusion equation in 1D
Nonlinear diffusion equation in 1D
Exercises
Exercise B.1: Truncation error of a weighted mean
Exercise B.2: Simulate the error of a weighted mean
Exercise B.3: Verify a truncation error formula
Problem B.4: Truncation error of the Backward Euler scheme
Exercise B.5: Empirical estimation of truncation errors
Exercise B.6: Correction term for a Backward Euler scheme
Problem B.7: Verify the effect of correction terms
Problem B.8: Truncation error of the Crank-Nicolson scheme
Problem B.9: Truncation error of \( u'=f(u,t) \)
Exercise B.10: Truncation error of \( [D_t D_tu]^n \)
Exercise B.11: Investigate the impact of approximating \( u'(0) \)
Problem B.12: Investigate the accuracy of a simplified scheme
Appendix: Software engineering; wave equation model
A 1D wave equation simulator
Mathematical model
Numerical discretization
A solver function
Saving large arrays in files
Using savez to store arrays in files
Using joblib to store arrays in files
Using a hash to create a file or directory name
Software for the 1D wave equation
Making hash strings from input data
Avoiding rerunning previously run cases
Verification
Programming the solver with classes
Class Problem
Class Mesh
Class Function
Class Solver
Migrating loops to Cython
Declaring variables and annotating the code
Visual inspection of the C translation
Building the extension module
Calling the Cython function from Python
Migrating loops to Fortran
The Fortran subroutine
Building the Fortran module with f2py
How to avoid array copying
Migrating loops to C via Cython
Translating index pairs to single indices
The complete C code
The Cython interface file
Building the extension module
Migrating loops to C via f2py
Migrating loops to C++ via f2py
Exercises
Exercise C.1: Explore computational efficiency of numpy.sum versus built-in sum
Exercise C.2: Make an improved numpy.savez function
Exercise C.3: Visualize the impact of the Courant number
Exercise C.4: Visualize the impact of the resolution
References