$$ \newcommand{\uex}{{u_{\small\mbox{e}}}} \newcommand{\uexd}[1]{{u_{\small\mbox{e}, #1}}} \newcommand{\vex}{{v_{\small\mbox{e}}}} \newcommand{\vexd}[1]{{v_{\small\mbox{e}, #1}}} \newcommand{\Aex}{{A_{\small\mbox{e}}}} \newcommand{\half}{\frac{1}{2}} \newcommand{\halfi}{{1/2}} \newcommand{\tp}{\thinspace .} \newcommand{\Ddt}[1]{\frac{D #1}{dt}} \newcommand{\E}[1]{\hbox{E}\lbrack #1 \rbrack} \newcommand{\Var}[1]{\hbox{Var}\lbrack #1 \rbrack} \newcommand{\Std}[1]{\hbox{Std}\lbrack #1 \rbrack} \newcommand{\xpoint}{\boldsymbol{x}} \newcommand{\normalvec}{\boldsymbol{n}} \newcommand{\Oof}[1]{\mathcal{O}(#1)} \newcommand{\x}{\boldsymbol{x}} \newcommand{\X}{\boldsymbol{X}} \renewcommand{\u}{\boldsymbol{u}} \renewcommand{\v}{\boldsymbol{v}} \newcommand{\w}{\boldsymbol{w}} \newcommand{\V}{\boldsymbol{V}} \newcommand{\e}{\boldsymbol{e}} \newcommand{\f}{\boldsymbol{f}} \newcommand{\F}{\boldsymbol{F}} \newcommand{\stress}{\boldsymbol{\sigma}} \newcommand{\strain}{\boldsymbol{\varepsilon}} \newcommand{\stressc}{{\sigma}} \newcommand{\strainc}{{\varepsilon}} \newcommand{\I}{\boldsymbol{I}} \newcommand{\T}{\boldsymbol{T}} \newcommand{\dfc}{\alpha} % diffusion coefficient \newcommand{\ii}{\boldsymbol{i}} \newcommand{\jj}{\boldsymbol{j}} \newcommand{\kk}{\boldsymbol{k}} \newcommand{\ir}{\boldsymbol{i}_r} \newcommand{\ith}{\boldsymbol{i}_{\theta}} \newcommand{\iz}{\boldsymbol{i}_z} \newcommand{\Ix}{\mathcal{I}_x} \newcommand{\Iy}{\mathcal{I}_y} \newcommand{\Iz}{\mathcal{I}_z} \newcommand{\It}{\mathcal{I}_t} \newcommand{\If}{\mathcal{I}_s} % for FEM \newcommand{\Ifd}{{I_d}} % for FEM \newcommand{\Ifb}{{I_b}} % for FEM \newcommand{\setb}[1]{#1^0} % set begin \newcommand{\sete}[1]{#1^{-1}} % set end \newcommand{\setl}[1]{#1^-} \newcommand{\setr}[1]{#1^+} \newcommand{\seti}[1]{#1^i} \newcommand{\sequencei}[1]{\left\{ {#1}_i \right\}_{i\in\If}} \newcommand{\basphi}{\varphi} \newcommand{\baspsi}{\psi} \newcommand{\refphi}{\tilde\basphi} \newcommand{\psib}{\boldsymbol{\psi}} \newcommand{\sinL}[1]{\sin\left((#1+1)\pi\frac{x}{L}\right)} \newcommand{\xno}[1]{x_{#1}} \newcommand{\Xno}[1]{X_{(#1)}} \newcommand{\yno}[1]{y_{#1}} \newcommand{\Yno}[1]{Y_{(#1)}} \newcommand{\xdno}[1]{\boldsymbol{x}_{#1}} \newcommand{\dX}{\, \mathrm{d}X} \newcommand{\dx}{\, \mathrm{d}x} \newcommand{\ds}{\, \mathrm{d}s} \newcommand{\Real}{\mathbb{R}} \newcommand{\Integerp}{\mathbb{N}} \newcommand{\Integer}{\mathbb{Z}} $$ next

Finite difference methods for wave motion

Finite difference methods for wave motion

Hans Petter Langtangen [1, 2]

[1] Center for Biomedical Computing, Simula Research Laboratory
[2] Department of Informatics, University of Oslo

Dec 12, 2013

This is still a preliminary version.

Table of contents

Simulation of waves on a string
      Discretizing the domain
            Uniform meshes
      The discrete solution
      Fulfilling the equation at the mesh points
      Replacing derivatives by finite differences
            Algebraic version of the PDE
            Algebraic version of the initial conditions
      Formulating a recursive algorithm
      Sketch of an implementation
Verification
      A slightly generalized model problem
      Using an analytical solution of physical significance
      Manufactured solution
      Constructing an exact solution of the discrete equations
Implementation
      Making a solver function
      Verification: exact quadratic solution
      Visualization: animating the solution
            Visualization via SciTools
            Making movie files
            Skipping frames for animation speed
            Visualization via Matplotlib
      Running a case
      The benefits of scaling
Vectorization
      Operations on slices of arrays
      Finite difference schemes expressed as slices
      Verification
      Efficiency measurements
Exercises
      Exercise 1: Simulate a standing wave
            Remarks
      Exercise 2: Add storage of solution in a user action function
      Exercise 3: Use a class for the user action function
      Exercise 4: Compare several Courant numbers in one movie
      Project 5: Calculus with 1D mesh functions
Generalization: reflecting boundaries
      Neumann boundary condition
      Discretization of derivatives at the boundary
      Implementation of Neumann conditions
      Index set notation
      Alternative implementation via ghost cells
            Idea
            Implementation
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 model PDE 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 6: Find the analytical solution to a damped wave equation
      Problem 7: Explore symmetry boundary conditions
      Exercise 8: Send pulse waves through a layered medium
      Exercise 9: Compare discretizations of a Neumann condition
Analysis of the difference equations
      Properties of the solution of the wave equation
      More precise definition of Fourier representations
      Stability
            Preliminary results
            Numerical wave propagation
      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
            Discretizing the PDEs
            Handling boundary conditions where is \( u \) known
            Discretizing the Neumann condition
Implementation
      Scalar computations
            Domain and mesh
            Solution arrays
            Index sets
            Computing the solution
      Vectorized computations
      Verification
            Testing a quadratic solution
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
            Efficiency
Migrating loops to Fortran
      The Fortran subroutine
      Building the Fortran module with f2py
      How to avoid array copying
            Efficiency
Migrating loops to C via Cython
      Translating index pairs to single indices
      The complete C code
      The Cython interface file
      Building the extension module
            Efficiency
Migrating loops to C via f2py
      Migrating loops to C++ via f2py
Using classes to implement a simulator
Exercises
      Exercise 10: Check that a solution fulfills the discrete model
      Project 11: Calculus with 2D/3D mesh functions
      Exercise 12: Implement Neumann conditions in 2D
      Exercise 13: Test the efficiency of compiled loops in 3D
Applications of wave equations
      Waves on a string
            Damping
            External forcing
            Modeling the tension via springs
      Waves on a membrane
      Elastic waves in a rod
      The acoustic model for seismic waves
            Anisotropy
      Sound waves in liquids and gases
      Spherical waves
      The linear shallow water equations
            Wind drag on the surface
            Bottom drag
            Effect of the Earth's rotation
      Waves in blood vessels
      Electromagnetic waves
Exercises
      Exercise 14: Simulate waves on a non-homogeneous string
      Exercise 15: Simulate damped waves on a string
      Exercise 16: Simulate elastic waves in a rod
      Exercise 17: Simulate spherical waves
      Exercise 18: Explain why numerical noise occurs
      Exercise 19: Investigate harmonic averaging in a 1D model
      Problem 20: Implement open boundary conditions
      Problem 21: Earthquake-generated tsunami over a subsea hill
      Problem 22: Earthquake-generated tsunami over a 3D hill
      Problem 23: Investigate Matplotlib for visualization
      Problem 24: Investigate visualization packages
      Problem 25: Implement loops in compiled languages
      Exercise 26: Simulate seismic waves in 2D
      Project 27: Model 3D acoustic waves in a room
      Project 28: Solve a 1D transport equation
      Problem 29: General analytical solution of a 1D damped wave equation
      Problem 30: General analytical solution of a 2D damped wave equation

next