Introduction to Finite Differences
and
Spring 2021
Finite Elements in the Geosciences
Instructors:
More Information:
Over the last
several decades, computing has played an
increasingly important role in the earth sciences. This trend has developed
along with the rapid increase in computational power
available to scientists, as well as the increasing
availability of numerical codes that can be used to
solve problems. However, to use these numerical
codes effectively, scientists must have a basic
understanding of how these codes work. This class
will introduce students to basic numerical
approaches to solving problems in the earth
sciences. We will focus primarily on finite
difference and finite element methods. Example
problems will be insprired by
the solid earth geopysics, but will be general enough
that the understanding gained from this class can
be used in broader fields of the earth sciences.
Required Texts:
Recommended Texts:
Class Meetings:
Homework: To catch up to speed on programming in Python, we
will read the first several chapters of S20 and/or LL20. Please tell
us which book you find more useful for this introduction.
Summary: We decided to proceed using LL20,
starting with Chapter 6. The earlier chapters of LL20 and the
entirety of S20 are useful references for implementing Python code,
and it would be useful to read them and implement some of the
examples. We do not plan to go through these chapters in detail,
unless we need to learn something specific from them.
Homework: For next time, read LL20 Chapter 6:
Computing Integrals and Testing Code. Carefully and write your own
Python codes to implement the examples in this chapter. For a
specific exercise, compute the integral of the sine function from 0
to pi (as in exercise 6.4a on page 169) using your own code to
implement the trapezoidal and midpoint methods. Make your own
version of the table on page 145, showing how the solution changes
as you increase n. For a bonus application, try implementing the
monte carlo method to solve this integral. Be prepared to share your
code, and to demonstrate how it works, with the group (e.g., using
share screen on Zoom).
Read Ahead: Next time, we will continue with LL20
Chapters 7 and 8. It would be a good idea to familiarize yourself
with these chapters in advance of our next meeting.
Summary: We talked about methods for numerical
integration (LL20 Chapter 6), using trapezoid rule and the midpoint method. We
discussed the monte carlo method, and generally talked about
programming techniques.
Homework: For next time, read LL20 Sections 8.1-8.3.
These sections give an introduction to writing code to solve
ordinary differential equations, with applications to exponential
problems. Please write code to implement the problems discussed in
these sections. Exercises 8.1-8.3 on pages 273-274 will help with
this.
Summary: We went through LL20, sections 8.1 to 8.3
and discussed some examples from the exercises assigned.
Homework: For next time, read sections 8.4 and 8.5 of LL20
(osciallating problems), and do exercises 8.9 and 8.11.
Summary: We went through LL20, sections 8.4 to
8.5. We reviewed the differences between the Forward Euler method and
more advanced finite difference methods such as Euler-Cromer (first
order) or Runge-Kutta (second or even fourth order) methods. We tested
our codes to see how the solution is affect by changing resolution.
Homework: For next time, read chapter 9 of LL20
(partial differential equations) and prepare a coding solution for exercise
9.2b.
Summary: We went through LL20, chapter 9. We
reviewed the diffusion equation, and how it is implemented in
finite-differences. We payed special attention to the boundary
conditions, including how they are implemented in the finite
difference solution. We discussed the temporal resolution (time step
size) that is necessary to obtain a stable solution, for a given spatial
resolution.
Homework: For next time, we begin the LL17 book,
which goes into more detail about finite difference methods. We start
with Chapter 1 (vibrational problems), sections 1.1 to 1.7. The book
goes into more detail than we need about some of the methods. Thus,
you can pass through the following subsections relatively quickly:
1.3.2-1.3.5, 1.4.6, 1.5.6, 1.74.
We will not assign an exercise this week, so we can focus on the
methods. However, it would be good to familiarize yourself with some
of the codes associated with this chapter, and try to get some of them
to work. Here is the website that provides some code for this book: [LL17
codes]
Note that the relevant code for this chapter can be found using the
"vib" link.
Summary: We went through LL17, chapter 1. We
discussed the stability or instability of finite difference solutions
to the vibrational ODE, how we can determine if a solution is stable,
and the conditions that make a solution unstable. We discussed the
different solution schemes (e.g., Centered Finite Difference, Forward
Euler, Backward Euler, Crank-Nicolson, Runge-Kutta, Forward-Backward
Euler-Cromer) and their implementation.
Homework: For next time, we continue with the LL17 book,
moving into solutions to the wave equation in Chapter 2. In
particular, we will discuss sections 2.1 to 2.4 in detail, but we will
omit sections 2.2.3, 2.2.4, 2.3.4, and 2.4.4. To practice
implementation, let's do exercises 2.1 and 2.4 (in section 2.5). Also,
let's read sections 2.6.1 to 2.6.3 (about reflecting boundary
conditions), section 2.7 (about variable wave velocity), and section
2.11 (about FD methods in 2D and 3D).
Summary: We went through LL17, chapter 2. We
discussed using a stencil to implement a finite difference solution to
the wave equation in 1 dimension (and also in 2 dimensions). The
discussed the importance of the Courant condition, which specifies the
relationship between spatial and temporal resolutions, and leads to an
unstable solution if it is violated. We implented standing waves using
a finite difference scheme, and discussed implementation of boundary conditions.
Homework: For next time, we will finish our reading
of the LL17 book, looking at the diffusion equation (Chapter 3) and
advection (Chapter 4). For Chapter 3, let's read sections 3.1 and 3.2,
as well as 3.5.1 to 3.3.5 (heterogeneous media) and 3.6.1-3.6.2 (two
dimensions). For Chapter 4, we will read sections 4.1.1 to 4.1.4.
Summary: We went through LL17, sections 3.1-3.2,
3.5, 3.6, and 4.1. In doing so we discussed the diffusion equation and
in particular the difference between the explicity solution (which we
used before for the wave equation) and the implicity solution. The
implicity solution is useful for the diffusion equation because it
allows for longer timesteps for a given spatial resolution (larger F
parameter). However, the implicit solution requires inversion of a
generally sparse or banded matrix, and there are efficient ways to do
this. The advection equation (section 4.1) is most readily solved with
a simple upwind differences scheme.
Homework: For next time, we start with finite
elements using the LM19 book. Let's read chapter 1, and chapter 2
sections 2.1-2.4, as well as a quick look at section 2.5. We have a
three-part exercise: (a) Write your own numerical integration program
(as we did early in the semester) to reproduce the example in section
2.2.3 and plot the result (you should get something similar to figure
2.2). (b) Now apply this same numerical integration method to
approximate the function y=exp(x^2) using a linear function in the
span [0,2]. (c) Finally, find the analytical expressions for the 3
Lagrange polynomials of order 2 for an element with points at [0.0,
0.5, 1.0], and plot them.
Summary: We went through chapters 1 and 2 of LM19,
and sections 2.1-2.4 in particular. We discussed how to approximate
vectors using basis vectors (easier for orthogonal vectors)
and we expanded on this idea to include approximation of functions
using (preferably orthogonal) basis functions. We examined the example
of the least squares method, which is related to the Taylor series expansion of polynomial
functions. Better sets of basis functions (because they have a higher
degree or orthogonality) are the sin and cos functions (e.g., Fourier
series), as well as the Lagrange polynomials (good for interpolation
problems) and Bernstein polynomials. We set up the formalism for
expressing any given function as a sum of a set of basis
functions, where each basis function is multiplied by a given
coefficient. We discussed the method for computing these coefficients,
and used the exercise in section 2.2.3 as an example.
Homework: For next time, we will continue with the LM19
book, moving into finite element basis functions in chapter 3. Let's
read sections 3.1, 3.2, and 3.5, and look quickly at section
3.6. Let's do exercises 3.1 and 3.3 (which can be done on paper), and
exercise 3.5 (it is sufficient to just map out what should be programmed) in section 3.8.
Summary: This week we spend most of our time
discussing section 3.1 of LM19. This discussion involved defining the
finite element basis functions, which are simply linear functions for
first-order functions (P1), and quadratics for 2nd order (P2). These
functions are defined only locally for each element (each element has 2
or more nodes associated with it). We how we can take integrals of
combinations of these P1 or P2 functions to obtain "local" versions of
the A matrix and b vector. These local versions can then be assembled
into the "global" A matrix and b vector. We looked at implementation
(section 3.2), and extensions into 2 and 3 dimensions (3.6)
Homework: For next time, we will continue with
LM19, reading sections 4.1 (approximating differential equations), 5.1
(computing with finite elements) and 6.1 (time-dependence). We will
not assign an exercise this week - our goal is to understand the
reading.
Summary: We fininshed our reading of LM19, and
discussed sections 4.1, 5.1, and 6.1, extending our analysis of basis
functions to the solutions of differential equations. We first
discussed how we can use basis functions generally to approximate
solutions to differential equations (sections 4.1), including
integration by parts to convert the differential equation into a "weak form"
that allows approximation of the solution using basis functions are not
second order differentiable. Next, we applied this understanding to
the finite element basis function (section 5.1) and discussed how to
construct a matrix system for P1 elements. Finally, we discussed how
to set up finite differences in time to solve time-dependent problems
with finite element solutions in space (section 6.1). This involves
constructing both the mass matrix M and the stiffness matrix K. Both
of these usually only need to be constructed once because they are
dependent on the finite element mesh and physical properties of the
material (e.g., viscosity or a diffusivity constant), which do not
change with time for most problems.
Homework: We will meet next week to discuss the
projects, so it would be helpful to get started on them this week so we
can answer questions. We also discussed how it would be useful to do an
exercise about constructing the mass and stiffness matrices. Here is
an exercise:
Consider the linear system described in equation (6.16), which solves
for time-dependent diffusion. Let's assume a diffusion coefficient of
1 mm^2/s and a domain width of 5 mm. First, break the domain into 5
equal-sized elements (with nodes at 0.0, 1.0, 2.0, 3.0, 4.0, and 5.0) and
construct M and K. Next, let's assume that we need higher resolution
on the left side, so we shift the nodes slightly (now nodes are at 0.0,
0.5, 1.2, 2.2, 3.5, and 5.0). Construct the M and K matrices for this
non-uniform mesh.
Project writeups are due by 13.00 on Friday, May 28.