Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Define mass and stiffness matrix in the glossary. #14826

Merged
merged 3 commits into from
Feb 25, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
120 changes: 119 additions & 1 deletion doc/doxygen/headers/glossary.h
Original file line number Diff line number Diff line change
Expand Up @@ -77,7 +77,7 @@
\end{array}\right),
* @f}
* where $U,P$ are the values of velocity and pressure degrees of freedom,
* respectively, $M$ is the mass matrix on the velocity space, $B$ corresponds to
* respectively, $M$ is the @ref GlossMassMatrix "mass matrix" on the velocity space, $B$ corresponds to
* the negative divergence operator, and $B^T$ is its transpose and corresponds
* to the negative gradient.
*
Expand Down Expand Up @@ -1322,6 +1322,57 @@
* @see @ref manifold "The module on Manifolds"
*
*
* <dt class="glossary">@anchor GlossMassMatrix <b>Mass matrix</b></dt>
* <dd>The "mass matrix" is a matrix of the form
* @f{align*}{
* M_{ij} = \int_\Omega \varphi_i(\mathbf x) \varphi_j(\mathbf x)\; dx,
* @f}
* possibly with a coefficient inside the integral, and
* where $\varphi_i(\mathbf x)$ are the shape functions of a finite element.
* The origin of the term refers to the fact that in structural mechanics
* (where the finite element method originated), one often starts from the
* elastodynamics (wave) equation
* @f{align*}{
* \rho \frac{\partial^2 u}{\partial t^2}
* -\nabla \cdot C \nabla u = f.
* @f}
* If one multiplies this equation by a test function $\varphi_i$,
* integrates over $\Omega$, and then discretizes by the substitution
* $u(\mathbf x,t) \to u_h(\mathbf x)=\sum_j U_j(t) \varphi_j(\mathbf x)$,
* then the first term above results in
* @f{align*}{
* \sum_j \left[\int_\Omega \rho \varphi_i \varphi_j \right]
* \frac{\partial^2 U_j(t)}{\partial t^2}
* @f}
* which can be written as
* @f{align*}{
* M
* \frac{\partial^2 U(t)}{\partial t^2}
* @f}
* where
* @f{align*}{
* M_{ij} = \int_\Omega \rho(\mathbf x)\varphi_i(\mathbf x) \varphi_j(\mathbf x)\; dx.
* @f}
* Since the matrix entries are a (weighted) integral over a mass density, they
* have the units of "mass", giving the "mass matrix" its name.
*
* In mathematics, where we often consider non-dimensionalized equations, we
* end up with the case $\rho=1$, and as a consequence the matrix without
* the coefficient,
* @f{align*}{
* M_{ij} = \int_\Omega \varphi_i(\mathbf x) \varphi_j(\mathbf x)\; dx,
* @f}
* also carries the name "mass matrix".
*
* The mass matrix is almost always written with the symbol $M$. See, for example,
* step-23, step-26, and a number of the other time dependent equations solved by
* tutorial programs.
*
* See also the @ref GlossStiffnessMatrix "stiffness matrix"
* for a related case.
* </dt>
*
*
* <dt class="glossary">@anchor GlossMaterialId <b>Material id</b></dt>
* <dd>Each cell of a triangulation has associated with it a property called
* "material id". It is commonly used in problems with heterogeneous
Expand Down Expand Up @@ -1747,6 +1798,73 @@
* </dd>
*
*
* <dt class="glossary">@anchor GlossStiffnessMatrix <b>Stiffness matrix</b></dt>
* <dd>The "stiffness matrix" is a matrix of the form
* @f{align*}{
* A_{ij} = \int_\Omega \nabla\varphi_i(\mathbf x)
* \cdot \nabla\varphi_j(\mathbf x)\; dx,
* @f}
* possibly with a coefficient inside the integral, and
* where $\varphi_i(\mathbf x)$ are the shape functions of a finite element.
* The term is also used for variations of the case above, for example
* replacing the gradient by the symmetric gradient in the case where
* the solution variable is vector-valued (e.g., in elasticity, or the
* Stokes equations). The key feature is that in the integral, first
* derivatives are applied to both the test and trial functions,
* $\varphi_i,\varphi_j$.
*
* The origin of the term refers to the fact that in structural mechanics
* (where the finite element method originated), one often starts from the
* elastostatics equation
* @f{align*}{
* -\nabla \cdot C \nabla u = f.
* @f}
* In this equation, $C$ is the stress-strain tensor that, informally
* speaking, relates how much force one has to apply to obtain a
* unit displacement. In other words, it encodes the "stiffness" of
* the material: A large $C$, i.e., a large stiffness, means a large
* required force for a desired displacement and the other way around.
*
* If one multiplies this equation by a test function $\varphi_i$,
* integrates over $\Omega$, and then discretizes by the substitution
* $u(\mathbf x,t) \to u_h(\mathbf x)=\sum_j U_j(t) \varphi_j(\mathbf x)$,
* then after integration by parts one ends up with
* @f{align*}{
* \sum_j \left[\int_\Omega \nabla \varphi_i \cdot C \varphi_j \right]
* U_j
* @f}
* which can be written as
* @f{align*}{
* AU
* @f}
* where
* @f{align*}{
* A_{ij} = \int_\Omega \nabla\varphi_i(\mathbf x) \cdot C \nabla \varphi_j(\mathbf x)\; dx.
* @f}
* Since the matrix entries are (weighted) integrals of the stiffness
* coefficient, the resulting matrix is called the "stiffness matrix".
*
* In mathematics, where we often consider non-dimensionalized equations,
* we end up with the case $C=1$, and as a consequence the matrix without
* the coefficient,
* @f{align*}{
* A_{ij} = \int_\Omega \nabla\varphi_i(\mathbf x) \cdot \nabla\varphi_j(\mathbf x)\; dx
* @f}
* which corresponds to the Laplace or Poisson equation,
* @f{align*}{
* -\Delta u = f,
* @f}
* also carries the name "stiffness matrix".
*
* The stiffness matrix is almost always denotes by the symbol $A$. See, for example,
* step-4, step-6, as well as a number of the time dependent equations considered in
* programs such as step-23 or step-26.
*
* See also the @ref GlossStiffnessMatrix "stiffness matrix"
* for a related case.
* </dt>
*
*
* <dt class="glossary">@anchor GlossSubdomainId <b>Subdomain id</b></dt>
* <dd>Each cell of a triangulation has associated with it a property called
* the "subdomain id" that can be queried using a call like
Expand Down
2 changes: 1 addition & 1 deletion doc/doxygen/headers/vector_valued.h
Original file line number Diff line number Diff line change
Expand Up @@ -801,7 +801,7 @@
\begin{array}{cc} M & B \\ B^T & 0 \end{array}
\right),
@f}
* where $M$ represents the mass matrix that results from discretizing the
* where $M$ represents the @ref GlossMassMatrix "mass matrix" that results from discretizing the
* identity operator $\mathbf 1$ and $B$ the equivalent of the gradient
* operator.
*
Expand Down
2 changes: 1 addition & 1 deletion examples/doxygen/step_3_mixed.cc
Original file line number Diff line number Diff line change
Expand Up @@ -278,7 +278,7 @@ void Step3::setup_system()
// cell - a reference to the right FEValues object (constructed
// with the correct mapping object, finite element, and quadrature rule),
// which can be used as usual to compute the cell integrals.
// - The cell-local stiffness matrix and the right-hand-side vector have
// - The cell-local @ref GlossStiffnessMatrix "stiffness matrix" and the right-hand-side vector have
// different sizes depending on the cell type (6 DoFs vs. 9 DoFs) so that
// they might need to be resized for each cell.
//
Expand Down
2 changes: 1 addition & 1 deletion examples/step-11/step-11.cc
Original file line number Diff line number Diff line change
Expand Up @@ -266,7 +266,7 @@ namespace Step11
// That's quite simple, right?
//
// Two remarks are in order, though: First, these functions are used in a
// lot of contexts. Maybe you want to create a Laplace or mass matrix for
// lot of contexts. Maybe you want to create a Laplace or @ref GlossMassMatrix "mass matrix" for
// a vector values finite element; or you want to use the default Q1
// mapping; or you want to assembled the matrix with a coefficient in the
// Laplace operator. For this reason, there are quite a large number of
Expand Down
2 changes: 1 addition & 1 deletion examples/step-20/doc/intro.dox
Original file line number Diff line number Diff line change
Expand Up @@ -476,7 +476,7 @@ then we can subdivide the linear system $Ax=b$ into the following blocks:
\end{array}\right),
@f}
where $U,P$ are the values of velocity and pressure degrees of freedom,
respectively, $M$ is the mass matrix on the velocity space, $B^T$ corresponds to
respectively, $M$ is the @ref GlossMassMatrix "mass matrix" on the velocity space, $B^T$ corresponds to
the negative divergence operator, and $B$ is its transpose and corresponds
to the gradient.

Expand Down
2 changes: 1 addition & 1 deletion examples/step-20/step-20.cc
Original file line number Diff line number Diff line change
Expand Up @@ -634,7 +634,7 @@ namespace Step20
// After we have the pressure, we can compute the velocity. The equation
// reads $MU=-BP+F$, and we solve it by first computing the right hand
// side, and then multiplying it with the object that represents the
// inverse of the mass matrix:
// inverse of the @ref GlossMassMatrix "mass matrix":
U = op_M_inv * (F - op_B * P);
}

Expand Down
4 changes: 2 additions & 2 deletions examples/step-22/doc/intro.dox
Original file line number Diff line number Diff line change
Expand Up @@ -448,7 +448,7 @@ elaborate preconditioning scheme for the Schur complement.

<li>
Second, every time we multiplied with $B^TM^{-1}B$ we had to solve with the
mass matrix $M$. This wasn't particularly difficult, however, since the mass
@ref GlossMassMatrix "mass matrix" $M$. This wasn't particularly difficult, however, since the mass
matrix is always well conditioned and so simple to invert using CG and a
little bit of preconditioning.
</ol>
Expand Down Expand Up @@ -484,7 +484,7 @@ this Schur complement in a straightforward way with the CG method:
using no preconditioner, the condition number of the Schur complement matrix
depends on the size ratios of the largest to the smallest cells, and one still
needs on the order of 50-100 CG iterations. However, there is a simple cure:
precondition with the mass matrix on the pressure space and we get down to a
precondition with the @ref GlossMassMatrix "mass matrix" on the pressure space and we get down to a
number between 5-15 CG iterations, pretty much independently of the structure
of the mesh (take a look at the <a href="#Results">results section</a> of this
program to see that indeed the number of CG iterations does not change as we
Expand Down
2 changes: 1 addition & 1 deletion examples/step-22/doc/results.dox
Original file line number Diff line number Diff line change
Expand Up @@ -224,7 +224,7 @@ As explained during the generation of the sparsity pattern, it is
important to have the numbering of degrees of freedom in mind when
using preconditioners like incomplete LU decompositions. This is most
conveniently visualized using the distribution of nonzero elements in
the stiffness matrix.
the @ref GlossStiffnessMatrix "stiffness matrix".

If we don't do anything special to renumber degrees of freedom (i.e.,
without using DoFRenumbering::Cuthill_McKee, but with using
Expand Down
2 changes: 1 addition & 1 deletion examples/step-22/step-22.cc
Original file line number Diff line number Diff line change
Expand Up @@ -815,7 +815,7 @@ namespace Step22
SolverCG<Vector<double>> cg(solver_control);

// Now to the preconditioner to the Schur complement. As explained in
// the introduction, the preconditioning is done by a mass matrix in the
// the introduction, the preconditioning is done by a @ref GlossMassMatrix "mass matrix" in the
// pressure variable.
//
// Actually, the solver needs to have the preconditioner in the form
Expand Down
2 changes: 1 addition & 1 deletion examples/step-23/step-23.cc
Original file line number Diff line number Diff line change
Expand Up @@ -88,7 +88,7 @@ namespace Step23
// Next comes the declaration of the main class. It's public interface of
// functions is like in most of the other tutorial programs. Worth
// mentioning is that we now have to store four matrices instead of one: the
// mass matrix $M$, the Laplace matrix $A$, the matrix $M+k^2\theta^2A$ used
// @ref GlossMassMatrix "mass matrix" $M$, the Laplace matrix $A$, the matrix $M+k^2\theta^2A$ used
// for solving for $U^n$, and a copy of the mass matrix with boundary
// conditions applied used for solving for $V^n$. Note that it is a bit
// wasteful to have an additional copy of the mass matrix around. We will
Expand Down
2 changes: 1 addition & 1 deletion examples/step-24/step-24.cc
Original file line number Diff line number Diff line change
Expand Up @@ -95,7 +95,7 @@ namespace Step24
unsigned int timestep_number;
const double theta;

// Here's what's new: first, we need that boundary mass matrix $B$ that
// Here's what's new: first, we need that boundary @ref GlossMassMatrix "mass matrix" $B$ that
// came out of the absorbing boundary condition. Likewise, since this
// time we consider a realistic medium, we must have a measure of the
// wave speed $c_0$ that will enter all the formulas with the Laplace
Expand Down
2 changes: 1 addition & 1 deletion examples/step-25/doc/intro.dox
Original file line number Diff line number Diff line change
Expand Up @@ -166,7 +166,7 @@ Above, the matrix $F_h'(\cdot)$ and the vector $F_h(\cdot)$ denote the discrete
\f}
Again, note that the first matrix equation above is, in fact, the
definition of an iterative procedure, so it is solved multiple times
until a stopping criterion is met. Moreover, $M$ is the mass matrix,
until a stopping criterion is met. Moreover, $M$ is the @ref GlossMassMatrix "mass matrix",
i.e. $M_{ij} = \left( \varphi_i,\varphi_j \right)_{\Omega}$, $A$ is
the Laplace matrix, i.e. $A_{ij} = \left( \nabla \varphi_i, \nabla
\varphi_j \right)_{\Omega}$, $S$ is the nonlinear term in the
Expand Down
2 changes: 1 addition & 1 deletion examples/step-26/doc/intro.dox
Original file line number Diff line number Diff line change
Expand Up @@ -75,7 +75,7 @@ $k_n$:
\theta F^n
\right],
@f}
where $M$ is the mass matrix and $A$ is the stiffness matrix that results from
where $M$ is the @ref GlossMassMatrix "mass matrix" and $A$ is the @ref GlossStiffnessMatrix "stiffness matrix" that results from
discretizing the Laplacian. Bringing all known quantities to the right hand
side yields the linear system we have to solve in every step:
@f{align*}
Expand Down
2 changes: 1 addition & 1 deletion examples/step-26/doc/results.dox
Original file line number Diff line number Diff line change
Expand Up @@ -182,7 +182,7 @@ translated it to the following ones:
&
\theta k \left| a_{ij} \right| &\geq m_{ij},\qquad j \neq i,
@f}
where $M = m_{ij}$ denotes the mass matrix and $A = a_{ij}$ the stiffness
where $M = m_{ij}$ denotes the @ref GlossMassMatrix "mass matrix" and $A = a_{ij}$ the stiffness
matrix with $a_{ij} \leq 0$ for $j \neq i$, respectively. With
$a_{ij} \leq 0$, we can formulate bounds for the global time step $k$ as
follows:
Expand Down
2 changes: 1 addition & 1 deletion examples/step-28/doc/intro.dox
Original file line number Diff line number Diff line change
Expand Up @@ -469,7 +469,7 @@ these three cases, as follows:
= (B_c M_{K_c})^{ij} \phi_{g'}^j,
@f}
where $M_{K_c}^{lj}=\int_{K_c} f(x) \varphi_{g'}^l(x) \varphi_{g'}^j(x)$ is
the weighted mass matrix on child $c$ of cell $K$.
the weighted @ref GlossMassMatrix "mass matrix" on child $c$ of cell $K$.

The next question is what happens if a child $K_c$ of $K$ is not
active. Then, we have to apply the process recursively, i.e. we have to
Expand Down
2 changes: 1 addition & 1 deletion examples/step-28/step-28.cc
Original file line number Diff line number Diff line change
Expand Up @@ -842,7 +842,7 @@ namespace Step28
{
// The first case is that both cells are no further refined. In that case,
// we can assemble the relevant terms (see the introduction). This
// involves assembling the mass matrix on the finer of the two cells (in
// involves assembling the @ref GlossMassMatrix "mass matrix" on the finer of the two cells (in
// fact there are two mass matrices with different coefficients, one for
// the fission distribution cross section $\chi_g\nu\Sigma_{f,g'}$ and one
// for the scattering cross section $\Sigma_{s,g'\to g}$). This is
Expand Down
6 changes: 3 additions & 3 deletions examples/step-31/doc/intro.dox
Original file line number Diff line number Diff line change
Expand Up @@ -549,7 +549,7 @@ this yields for the simpler case of uniform time steps of size $k$:
On the left side of this equation remains the term from the time
derivative and the original (physical) diffusion which we treat
implicitly (this is actually a nice term: the matrices that result
from the left hand side are the mass matrix and a multiple of the
from the left hand side are the @ref GlossMassMatrix "mass matrix" and a multiple of the
Laplace matrix &mdash; both are positive definite and if the time step
size $k$ is small, the sum is simple to invert). On the right hand
side, the terms in the first line result from the time derivative; in
Expand Down Expand Up @@ -710,7 +710,7 @@ S^{-1}$ by a CG solver applied to the mass matrix on the pressure
space was a good choice. In a small deviation from step-22, we
here have a coefficient $\eta$ in the momentum equation, and by the same
derivation as there we should arrive at the conclusion that it is the weighted
mass matrix with entries $\tilde S_{ij}=(\eta^{-1}\varphi_i,\varphi_j)$ that
@ref GlossMassMatrix "mass matrix" with entries $\tilde S_{ij}=(\eta^{-1}\varphi_i,\varphi_j)$ that
we should be using.

It was more complicated to come up with a good replacement $\tilde
Expand Down Expand Up @@ -830,7 +830,7 @@ scheme and the current and previous time step. This being the sum of a
symmetric positive definite and a symmetric positive semidefinite matrix, the
result is also symmetric positive definite. Furthermore, $\frac\beta\alpha$ is
a number proportional to the time step, and so becomes small whenever the mesh
is fine, damping the effect of the then ill-conditioned stiffness matrix.
is fine, damping the effect of the then ill-conditioned @ref GlossStiffnessMatrix "stiffness matrix".

As a consequence, inverting this matrix with the Conjugate Gradient algorithm,
using a simple preconditioner, is trivial and very cheap compared to inverting
Expand Down
2 changes: 1 addition & 1 deletion examples/step-31/step-31.cc
Original file line number Diff line number Diff line change
Expand Up @@ -959,7 +959,7 @@ namespace Step31
}

// The creation of the temperature matrix (or, rather, matrices, since we
// provide a temperature mass matrix and a temperature stiffness matrix,
// provide a temperature mass matrix and a temperature @ref GlossStiffnessMatrix "stiffness matrix",
// that will be added together for time discretization) follows the
// generation of the Stokes matrix &ndash; except that it is much easier
// here since we do not need to take care of any blocks or coupling
Expand Down
2 changes: 1 addition & 1 deletion examples/step-32/doc/intro.dox
Original file line number Diff line number Diff line change
Expand Up @@ -323,7 +323,7 @@ system matrix $M$ in block form such that the product matrix
is of a form that Krylov-based iterative solvers like GMRES can solve in a
few iterations. We then replaced the exact inverse of $A$ by the action
of an AMG preconditioner $\tilde{A}$ based on a vector Laplace matrix,
approximated the Schur complement $S = B A^{-1} B^T$ by a mass matrix $M_p$
approximated the Schur complement $S = B A^{-1} B^T$ by a @ref GlossMassMatrix "mass matrix" $M_p$
on the pressure space and wrote an <tt>InverseMatrix</tt> class for
implementing the action of $M_p^{-1}\approx S^{-1}$ on vectors. In the
InverseMatrix class, we used a CG solve with an incomplete Cholesky (IC)
Expand Down
2 changes: 1 addition & 1 deletion examples/step-32/step-32.cc
Original file line number Diff line number Diff line change
Expand Up @@ -2457,7 +2457,7 @@ namespace Step32
// having inhomogeneous boundary conditions, by just making a right hand
// side at this point (compare the comments for the <code>project()</code>
// function above): We create some matrix columns with exactly the values
// that would be entered for the temperature stiffness matrix, in case we
// that would be entered for the temperature @ref GlossStiffnessMatrix "stiffness matrix", in case we
// have inhomogeneously constrained dofs. That will account for the correct
// balance of the right hand side vector with the matrix system of
// temperature.
Expand Down
4 changes: 2 additions & 2 deletions examples/step-36/doc/intro.dox
Original file line number Diff line number Diff line change
Expand Up @@ -83,8 +83,8 @@ continuous eigenvalue problem: @f[ \sum_j [(\nabla\phi_i,
\varepsilon_h \sum_j (\phi_i, \phi_j) \tilde{\psi}_j\quad. @f] In
matrix and vector notation, this equation then reads: @f[ A
\tilde{\Psi} = \varepsilon_h M \tilde{\Psi} \quad, @f] where $A$ is
the stiffness matrix arising from the differential operator $L$, and
$M$ is the mass matrix. The solution to the eigenvalue problem is an
the @ref GlossStiffnessMatrix "stiffness matrix" arising from the differential operator $L$, and
$M$ is the @ref GlossMassMatrix "mass matrix". The solution to the eigenvalue problem is an
eigenspectrum $\varepsilon_{h,\ell}$, with associated eigenfunctions
$\Psi_\ell=\sum_j \phi_j\tilde{\psi}_j$.

Expand Down
2 changes: 1 addition & 1 deletion examples/step-36/step-36.cc
Original file line number Diff line number Diff line change
Expand Up @@ -86,7 +86,7 @@ namespace Step36
DoFHandler<dim> dof_handler;

// With these exceptions: For our eigenvalue problem, we need both a
// stiffness matrix for the left hand side as well as a mass matrix for
// @ref GlossStiffnessMatrix "stiffness matrix" for the left hand side as well as a @ref GlossMassMatrix "mass matrix" for
// the right hand side. We also need not just one solution function, but a
// whole set of these for the eigenfunctions we want to compute, along
// with the corresponding eigenvalues:
Expand Down
2 changes: 1 addition & 1 deletion examples/step-37/step-37.cc
Original file line number Diff line number Diff line change
Expand Up @@ -191,7 +191,7 @@ namespace Step37
// parameter as is done in the MatrixFreeOperators::LaplaceOperator class.
//
// As a sidenote, if we implemented several different operations on the same
// grid and degrees of freedom (like a mass matrix and a Laplace matrix), we
// grid and degrees of freedom (like a @ref GlossMassMatrix "mass matrix" and a Laplace matrix), we
// would define two classes like the current one for each of the operators
// (derived from the MatrixFreeOperators::Base class), and let both of them
// refer to the same MatrixFree data cache from the general problem
Expand Down
2 changes: 1 addition & 1 deletion examples/step-41/doc/intro.dox
Original file line number Diff line number Diff line change
Expand Up @@ -185,7 +185,7 @@ get the equations
&[BU-G]_i \geq 0, \quad \Lambda_i \leq 0,\quad \Lambda_i[BU-G]_i = 0
\qquad \forall i.&
@f}
where $B$ is the mass matrix on the chosen finite element space and the
where $B$ is the @ref GlossMassMatrix "mass matrix" on the chosen finite element space and the
indices $i$ above are for all degrees of freedom in the set $\cal S$ of degrees of
freedom located in the interior of the domain
(we have Dirichlet conditions on the perimeter). However, we
Expand Down