Merge lp:~hnarayanan/fenics-book/proofreading into lp:fenics-book

[Harish Narayanan]

[lp:~hnarayanan/fenics-book/proofreading]

639. By Harish Narayanan on 2011-05-13

merge

640. By Harish Narayanan on 2011-05-13

Some more fixes for [lopes]

641. By Harish Narayanan on 2011-05-13

Final fixes for [lopes]

642. By Harish Narayanan on 2011-05-13

Standardised backmatter a little bit and made it look more like the rest of the chapters

643. By Harish Narayanan on 2011-05-13

Updated STATUS

[Anders Logg (logg)]

Thanks!

--
Anders

On Mon, May 16, 2011 at 07:33:34AM -0000, Harish Narayanan wrote:
> You have been requested to review the proposed merge of lp:~hnarayanan/fenics-book/proofreading into lp:fenics-book.
>
> For more details, see:
> https://code.launchpad.net/~hnarayanan/fenics-book/proofreading/+merge/60776
>
>
>

> === modified file 'STATUS'
> --- STATUS 2011-05-13 08:16:02 +0000
> +++ STATUS 2011-05-13 10:11:49 +0000
> @@ -95,7 +95,7 @@
> Nuno D. Lopes <email address hidden> 6 [kirby-5], APPLIED
> Kent-Andre Mardal <email address hidden> 13, 14 [kirby-3], [kirby-2] SUBMITTED (AL will fix), SUBMITTED (RCK will fix)
> Mikael Mortensen <email address hidden> 32 [nikbakht]
> -Harish Narayanan <email address hidden> 8, 27 [oelgaard-2], [lopes] SUBMITTED, WAITING
> +Harish Narayanan <email address hidden> 8, 27 [oelgaard-2], [lopes] APPLIED, SUBMITTED
> Murtazo Nazarov <email address hidden> 21 [terrel]
> Mehdi Nikbakht <email address hidden> 36 [lezar] APPLIED
> Johannes Ring <email address hidden> 7, 33 [logg-3], [schroll] APPLIED, APPLIED
>
> === modified file 'backmatter/authors.tex'
> --- backmatter/authors.tex 2011-04-29 11:21:33 +0000
> +++ backmatter/authors.tex 2011-05-13 10:11:49 +0000
> @@ -1,5 +1,7 @@
> \chapter*{List of authors}
> \addcontentsline{toc}{chapter}{List of authors}
> +\addtocounter{chapter}{1}
> +\vspace{-1cm}
>
> \emph{The following authors have contributed to this book.}
>
>
> === modified file 'backmatter/fdl-1.3.tex'
> --- backmatter/fdl-1.3.tex 2011-05-10 18:26:12 +0000
> +++ backmatter/fdl-1.3.tex 2011-05-13 10:11:49 +0000
> @@ -1,25 +1,20 @@
> -\chapter*{\rlap{GNU Free Documentation License}}
> -\phantomsection % so hyperref creates bookmarks
> +\chapter*{GNU Free Documentation License}
> \addcontentsline{toc}{chapter}{GNU Free Documentation License}
> -%\label{label_fdl}
> -
> - \begin{center}
> -
> - Version 1.3, 3 November 2008
> -
> -
> - Copyright \copyright{} 2000, 2001, 2002, 2007, 2008 Free Software Foundation, Inc.
> -
> - \bigskip
> -
> - <http://fsf.org/>
> -
> - \bigskip
> -\end{center}
> -
> - Everyone is permitted to copy and distribute verbatim copies
> - of this license document, but changing it is not allowed.
> -
> +\addtocounter{chapter}{1}
> +\vspace{-1cm}
> +
> +Version 1.3, 3 November 2008
> +
> +Copyright \copyright{} 2000, 2001, 2002, 2007, 2008 Free Software Foundation, Inc.
> +
> +\bigskip
> +
> + <http://fsf.org/>
> +
> +\bigskip
> +
> +Everyone is permitted to copy and distribute verbatim copies of this
> +license document, but changing it is not allowed.
>
> \def\thesection{\arabic{section}.}
> \setcounter{section}{-1}
>
> === modified file 'chapters/lopes/chapter.tex'
> --- chapters/lopes/chapter.tex 2011-04-29 09:32:01 +0000
> +++ chapters/lopes/chapter.tex 2011-05-13 10:11:49 +0000
> @@ -3,41 +3,42 @@
> {Nuno D. Lopes, Pedro Jorge Silva Pereira and Lu{\'\i}s Trabucho}
> {lopes}
>
> -The main motivation of this work is the implementation of a
> -general finite element solver for ...

=== modified file 'STATUS'
--- STATUS	2011-05-13 08:16:02 +0000
+++ STATUS	2011-05-13 10:11:49 +0000
@@ -95,7 +95,7 @@
 Nuno D. Lopes <[e-mail]>             6       [kirby-5],                 APPLIED
 Kent-Andre Mardal <[e-mail]>         13, 14  [kirby-3], [kirby-2]       SUBMITTED (AL will fix), SUBMITTED (RCK will fix)
 Mikael Mortensen <[e-mail]>     32      [nikbakht]
-Harish Narayanan <[e-mail]>            8, 27   [oelgaard-2], [lopes]      SUBMITTED, WAITING
+Harish Narayanan <[e-mail]>            8, 27   [oelgaard-2], [lopes]      APPLIED, SUBMITTED
 Murtazo Nazarov <[e-mail]>           21      [terrel]
 Mehdi Nikbakht <[e-mail]>         36      [lezar]                    APPLIED
 Johannes Ring <[e-mail]>              7, 33   [logg-3], [schroll]        APPLIED, APPLIED

=== modified file 'backmatter/authors.tex'
--- backmatter/authors.tex	2011-04-29 11:21:33 +0000
+++ backmatter/authors.tex	2011-05-13 10:11:49 +0000
@@ -1,5 +1,7 @@
 \chapter*{List of authors}
 \addcontentsline{toc}{chapter}{List of authors}
+\addtocounter{chapter}{1}
+\vspace{-1cm}
 
 \emph{The following authors have contributed to this book.}
 

=== modified file 'backmatter/fdl-1.3.tex'
--- backmatter/fdl-1.3.tex	2011-05-10 18:26:12 +0000
+++ backmatter/fdl-1.3.tex	2011-05-13 10:11:49 +0000
@@ -1,25 +1,20 @@
-\chapter*{\rlap{GNU Free Documentation License}}
-\phantomsection  % so hyperref creates bookmarks
+\chapter*{GNU Free Documentation License}
 \addcontentsline{toc}{chapter}{GNU Free Documentation License}
-%\label{label_fdl}
-
- \begin{center}
-
-       Version 1.3, 3 November 2008
-
-
- Copyright \copyright{} 2000, 2001, 2002, 2007, 2008  Free Software Foundation, Inc.
-
- \bigskip
-
-     <http://fsf.org/>
-
- \bigskip
-\end{center}
-
- Everyone is permitted to copy and distribute verbatim copies
- of this license document, but changing it is not allowed.
-
+\addtocounter{chapter}{1}
+\vspace{-1cm}
+
+Version 1.3, 3 November 2008
+
+Copyright \copyright{} 2000, 2001, 2002, 2007, 2008  Free Software Foundation, Inc.
+
+\bigskip
+
+ <http://fsf.org/>
+
+\bigskip
+
+Everyone is permitted to copy and distribute verbatim copies of this
+license document, but changing it is not allowed.
 
 \def\thesection{\arabic{section}.}
 \setcounter{section}{-1}

=== modified file 'chapters/lopes/chapter.tex'
--- chapters/lopes/chapter.tex	2011-04-29 09:32:01 +0000
+++ chapters/lopes/chapter.tex	2011-05-13 10:11:49 +0000
@@ -3,41 +3,42 @@
               {Nuno D. Lopes, Pedro Jorge Silva Pereira and Lu{\'\i}s Trabucho}
               {lopes}
 
-The main motivation of this work is the implementation of a
-general finite element solver for some of the improved Boussinesq
-models\index{Boussinesq models}. Here, we use an extension of the model
-proposed by~\citet{ZhaoTengCheng2004} to investigate the behavior
-of surface water waves. The equations in this model do not contain
-spatial derivatives with an order higher than 2. Some effects like
-energy dissipation and wave generation by natural phenomena or external
-physical mechanisms are also included.  As a consequence, some modified
-dispersion relations are derived for this extended model. A matrix-based
-linear stability analysis of the proposed model is presented.  It is
-shown that this model is robust with respect to instabilities related
-to steep bottom gradients.
+The main motivation of this work is the implementation of a general
+finite element solver for some of the improved Boussinesq
+models\index{Boussinesq models}. Here, we use an extension of the
+model proposed by~\citet{ZhaoTengCheng2004} to investigate the
+behavior of surface water waves. The equations in this model do not
+contain spatial derivatives with an order higher than 2. Some effects
+like energy dissipation and wave generation by natural phenomena or
+external physical mechanisms are also included.  As a consequence,
+some modified dispersion relations are derived for this extended
+model. A matrix-based linear stability analysis of the proposed model
+is presented.  It is shown that this model is robust with respect to
+instabilities related to steep bottom gradients.
 
 \section{Overview}
 
-The \fenics project, via \dolfin, \ufl and \ffc, provides good technical
-and scientific support for the implementation of large scale industrial
-models based on the finite element method. Specifically, all the finite
-element matrices and vectors are automatically generated and assembled
-by \dolfin\ and \ffc, directly from the variational formulation of the
-problem which is declared using \ufl. Moreover, \dolfin\ provides a user
-friendly interface for the libraries needed to solve the finite element
-system of equations.
+The \fenics project, via \dolfin, \ufl and \ffc, provides good
+technical and scientific support for the implementation of large scale
+industrial models based on the finite element method. Specifically,
+all the finite element matrices and vectors are automatically
+generated and assembled by \dolfin\ and \ffc, directly from the
+variational formulation of the problem which is declared using
+\ufl. Moreover, \dolfin\ provides a user friendly interface for the
+libraries needed to solve the finite element system of equations.
 
-Numerical implementation of Boussinesq equations goes back to
-the works of \citet{Peregrine1967} and \citet{Wu1981}, and later
-by the development of improved dispersion characteristics (see,
-e.g., \citet{MadsenEtAl1991,Nwogu1993,ChenLiu1994} as well as
+Numerical implementation of Boussinesq equations goes back to the
+works of \citet{Peregrine1967} and \citet{Wu1981}, and later by the
+development of improved dispersion characteristics (see, e.g.,
+\citet{MadsenEtAl1991,Nwogu1993,ChenLiu1994} as well as
 \citet{BejiNadaoka1996}).
 
-We implement a solver for some of the Boussinesq type systems to model the
-evolution of surface water waves\index{water waves} in a variable depth
-seabed.  This type of model is used, for instance, in harbor simulation
-(see Figure~\ref{fig:lopes:harbor} for an example of a standard harbor),
-tsunami generation and wave propagation as well as in coastal dynamics.
+We implement a solver for some of the Boussinesq type systems to model
+the evolution of surface water waves\index{water waves} in a variable
+depth seabed.  This type of model is used, for instance, in harbor
+simulation (see Figure~\ref{fig:lopes:harbor} for an example of a
+standard harbor), tsunami generation and wave propagation as well as
+in coastal dynamics.
 
 \begin{figure}
   \begin{center}
@@ -47,26 +48,28 @@
   \label{fig:lopes:harbor}
 \end{figure}
 
-In Section~\ref{sec:lopes:dolfwave}, we begin by describing the DOLFWAVE
-application which is a \fenics based application for the simulation of
-surface water waves (see \url{http://ptmat.fc.ul.pt/~ndl/}).
+In Section~\ref{sec:lopes:dolfwave}, we begin by describing the
+DOLFWAVE application which is a \fenics based application for the
+simulation of surface water waves (see
+\url{http://ptmat.fc.ul.pt/~ndl/}).
 
 \editornote{URL should be archival.}
 
-The governing equations for surface water waves are presented
-in Section~\ref{sec:lopes:modelderivation}. From these equations
+The governing equations for surface water waves are presented in
+Section~\ref{sec:lopes:modelderivation}. From these equations
 different types of models can be derived. There are several Boussinesq
 models and some of the most widely used are those based on the wave
 surface Elevation and horizontal Velocities formulation (BEV) (see,
 e.g., \citet{WalkleyBerzins2002}, \citet{WooLiu2004a} as well as
 \citet{WooLiu2004b} for finite element discretizations of BEV models).
 However, we only consider the wave surface Elevation and velocity
-Potential (BEP) formulation (see, e.g., ~\citet{LangtangenPedersen1998}
-for a finite element discretization of a BEP model).  Thus, the number
-of unknowns is reduced from five (the three velocity components, the
-pressure and the wave surface elevation) in the BEV models to three (the
-velocity potential, the pressure and the wave surface elevation) in the
-BEP models. Two different types of BEP models are taken into account:
+Potential (BEP) formulation (see, e.g.,
+~\citet{LangtangenPedersen1998} for a finite element discretization of
+a BEP model).  Thus, the number of unknowns is reduced from five (the
+three velocity components, the pressure and the wave surface
+elevation) in the BEV models to three (the velocity potential, the
+pressure and the wave surface elevation) in the BEP models. Two
+different types of BEP models are taken into account:
 %%
 \begin{enumerate}
 \item \label{lopes:i}  a standard  model containing sixth-order
@@ -76,14 +79,14 @@
 \end{enumerate}
 
 A standard technique is used in order to derive the Boussinesq-type
-model mentioned in~\ref{lopes:i}.  In the subsequent sections, only the
-ZTC-type model is considered.  Note that these two models are
+model mentioned in~\ref{lopes:i}.  In the subsequent sections, only
+the ZTC-type model is considered.  Note that these two models are
 complemented with some extra terms, due to the inclusion of effects
 like energy dissipation and wave generation by moving an impermeable
 bottom or using a source function.
 
-An important characteristic of the extended ZTC
-model, including dissipative effects, is presented in
+An important characteristic of the extended ZTC model, including
+dissipative effects, is presented in
 Section~\ref{sec:lopes:dispersionproperties}, namely, the dispersion
 relation.
 
@@ -165,7 +168,6 @@
 model, specifically C-grid is implemented for the spatial
 discretization and a leap-frog scheme is used for the time stepping.
 
-
 \section{DOLFWAVE}
 \label{sec:lopes:dolfwave}
 
@@ -179,8 +181,8 @@
 We have already implemented solvers for the following cases:
 %%
 \begin{enumerate}
-\item \label{lopes:list_i} Shallow water wave models for unidirectional
-  long waves in one horizontal dimension;
+\item \label{lopes:list_i} Shallow water wave models for
+  unidirectional long waves in one horizontal dimension;
 \item  Boussinesq-type models for moderately long
   waves with small amplitude in shallow water.
 \end{enumerate}
@@ -230,7 +232,7 @@
 We use \ufl\ for the declaration of the finite element discretization
 of the variational forms related to the models mentioned above (see
 the \ufl\ form files in the following directories of the DOLFWAVE code
-tree: \emp{dolfwave/src/1hd1sforms}, \emp{dolfwave/src/1hdforms} and
+tree: \emp{dolfwave/src/1hd1sforms}, \emp{dolfwave/src/1hdforms} and\\
 \emp{dolfwave/src/2hdforms}).  These files are compiled using \ffc to
 generate the C++ code of the finite element discretization of the
 variational forms (see Chapter~\ref{chap:alnes-1}).  DOLFWAVE is based
@@ -242,12 +244,11 @@
 the C++ code for boundary conditions and source functions are
 included. Scripts for visualization and data analysis are also part of
 the application. The Xd3d post-processor is used in some cases (see
-\url{http://www.cmap.polytechnique.fr/~jouve/xd3d/}).  DOLFWAVE
-has a large number of demos covering all the implemented models (see
+\url{http://www.cmap.polytechnique.fr/~jouve/xd3d/}).  DOLFWAVE has a
+large number of demos covering all the implemented models (see
 \emp{dolfwave/demo}). Different physical effects are illustrated.  All
 the numerical examples in this work are included in the demos.
 
-
 \section{Model derivation}
 \label{sec:lopes:modelderivation}
 
@@ -387,7 +388,7 @@
 dispersive terms.  For simplicity, in what follows, we drop the prime
 notation.
 
-The Boussinesq approach consists in reducing a $3$D problem to a $2$D
+The Boussinesq approach consists of reducing a $3$D problem to a $2$D
 one.  This may be accomplished by expanding the velocity potential in
 a Taylor power series in terms of $z$.  Using Laplace's equation, in a
 dimensionless form, we can obtain the following expression for the
@@ -419,19 +420,19 @@
     and}\ \varepsilon<1.
 \end{equation}
 %%
-Note that the Ursell number is defined by
-$U_r = \varepsilon/\mu^2$ and plays a central role in deciding
-the choice of approximations which correspond to very different
-physics.  The regime of weakly nonlinear, small amplitude and
-moderately long waves in shallow water is characterized by
-\eqref{eq:lopes:ursell} $\left(O(\mu^2)=O(\varepsilon)\right.$; that
-is, $\left.\frac{H^2}{L^2}\sim\frac{A}{H} \right)$. Boussinesq
-equations account for the effects of nonlinearity $\varepsilon$ and
-dispersion $\mu^2$ to the leading order.  When $\varepsilon\gg \mu^2$
-, they reduce to the Airy equations.  When $\varepsilon\ll \mu^2$ they
-reduce to the linearized approximation with weak dispersion.  Finally,
-if we assume that $\varepsilon\rightarrow 0$ and $\mu^2\rightarrow 0$,
-the classical linearized wave equation is obtained.
+Note that the Ursell number is defined by $U_r = \varepsilon/\mu^2$
+and plays a central role in deciding the choice of approximations
+which correspond to very different physics.  The regime of weakly
+nonlinear, small amplitude and moderately long waves in shallow water
+is characterized by \eqref{eq:lopes:ursell}
+$\left(O(\mu^2)=O(\varepsilon)\right.$; that is,
+$\left.\frac{H^2}{L^2}\sim\frac{A}{H} \right)$. Boussinesq equations
+account for the effects of nonlinearity $\varepsilon$ and dispersion
+$\mu^2$ to the leading order.  When $\varepsilon\gg \mu^2$ , they
+reduce to the Airy equations.  When $\varepsilon\ll \mu^2$ they reduce
+to the linearized approximation with weak dispersion.  Finally, if we
+assume that $\varepsilon\rightarrow 0$ and $\mu^2\rightarrow 0$, the
+classical linearized wave equation is obtained.
 
 A sixth-order spatial derivative model is obtained if $\phi_1$ is
 expanded in terms of $\phi_0$ and all terms up to $O(\mu^8)$ are
@@ -667,7 +668,7 @@
     several models.}
   \label{fig:lopes:dispersion}
 \end{figure}
-%%
+
 From Figure~\ref{fig:lopes:dispersion}, we can also see that
 these two dissipative models admit critical wave numbers $k_1$ and
 $k_2$, such that the positive part of $\displaystyle
@@ -1074,7 +1075,7 @@
 reference for comparison.  For both models, full reflective boundary
 conditions are considered.
 
-We start by assuming  a separated  solution of the form
+We start by assuming a separated solution of the form
 %%
 \begin{equation}
 \label{eq:lopes:exp}

=== modified file 'chapters/narayanan/chapter.tex'
--- chapters/narayanan/chapter.tex	2011-05-03 21:59:00 +0000
+++ chapters/narayanan/chapter.tex	2011-05-13 10:11:49 +0000
@@ -29,8 +29,8 @@
 mechanical response of many polymeric and biological materials. Such
 materials are capable of undergoing finite deformation, and their
 material response is often characterized by complex, nonlinear
-constitutive relationships. (See, for example, \citep{Holzapfel2000}
-and \citep{TruesdellNoll1965} and the references within for several
+constitutive relationships. (See, for example, \citet{Holzapfel2000}
+and \citet{TruesdellNoll1965} and the references within for several
 examples.) Because of these difficulties, predicting the response of
 arbitrary structures composed of such materials to arbitrary loads
 requires numerical computation, usually based on the finite element
@@ -91,13 +91,14 @@
 are parametrized by reference positions. This is commonly termed the
 {\em material} or {\em Lagrangian} description.
 
-In its most basic terms, the {\em deformation} of the body over a time
+In its most basic terms, the {\em deformation} of the body over time
 $t \in [0, T]$ is a sufficiently smooth bijective map $\Bvarphi:
 \overline{\Omega} \times [0,T] \rightarrow \mathbb{R}^{2, 3}$, where
 $\overline{\Omega} :=\ \overline{\Omega\cup\partial\Omega}$ and
 $\partial\Omega$ is the boundary of $\Omega$. The restrictions on the
-map ensure that the motion it describes is physical (e.g., disallowing
-the interpenetration of matter or the formation of cracks). From this,
+map ensure that the motion it describes is physical and within the
+range of applicability of the theory (e.g., disallowing the
+interpenetration of matter or the formation of cracks). From this map,
 we can construct the {\em displacement field},
 %%
 \begin{equation}
@@ -159,16 +160,16 @@
   \label{eq:narayanan:balance-of-momentum}
 \end{equation}
 %%
-where $\rho$ is the {\em reference density} of the body,
-$\bP$ is the {\em first Piola--Kirchhoff stress tensor}, $\mathrm{Div}
-(\cdot)$ is the divergence operator and $\bB$ is the body force per
-unit volume. Along with \eqref{eq:narayanan:balance-of%
--momentum}, we have initial conditions $\bu(\bX, 0) = \bu_{0}(\bX)$
-and $\frac{\partial\bu}{\partial t}(\bX, 0) = \bv_{0}(\bX)$ in
-$\Omega$, and boundary conditions $\bu(\bX, t) = \bg(\bX, t)$ on
-$\partial\Omega_{\mathrm{D}}$ and $\bP \bN = \bT$ on
-$\partial\Omega_{\mathrm{N}}$. Here, $\bN$ is the outward normal
-on the boundary.
+where $\rho$ is the {\em reference density} of the body, $\bP$ is the
+{\em first Piola--Kirchhoff stress tensor}, $\mathrm{Div} (\cdot)$ is
+the divergence operator and $\bB$ is the body force per unit
+volume. Along with \eqref{eq:narayanan:balance-of-momentum}, we have
+initial conditions $\bu(\bX, 0) = \bu_{0}(\bX)$ and
+$\frac{\partial\bu}{\partial t}(\bX, 0) = \bv_{0}(\bX)$ in $\Omega$,
+and boundary conditions $\bu(\bX, t) = \bg(\bX, t)$ on
+$\partial\Omega_{\mathrm{D}}$ and $\bP(\bX, t) \bN(\bX) = \bT(\bX, t)$
+on $\partial\Omega_{\mathrm{N}}$. Here, $\bN(X)$ is the outward normal
+on the boundary at the point $\bX$.
 
 We focus on the balance of linear momentum because, in a continuous
 sense, the other fundamental balance principles that materials must
@@ -190,14 +191,14 @@
 
 \subsection{Accounting for different materials}
 
-It is important to reiterate that
-\eqref{eq:narayanan:balance-of-momentum} is valid for all
-materials. In order to differentiate between different materials and
-to characterize their specific mechanical responses, the theory turns
-to {\em constitutive relationships}, which are models for describing
-the real mechanical behavior of matter. In the case of nonlinear
-elastic (or {\em hyperelastic}) materials, this description is posed
-in the form of a stress-strain relationship through an objective and
+It is important to reiterate that \eqref{eq:narayanan:balance%
+-of-momentum} is valid for all materials. In order to differentiate
+between different materials and to characterize their specific
+mechanical responses, the theory turns to {\em constitutive
+relationships}, which are models for describing the real mechanical
+behavior of matter. In the case of nonlinear elastic (or {\em
+hyperelastic}) materials, this description is usually posed in the
+form of a stress-strain relationship through an objective and
 frame-indifferent Helmholtz free energy function called the {\em
 strain energy function}, $\psi$. This is an energy defined per unit
 reference volume and is solely a function of the local {\em strain
@@ -216,8 +217,10 @@
 \\
 Deformation gradient &  $\bF = \bone + \mathrm{Grad}(\bu)$
 \\
-Right Cauchy--Green tensor & $\bC = \bF^{\mathrm{T}} \bF$\\
-Green--Lagrange strain tensor & $\bE = \frac{1}{2} \left(\bC - \bone\right)$
+Right Cauchy--Green tensor & $\bC = \bF^{\mathrm{T}} \bF$
+\\
+Green--Lagrange strain tensor & $\bE = \frac{1}{2} \left(\bC -
+  \bone\right)$
 \\
 Left Cauchy--Green tensor & $\bb = \bF \bF^{\mathrm{T}}$
 \\
@@ -300,7 +303,7 @@
 following {\em constitutive relationship}:
 %%
 \begin{equation}
-\bS = \bF^{-1} \frac{\partial \psi(\bF)}{\partial \bF}.
+  \bS = \bF^{-1} \frac{\partial \psi(\bF)}{\partial \bF}.
 \end{equation}
 %%
 The second Piola--Kirchhoff stress tensor is related to the first
@@ -309,25 +312,25 @@
 
 As already mentioned, the strain energy function can be posed in
 equivalent forms in terms of different strain measures. (Again, the
-reader is directed to classical texts to motivate this.) In order to
-then arrive at the second Piola--Kirchhoff stress tensor, we turn to
-the chain rule of differentiation. For example,
+interested reader is directed to classical texts to motivate this.) In
+order to then arrive at the second Piola--Kirchhoff stress tensor, we
+turn to the chain rule of differentiation. For example,
 %%
 \begin{equation}
-\begin{split}
-  \bS & = 2 \frac{\partial \psi(\bC)}{\partial \bC} = \frac{\partial
-    \psi(\bE)}{\partial \bE}\\
-      & = 2\left[\left(\frac{\partial \psi (I_{1}, I_{2},
-            I_{3})}{\partial I_{1}} + I_{1} \frac{\partial \psi (I_{1}, I_{2},
-            I_{3})}{\partial I_{2}} \right) \bone - \frac{\partial \psi (I_{1},
-          I_{2}, I_{3})}{\partial I_{2}} \bC + I_{3} \frac{\partial \psi (I_{1},
-          I_{2}, I_{3})}{\partial I_{3}} \bC^{-1} \right]\\
-      & = \sum_{A = 1}^{3}
-      \frac{1}{\lambda_{A}} \frac{\partial \psi(\lambda_{1},
-        \lambda_{2}, \lambda_{3})}{\partial
-        \lambda_{A}} \bN^{A}\otimes\bN^{A} = \ldots
-\end{split}
-\label{eq:narayanan:secondpk}
+  \begin{split}
+    \bS & = 2 \frac{\partial \psi(\bC)}{\partial \bC} =
+    \frac{\partial \psi(\bE)}{\partial \bE}\\
+    & = 2\left[\left(\frac{\partial \psi (I_{1}, I_{2},
+          I_{3})}{\partial I_{1}} + I_{1} \frac{\partial \psi (I_{1},
+          I_{2}, I_{3})}{\partial I_{2}} \right) \bone -
+      \frac{\partial \psi (I_{1}, I_{2}, I_{3})}{\partial I_{2}} \bC +
+      I_{3} \frac{\partial \psi (I_{1}, I_{2}, I_{3})}{\partial I_{3}}
+      \bC^{-1} \right]\\
+    & = \sum_{A = 1}^{3} \frac{1}{\lambda_{A}} \frac{\partial
+      \psi(\lambda_{1}, \lambda_{2}, \lambda_{3})}{\partial
+      \lambda_{A}} \bN^{A}\otimes\bN^{A} = \ldots
+  \end{split}
+  \label{eq:narayanan:secondpk}
 \end{equation}
 
 Using definitions such as the ones explicitly provided in
@@ -335,12 +338,12 @@
 Piola--Kirchhoff stress tensor from the strain energy function by
 suitably differentiating it with respect to the appropriate strain
 measure. This allows the user to easily specify material models in
-terms of each of the strain measures introduced in
-Table~\ref{tab:narayanan:straindefs}. The base class for all material
-models, \emp{MaterialModel}, encapsulates this functionality. The
-relevant method of this class is provided in
-Figure~\ref{code:narayanan:material_model_base.py}. The implementation
-relies heavily on the UFL \emp{diff} operator.
+terms of each of the strain measures introduced in Table~\ref{tab:%
+narayanan:straindefs}. The base class for all material models,
+\emp{MaterialModel}, encapsulates this functionality. The relevant
+method of this class is provided in Figure~\ref{code:narayanan:%
+material_model_base.py}. The implementation relies heavily on the UFL
+\emp{diff} operator.
 
 \begin{figure}
 \begin{python}
@@ -370,7 +373,8 @@
 \end{python}
 \label{code:narayanan:material_model_base.py}
 \caption{Partial listing of the method that suitably computes the
-  second Piola--Kirchhoff stress tensor based on the strain measure.}
+  second Piola--Kirchhoff stress tensor based on the chosen strain
+  measure.}
 \end{figure}
 
 The generality of the material model base class allows for the (almost
@@ -448,50 +452,49 @@
 \subsection{The finite element formulation of the balance of linear
   momentum}
 
-By taking the dot product of
-\eqref{eq:narayanan:balance-of-momentum} with a test function
-$\bv \in \hat{\mathrm{V}}$ and integrating over the reference domain
-and time, we have
+By taking the dot product of \eqref{eq:narayanan:balance-of-momentum}
+with a test function $\bv \in \hat{\mathrm{V}}$ and integrating over
+the reference domain and time, we have
 %%
 \begin{equation}
- \int_{0}^{T} \int_{\Omega} \rho \frac{\partial^{2}
+  \int_{0}^{T} \int_{\Omega} \rho \frac{\partial^{2}
     \bu}{\partial t^{2}} \cdot \bv \dx \dt\;
   = \int_{0}^{T}\int_{\Omega} \mathrm{Div}(\bP) \cdot \bv \dx \dt\;
   + \int_{0}^{T} \int_{\Omega} \bB \cdot \bv \dx \dt.
 \end{equation}
 %%
 Noting that the traction vector $\bT = \bP \bN$ on
-$\partial\Omega_{\mathrm{N}}$ ($\bN$ being the outward normal
-on the boundary) and that by definition
-$\bv|_{\partial\Omega_{\mathrm{D}}} = 0$, we apply the divergence
-theorem to arrive at the following weak form of the balance of linear
-momentum: \\*
-Find $\bu \in \mathrm{V}$, such that $\foralls \bv \in \hat{\mathrm{V}}$:
+$\partial\Omega_{\mathrm{N}}$ ($\bN$ being the outward normal on the
+boundary) and that by definition $\bv|_{\partial\Omega_{\mathrm{D}}} =
+0$, we apply the divergence theorem to arrive at the following weak
+form of the balance of linear momentum: \\*
+Find $\bu \in \mathrm{V}$, such that $\foralls \bv \in
+\hat{\mathrm{V}}$:
+%%
 \begin{equation}
   \int_{0}^{T} \int_{\Omega} \rho \frac{\partial^{2} \bu}{\partial
-t^{2}} \cdot \bv \dx \dt\; + \int_{0}^{T} \int_{\Omega}
-\bP:\mathrm{Grad}(\bv) \dx \dt\; = \int_{0}^{T} \int_{\Omega} \bB
-\cdot \bv \dx \dt\; + \int_{0}^{T} \int_{\partial\Omega_{\mathrm{N}}}
-\bT \cdot \bv \ds \dt,
+    t^{2}} \cdot \bv \dx \dt\; + \int_{0}^{T} \int_{\Omega}
+  \bP:\mathrm{Grad}(\bv) \dx \dt\; = \int_{0}^{T} \int_{\Omega} \bB
+  \cdot \bv \dx \dt\; + \int_{0}^{T} \int_{\partial\Omega_{\mathrm{N}}}
+  \bT \cdot \bv \ds \dt,
 \label{eq:narayanan:weakform}
 \end{equation}
 %%
 with initial conditions $\bu(\bX, 0) = \bu_{0}(\bX)$ and
-$\frac{\partial\bu}{\partial t}(\bX, 0) = \bv_{0}(\bX)$ in
-$\Omega$, and boundary conditions $\bu(\bX, t) = \bg(\bX, t)$ on
+$\frac{\partial\bu}{\partial t}(\bX, 0) = \bv_{0}(\bX)$ in $\Omega$,
+and boundary conditions $\bu(\bX, t) = \bg(\bX, t)$ on
 $\partial\Omega_{\mathrm{D}}$.
 
 The finite element formulation implemented in \twist{} follows the
-Galerkin approximation of the above weak
-form~\eqref{eq:narayanan:weakform}, by looking for solutions in a finite
-solution space  $\mathrm{V}_{h} \subset \mathrm{V}$ and allowing for
-test functions in a finite approximation of the test space
-$\hat{\mathrm{V}}_{h} \subset \hat{\mathrm{V}}$.\footnote{We now
-  note an inherent advantage in choosing the Lagrangian description in
-  formulating the theory. The fact that the integrals in
-  \eqref{eq:narayanan:weakform}, along with the various fields
-  and differential operators, are defined over the fixed domain
-  $\Omega$ means that one need not be concerned with the
+Galerkin approximation of the above weak form~\eqref{eq:narayanan:%
+weakform}, by looking for solutions in a finite solution space
+$\mathrm{V}_{h} \subset \mathrm{V}$ and allowing for test functions in
+a finite approximation of the test space $\hat{\mathrm{V}}_{h} \subset
+\hat{\mathrm{V}}$.\footnote{We now note an inherent advantage in
+  choosing the Lagrangian description in formulating the theory. The
+  fact that the integrals in \eqref{eq:narayanan:weakform}, along with
+  the various fields and differential operators, are defined over the
+  fixed domain $\Omega$ means that one need not be concerned with the
   complexity associated with calculations on a moving computational
   domain when implementing this formulation.}
 
@@ -507,34 +510,33 @@
 \label{eq:narayanan:staticweakform}
 \end{equation}
 %%
-Since \twist{} provides the necessary functionality to
-easily compute the first Piola--Kirchhoff stress tensor, $\bP$, given
-a displacement field, $\bu$, for arbitrary material models,
-\eqref{eq:narayanan:staticweakform} is just a nonlinear functional in
-terms of $\bu$. The automatic differentiation capabilities of
-UFL\footnote{An earlier chapter on UFL~(\ref{chap:alnes-1}) provides a
-detailed look at the capabilities of UFL, as well as insights into how
-it achieves its functionality. Even so, we note the following
-differentiation capabilities of UFL because of their pivotal relevance
-to this work:
+Since \twist{} provides the necessary functionality to easily compute
+the first Piola--Kirchhoff stress tensor, $\bP$, given a displacement
+field, $\bu$, for arbitrary material models, \eqref{eq:narayanan:%
+staticweakform} is just a nonlinear functional in terms of $\bu$. The
+automatic differentiation capabilities of UFL\footnote{An earlier
+  chapter on UFL~(\ref{chap:alnes-1}) provides a detailed look at the
+  capabilities of UFL, as well as insights into how it achieves its
+  functionality. Even so, we note the following differentiation
+  capabilities of UFL because of their pivotal relevance to this work:
 %%
-\begin{itemize}
-\item Computing spatial derivatives of fields, which allows for the
-  construction of differential operators such as such as
-  $\mathrm{Grad(\cdot)}$ or $\mathrm{Div(\cdot)}$:\\
-  \texttt{df\_i = Dx(f, i)}
-\item Differentiating arbitrary expressions with respect to variables
-  they are functions of:\\
-\texttt{g = variable(cos(cell.x[0]))}\\
-\texttt{f = exp(g**2)}\\
-\texttt{h = diff(f, g)}
-\item Differentiating forms with respect to coefficients of a discrete
-  function, allowing for automatic linearizations of nonlinear
-  variational forms:\\
-\texttt{a = derivative(L, w, u)}
-\end{itemize}} make
-this nonlinear form straightforward to implement, as evidenced by the
-code listing in Figure~\ref{code:narayanan:staticmomentumsolver}.
+  \begin{itemize}
+  \item Computing spatial derivatives of fields, which allows for the
+    construction of differential operators such as such as
+    $\mathrm{Grad(\cdot)}$ or $\mathrm{Div(\cdot)}$:\\
+    \texttt{df\_i = Dx(f, i)}
+  \item Differentiating arbitrary expressions with respect to variables
+    they are functions of:\\
+    \texttt{g = variable(cos(cell.x[0]))}\\
+    \texttt{f = exp(g**2)}\\
+    \texttt{h = diff(f, g)}
+  \item Differentiating forms with respect to coefficients of a discrete
+    function, allowing for automatic linearizations of nonlinear
+    variational forms:\\
+    \texttt{a = derivative(L, w, u)}
+  \end{itemize}} make this nonlinear form straightforward to
+implement, as evidenced by the code listing in Figure~\ref{code:%
+narayanan:staticmomentumsolver}.
 
 This listing provides the relevant section of the static balance of
 linear momentum solver class, \emp{StaticMomentumBalanceSolver}. The
@@ -597,15 +599,17 @@
 \br) \in \hat{\mathrm{V}}$:
 %%
 \begin{equation}
-\begin{split}
-  \int_{0}^{T} \int_{\Omega} \rho \frac{\partial \bw}{\partial t}
-\cdot \bv \dx \dt\; + \int_{0}^{T} \int_{\Omega}
-\bP:\mathrm{Grad}(\bv) \dx \dt\; &= \int_{0}^{T} \int_{\Omega} \bB
-\cdot \bv \dx \dt\; + \int_{0}^{T} \int_{\partial\Omega_{\mathrm{N}}}
-\bT \cdot \bv \ds \dt, \; \mathrm{and}\\ \int_{0}^{T} \int_{\Omega}
-\frac{\partial \bu}{\partial t} \cdot \br \dx \dt\; &= \int_{0}^{T}
-\int_{\Omega} \bw \cdot \br \dx \dt.
-\end{split}
+  \begin{split}
+    \int_{0}^{T} \int_{\Omega} \rho \frac{\partial \bw}{\partial t}
+    \cdot \bv \dx \dt\; + \int_{0}^{T} \int_{\Omega}
+    \bP:\mathrm{Grad}(\bv) \dx \dt\;
+    &= \int_{0}^{T} \int_{\Omega} \bB \cdot \bv \dx \dt\; +
+    \int_{0}^{T} \int_{\partial\Omega_{\mathrm{N}}} \bT \cdot \bv \ds
+    \dt, \; \mathrm{and}\\
+    \int_{0}^{T} \int_{\Omega} \frac{\partial \bu}{\partial t} \cdot
+    \br \dx \dt\;
+    &= \int_{0}^{T} \int_{\Omega} \bw \cdot \br \dx \dt.
+  \end{split}
 \label{eq:narayanan:cg1weakform}
 \end{equation}
 %%
@@ -620,14 +624,16 @@
 scheme:
 %%
 \begin{equation}
-\begin{split}
-  \int_{\Omega} \rho \frac{\left(\bw_{n+1} - \bw_{n}\right)}{\Delta t} \cdot \bv \dx\;
-  + \int_{\Omega} \bP(\bu_{\mathrm{mid}}):\mathrm{Grad}(\bv) \dx\;
-  &=  \int_{\Omega} \bB \cdot \bv \dx\;
-  +  \int_{\partial\Omega_{\mathrm{N}}} \bT \cdot \bv \ds,
-  \; \mathrm{and}\\
-    \int_{\Omega} \frac{\left(\bu_{n+1} - \bu_{n}\right)}{\Delta t} \cdot \br \dx\;
-  &=  \int_{\Omega} \bw_{\mathrm{mid}} \cdot \br \dx,
+  \begin{split}
+    \int_{\Omega} \rho \frac{\left(\bw_{n+1} - \bw_{n}\right)}{\Delta
+      t} \cdot \bv \dx\; + \int_{\Omega}
+    \bP(\bu_{\mathrm{mid}}):\mathrm{Grad}(\bv) \dx\;
+    &=  \int_{\Omega} \bB \cdot \bv \dx\;
+    +  \int_{\partial\Omega_{\mathrm{N}}} \bT \cdot \bv \ds,
+    \; \mathrm{and}\\
+    \int_{\Omega} \frac{\left(\bu_{n+1} - \bu_{n}\right)}{\Delta t}
+    \cdot \br \dx\;
+    &=  \int_{\Omega} \bw_{\mathrm{mid}} \cdot \br \dx,
 \end{split}
 \label{eq:narayanan:cg1weakformapprox}
 \end{equation}
@@ -690,7 +696,7 @@
     a = derivative(L, U, dU)
 \end{python}
 \caption{Relevant portion of the dynamic balance of linear momentum
-  balance solver using the CG$_{1}$ time-stepping scheme.}
+  solver using the CG$_{1}$ time-stepping scheme.}
 \label{code:narayanan:cg1}
 \end{figure}
 
@@ -701,7 +707,7 @@
   second example calculation.}  But it should also be noted that the mixed
 system that results from the formulation is computationally expensive
 and memory intensive as the number of variables being solved for have
-doubled.
+doubled due to the introduction of the velocity variable.
 
 \twist{} also provides a standard implementation of a finite
 difference time-stepping algorithm that is commonly used in the
@@ -727,9 +733,9 @@
 %%
 \begin{equation}
   \int_{\Omega} \rho a_{0} \cdot \bv \dx\; +  \int_{\Omega}
-\bP(u_{0}):\mathrm{Grad}(\bv) \dx\; - \int_{\Omega} \bB(X, 0)
-\cdot \bv \dx\; -  \int_{\partial\Omega_{\mathrm{N}}}
-\bT(X, 0) \cdot \bv \ds = 0.
+  \bP(u_{0}):\mathrm{Grad}(\bv) \dx\; - \int_{\Omega} \bB(X, 0)
+  \cdot \bv \dx\; -  \int_{\partial\Omega_{\mathrm{N}}}
+  \bT(X, 0) \cdot \bv \ds = 0.
 \end{equation}
 %%
 This provides the complete initial state ($\bu_{0}, \bv_{0}, \ba_{0}$)
@@ -738,15 +744,15 @@
 following definitions:
 %%
 \begin{equation}
-\begin{split}
-\bu_{n+1} & = \bu_{n} + \Delta t \bv_{n} + \Delta t^{2} \left[
-  \left(\frac{1}{2} - \beta \right) \ba_{n} + \beta \ba_{n+1}
-\right]\\
-\bv_{n+1} & = \bv_{n} + \Delta t \left[ (1-\gamma) \ba_{n} + \gamma
-  \ba_{n+1} \right]\\
-\bu_{n+\alpha} & = (1 - \alpha) \bu_{n} + \alpha \bu_{n+1}\\
-\bv_{n+\alpha} & = (1 - \alpha) \bv_{n} + \alpha \bv_{n+1}\\
-t_{n+\alpha} &= (1 - \alpha) t_{n} + \alpha t_{n + 1}
+  \begin{split}
+    \bu_{n+1} & = \bu_{n} + \Delta t \bv_{n} + \Delta t^{2} \left[
+      \left(\frac{1}{2} - \beta \right) \ba_{n} + \beta \ba_{n+1}
+    \right]\\
+    \bv_{n+1} & = \bv_{n} + \Delta t \left[ (1-\gamma) \ba_{n} +
+      \gamma \ba_{n+1} \right]\\
+    \bu_{n+\alpha} & = (1 - \alpha) \bu_{n} + \alpha \bu_{n+1}\\
+    \bv_{n+\alpha} & = (1 - \alpha) \bv_{n} + \alpha \bv_{n+1}\\
+    t_{n+\alpha} &= (1 - \alpha) t_{n} + \alpha t_{n + 1}
 \end{split}
 \label{eq:narayanan:hhtdefs}
 \end{equation}
@@ -756,23 +762,23 @@
 %%
 \begin{equation}
    \int_{\Omega} \rho a_{n+1} \cdot \bv \dx\; + \int_{\Omega} \bP(u_{n
-+ \alpha}):\mathrm{Grad}(\bv) \dx\; - \int_{\Omega} \bB(X, t_{n +
-\alpha}) \cdot \bv \dx\; - \int_{\partial\Omega_{\mathrm{N}}} \bT(X,
-t_{n + \alpha}) \cdot \bv \ds = 0,
+     + \alpha}):\mathrm{Grad}(\bv) \dx\; - \int_{\Omega} \bB(X, t_{n +
+     \alpha}) \cdot \bv \dx\; - \int_{\partial\Omega_{\mathrm{N}}}
+   \bT(X, t_{n + \alpha}) \cdot \bv \ds = 0,
 \label{eq:narayanan:newaccn}
 \end{equation}
 %%
-we can solve for the for the only unknown variable, the
-acceleration at the next step, $\ba_{n+1}$. The acceleration solution
-to \eqref{eq:narayanan:newaccn} is then used in the
-definitions~\eqref{eq:narayanan:hhtdefs} to update to new displacement
-and velocity values, and the problem is stepped through time.
+we can solve for the for the only unknown variable, the acceleration
+at the next step, $\ba_{n+1}$. The acceleration solution to
+\eqref{eq:narayanan:newaccn} is then used in the definitions%
+~\eqref{eq:narayanan:hhtdefs} to update to new displacement and
+velocity values, and the problem is stepped through time.
 
 We close this subsection on time-stepping algorithms with one usage
 detail pertaining to \twist. By default, when solving a dynamics
 problem, \twist{} assumes that the user wants to use the HHT
-method. In case one wants to override this behavior, they can do so
-by returning \emp{"CG(1)"} in the \emp{time\_stepping} method while
+method. In case one wants to override this behavior, they can do so by
+returning \emp{"CG(1)"} in the \emp{time\_stepping} method while
 specifying the problem. Figure~\ref{code:narayanan:dynamicrelease} is
 an example showing this.
 
@@ -781,12 +787,12 @@
 The algorithms discussed thus far serve primarily to explain the
 computational framework's inner working, and are not at the level at
 which the user usually interacts with \twist{} (unless they are
-interested in extending it). In practice, the functionality of \twist{}
-is exposed to the user through two primary problem definition classes:
-\emp{StaticHyperelasticity} and \emp{Hyperelasticity}. These classes
-reside in \emp{problem\_definitions.py}, and contain numerous methods
-for defining aspects of the nonlinear elasticity problem. As their
-names suggest, these are respectively used to describe static or
+interested in extending it). In practice, the functionality of
+\twist{} is exposed to the user through two primary problem definition
+classes: \emp{StaticHyperelasticity} and \emp{Hyperelasticity}. These
+classes reside in \emp{problem\_definitions.py}, and contain numerous
+methods for defining aspects of the nonlinear elasticity problem. As
+their names suggest, these are respectively used to describe static or
 dynamic problems in nonlinear elasticity.
 
 Over the course of the following examples, we will see how various
@@ -808,9 +814,9 @@
 3.8461$~N/m$^2$ and a spatially varying $\lambda = 5.8 x_{1} + 5.7 (1
 - x_{1})$~N/m$^2$. Here, $x_{1}$ is the first coordinate of the
 reference position, $\bX$.\footnote{The numerical parameters in this
-chapter have been arbitrarily chosen for illustration of the
-framework's use. They do not necessarily correspond to a real
-material.} In order to twist the cube, the face $x_{1} = 0$ is held
+  chapter have been arbitrarily chosen for illustration of the
+  framework's use. They do not necessarily correspond to a real
+  material.} In order to twist the cube, the face $x_{1} = 0$ is held
 fixed and the opposite face $x_{1} = 1$ is rotated 60 degrees using
 the Dirichlet condition defined in
 Figure~\ref{code:narayanan:statictwist}.
@@ -834,17 +840,17 @@
   functionality that it offers.}
 \end{figure}
 %%
-The problem is completely specified by defining
-relevant methods in the user-created class \emp{Twist} (see
+The problem is completely specified by defining relevant methods in
+the user-created class \emp{Twist} (see
 Figure~\ref{code:narayanan:statictwist}), which derives from the base
 class \emp{StaticHyperelasticity}. \twist{} only requires relevant
 methods to be provided, and for the current problem, this includes
-the computational domain, Dirichlet boundary conditions and material
-model. The methods are fairly self-explanatory, but the following points
-are to be noted. Firstly, \twist{} supports spatially-varying material
-parameters. Secondly, Dirichlet boundary conditions are posed in two
-parts: the conditions themselves, and the corresponding boundaries along
-which they act.
+those that define the computational domain, Dirichlet boundary
+conditions and material model. The methods are fairly
+self-explanatory, but the following points are to be noted. Firstly,
+\twist{} supports spatially-varying material parameters. Secondly,
+Dirichlet boundary conditions are posed in two parts: the conditions
+themselves, and the corresponding boundaries along which they act.
 
 \begin{figure}
 \begin{python}
@@ -993,14 +999,15 @@
 called, we see the relaxation of the pre-twisted cube. After initial
 unwinding of the twist, the body proceeds to twist in the opposite
 direction due to inertia. This process repeats itself, and snapshots
-of the displacement over the first 0.5~s are shown in
-Figure~\ref{fig:narayanan:releasedcube}. Figure~\ref{fig:narayanan:energies}
-highlights the energy conservation of the CG$_{1}$ numerical scheme used
-to time-step this problem by totaling the kinetic energy and
+of the displacement over the first 0.5~s are shown in Figure%
+~\ref{fig:narayanan:releasedcube}. Figure~\ref{fig:narayanan:energies}
+highlights the energy conservation of the CG$_{1}$ numerical scheme
+used to time-step this problem by totaling the kinetic energy and
 potential energy of the body over the course of the
 calculation. \twist{} provides this information through the methods
 \emp{kinetic\_energy(v)} and \emp{potential\_energy(u)}, where \emp{v}
-and \emp{u} are velocity and displacement fields respectively.
+and \emp{u} are the discrete velocity and displacement fields
+respectively.
 
 \begin{figure}
   \begin{center}