Minor demo and readme fixes (#98)

I was going through movement and saw minor things that needed to be
fixed:
- broken link in README.md
- typos and formatting issues in demos
- cross-referencing in demos

README.md

@@ -10,9 +10,9 @@ Movement is built for use with the [Firedrake finite element library](https://fi
 ## Documentation
 
 Movement has two sources of documentation:
-* The [website](mesh-adaptation.github.io), which includes long-form documentation, demos, and API documentation.
+* The [website](https://mesh-adaptation.github.io/movement/index.html), which includes long-form documentation, demos, and API documentation.
 * The [wiki](https://github.com/mesh-adaptation/mesh-adaptation-docs/wiki), which includes recommendations on how to interact with the codebase and development practices to be followed.
 
 ## Installation
 
-Forinstallation instructions, we refer to the [Wiki page](https://github.com/mesh-adaptation/mesh-adaptation-docs/wiki/Installation-Instructions).
+For installation instructions, we refer to the [Wiki page](https://github.com/mesh-adaptation/mesh-adaptation-docs/wiki/Installation-Instructions).

demos/laplacian_smoothing.py

@@ -32,7 +32,7 @@
 # i.e.,
 #
 # .. math::
-#    \frac1{h^2}(-v_{i-1} + 2v_i - v_{i+1})) = 0,
+#    \frac1{h^2}(-v_{i-1} + 2v_i - v_{i+1}) = 0,
 #
 # the left-hand side of which you might recognise as a finite difference approximation
 # of the second derivative, i.e., the Laplace operator.
@@ -61,13 +61,13 @@
 #    :align: center
 #
 # Suppose we wish to enforce a time-dependent velocity :math:`\mathbf{v}_f`
-# on the top boundary :math and see how the mesh responds. Consider the velocity
+# on the top boundary and see how the mesh responds. Consider the velocity
 #
 # .. math::
 #     \mathbf{v}_f(x,y,t) = \left[0, A\:\sin\left(\frac{2\pi t}T\right)\:\sin(\pi x)\right]
 #
 # acting only in the vertical direction, where :math:`A` is the amplitude and :math:`T`
-# is the time period. The displacement following a sinusoidal pattern along that
+# is the time period. The displacement follows a sinusoidal pattern along that
 # boundary. The boundary movement is also sinusoidal in time, such that it ramps up and
 # then reverses, with the analytical solution being that the final mesh coincides with
 # the initial mesh.
@@ -103,8 +103,8 @@ def boundary_velocity(x, t):
 #
 # To enforce this boundary velocity, we need to create a :class:`~.LaplacianSmoother`
 # instance and define a function for updating the boundary conditions. Since we are
-# going to enforce the velocity on the top boundary, we create a :class:`~.Function` to
-# represent the boundary condition values and pass this to a :class:`~.DirichletBC`
+# going to enforce the velocity on the top boundary, we create a :class:`~firedrake.function.Function` to
+# represent the boundary condition values and pass this to a :class:`~firedrake.bcs.DirichletBC`
 # object. We then define a function which updates it as time progresses. ::
 
 mover = LaplacianSmoother(mesh, timestep)
@@ -166,4 +166,8 @@ def update_boundary_velocity(t):
 print(f"l_infinity error: {linf_error:.3f} m")
 assert np.isclose(linf_error, 0.0)
 
+# .. code-block:: console
+#
+#    l_infinity error: 0.000 m
+#
 # This tutorial can be downloaded as a `Python script <laplacian_smoothing.py>`__.

demos/lineal_spring.py

@@ -122,8 +122,9 @@ def boundary_displacement(x, t):
 #    :figwidth: 60%
 #    :align: center
 #
-# To apply this boundary displacement, we need to create a :class:`~.SpringMover`
-# instance and define a function for updating the boundary conditions. ::
+# To apply this boundary displacement, we need to create a
+# :class:`~movement.spring.SpringMover` instance and define a function for updating the
+# boundary conditions. ::
 
 mover = SpringMover(mesh, timestep, method="lineal")
 top = Function(mover.coord_space)
@@ -182,6 +183,10 @@ def update_boundary_displacement(t):
 print(f"l_infinity error: {linf_error:.3f} m")
 assert linf_error < 0.002
 
+# .. code-block:: console
+#
+#    l_infinity error: 0.001 m
+#
 # Note that the mesh doesn't return to its original state quite as neatly with the lineal
 # spring method as it does with the Laplacian smoothing method. However, the result is
 # still very good (as can be seen from the plots above).
@@ -192,6 +197,12 @@ def update_boundary_displacement(t):
 print(f"Stiffness matrix shape: {K.shape}")
 print(f"Number of mesh vertices: {mesh.num_vertices()}")
 
+# .. code-block:: console
+#
+#    Stiffness matrix shape: (242, 242)
+#    Number of mesh vertices: 121
+#
+
 fig, axes = plt.subplots()
 axes.spy(K)
 plt.savefig("lineal_spring-sparsity.jpg")

demos/monge_ampere1.py

@@ -19,35 +19,37 @@
 # and physical meshes.
 #
 # Let :math:`\boldsymbol{\xi}` and :math:`\mathbf{x}` denote the coordinate fields in
-# the computational and physical domains. Skipping some of the details that can be
-# found in :cite:`McRae:2018`, of the possible mappings we choose one that takes the form
+# the computational and physical domains, respectively. Skipping some of the details
+# that can be found in :cite:`McRae:2018`, of the possible mappings we choose one that
+# takes the form
 #
 # .. math::
 #     \mathbf{x} = \boldsymbol{\xi}+\nabla_{\boldsymbol{\xi}}\phi,
 #
 # where :math:`\phi` is a convex potential function. Further, we choose the potential
-# such that is the solution of the Monge-Ampère type equation,
+# such that it is the solution of the Monge-Ampère type equation,
 #
 # .. math::
 #     m(\mathbf{x}) \det(\underline{\mathbf{I}}+\nabla_{\boldsymbol{\xi}}\nabla_{\boldsymbol{\xi}}\phi) = \theta,
 #
 # where :math:`m(\mathbf{x})` is a so-called *monitor function* and :math:`\theta` is a
 # strictly positive normalisation function. The monitor function is of key importance
 # because it specifies the desired *density* of the adapted mesh across the domain,
-# i.e., where resolution is focused. (Note that density is the reciprocal of area in 2D
-# or of volume in 3D.)
+# i.e., where resolution is focused. Note that density is the reciprocal of area in 2D
+# or of volume in 3D.
 #
-# We begin the example by importing from the namespaces of Firedrake and Movement.
+# We begin the example by importing from the namespaces of Firedrake and Movement. ::
 
-# To start with a simple example, consider a uniform mesh of the unit square. Feel free
-# to ignore the `"MOVEMENT_REGRESSION_TEST"`, as it is only used when this demo is run
-# in the test suite (to reduce its runtime). ::
 import os
 
 from firedrake import *
 
 from movement import *
 
+# To start with a simple example, consider a uniform mesh of the unit square. Feel free
+# to ignore the `"MOVEMENT_REGRESSION_TEST"`, as it is only used when this demo is run
+# in the test suite (to reduce its runtime). ::
+
 test = os.environ.get("MOVEMENT_REGRESSION_TEST")
 n = 10 if test else 20
 mesh = UnitSquareMesh(n, n)
@@ -91,7 +93,7 @@ def ring_monitor(mesh):
 
 
 # With an initial mesh and a monitor function, we are able to construct a
-# :class:`~.MongeAmpereMover` instance and adapt the mesh. By default, the Monge-Ampère
+# :class:`~movement.monge_ampere.MongeAmpereMover` instance and adapt the mesh. By default, the Monge-Ampère
 # equation is solved to a relative tolerance of :math:`10^{-8}`. However, for the
 # purposes of continuous integration testing, a tolerance of :math:`10^{-3}` is used
 # instead to further reduce the runtime. ::
@@ -134,7 +136,7 @@ def ring_monitor(mesh):
 # applied.
 #
 # In the `next demo <./monge_ampere_3d.py.html>`__, we will demonstrate
-# that the Monge Ampère method can also be applied in three dimensions.
+# that the Monge-Ampère method can also be applied in three dimensions.
 #
 # This tutorial can be dowloaded as a `Python script <monge_ampere1.py>`__.
 #

demos/monge_ampere_3d.py

@@ -15,7 +15,7 @@ def sinatan3(mesh):
     return 0.1 * sin(50 * x * z) + atan2(0.1, sin(5 * y) - 2 * x * z)
 
 
-# We will now try to use mesh movement to optimize the mesh such that it can
+# We will now try to use mesh movement to optimise the mesh such that it can
 # most accurately represent this function with limited resolution.
 # A good indicator for where resolution is required
 # is to look at the curvature of the function which can be expressed
@@ -27,13 +27,14 @@ def sinatan3(mesh):
 #    m = 1 + \alpha \frac{H(u_h):H(u_h)}{\max_{{\bf x}\in\Omega} H(u_h):H(u_h)}
 #
 # where :math:`:` indicates the inner product, i.e. :math:`\sqrt{H:H}` is the Frobenius norm
-# of :math:`H`. We have normalized such that the minimum of the monitor function is one (where
+# of :math:`H`. We have normalised such that the minimum of the monitor function is one (where
 # the error is zero), and its maximum is :math:`1 + \alpha` (where the curvature is maximal). This
-# means we can select the ratio between the largest and smallest cell volume  in the
+# means that we can select the ratio between the largest and smallest cell volume in the
 # moved mesh as :math:`1+\alpha`.
 #
 # As in the `previous Monge-Ampère demo <./monge_ampere1.py.html>`__, we use the
-# :class:`~.MongeAmpereMover` to perform the mesh movement based on this monitor. We need
+# :class:`~movement.monge_ampere.MongeAmpereMover` to perform the mesh movement based on
+# this monitor. We need
 # to provide the monitor as a callback function that takes the mesh as its
 # input. During the iterations of the mesh movement process the monitor will then
 # be re-evaluated in the (iteratively)
@@ -92,8 +93,8 @@ def monitor(mesh):
 #    :align: center
 #
 # In the `next demo <./monge_ampere_helmholtz.py.html>`__, we will demonstrate
-# how to optimize the mesh for the discretisation of a PDE with the aim to
-# minimize its discretisation error.
+# how to optimise the mesh for the discretisation of a PDE with the aim to
+# minimise its discretisation error.
 #
 # This tutorial can be dowloaded as a `Python script <monge_ampere_3d.py>`__.
 #

demos/monge_ampere_helmholtz.py

@@ -1,8 +1,8 @@
-# Using mesh movement to optimize the mesh for PDE solution
+# Using mesh movement to optimise the mesh for PDE solution
 # =========================================================
 #
-# In this demo will demonstrate how we might use mesh movement to obtain a mesh
-# that is optimized for solving a particular PDE. The general idea is that we
+# In this demo we will demonstrate how we might use mesh movement to obtain a mesh
+# that is optimised for solving a particular PDE. The general idea is that we
 # want to reduce numerical error by increasing resolution where the local error
 # is (expected to be) high and decrease it elsewhere.
 #
@@ -27,12 +27,12 @@
 #    f(x, y) &= -\nabla^2 u(x, y) + u(x, y)
 #        = \left[ -\frac{4(x-x_0)^2 + 4(y-y_0)^2}{w^4}
 #                + \frac{4}{w^2} + 1 \right]
-#          \exp\left(-\frac{(x-x_0)^2 +  (y-y_0)^2}{w^2}\right)
+#          \exp\left(-\frac{(x-x_0)^2 +  (y-y_0)^2}{w^2}\right),
 #
 # where :math:`(x_0, y_0)` is the centre of the Gaussian with width :math:`w`.
 # Note that here we *first* choose the solution :math:`u` after which we can
 # compute what rhs :math:`f` should be, by substitution in the PDE, in order
-# for :math:`u` to be the analytical solution. This so called Method of
+# for :math:`u` to be the analytical solution. This so-called Method of
 # Manufactured Solutions approach is an easy way to construct PDE
 # configurations for which we know the analytical solution, e.g. for testing
 # purposes.
@@ -101,7 +101,7 @@ def solve_helmholtz(mesh):
 #
 #    L2-norm error on initial mesh: 0.010233816824277465
 #
-# We will now try to use mesh movement to optimize the mesh to reduce
+# We will now try to use mesh movement to optimise the mesh to reduce
 # this numerical error. We use the same monitor function as
 # in the `previous Monge-Ampère demo <./monge_ampere_3d.py.html>`__
 # based on the norm of the Hessian of the solution.
@@ -167,7 +167,7 @@ def monitor(mesh):
 #   17   Min/Max 8.6775e-01   Residual 5.8585e-08   Variation (σ/μ) 9.9835e-01
 #   Converged in 17 iterations.
 #
-# Plotting the resulting mesh
+# Plotting the resulting mesh:
 
 from firedrake.pyplot import triplot
 
@@ -193,7 +193,7 @@ def monitor(mesh):
 #
 # Of course, in many practical problems we do not actually have access
 # to the exact solution. We can then use the Hessian of the numerical
-# solution in the monitor function. Calculating the Hessian we have to
+# solution in the monitor function. When calculating the Hessian we have to
 # be a bit careful however: since our numerical FEM solution :math:`u_h`
 # is a piecewise linear function, its first gradient results in a
 # piecewise *constant* function, a vector-valued constant