Quadrature docs

Author: pbrubeck

docs/source/element_list.py

@@ -1,5 +1,4 @@
 from finat.ufl.elementlist import ufl_elements
-# ~ from ufl.finiteelement.elementlist import ufl_elements
 from finat.element_factory import supported_elements
 import csv
 

docs/source/extruded-meshes.rst

@@ -304,7 +304,7 @@ supports a more general syntax:
 
     V = FunctionSpace(mesh, element)
 
-where ``element`` is a UFL :py:class:`~ufl.finiteelement.finiteelement.FiniteElement` object. This
+where ``element`` is a UFL :py:class:`~finat.ufl.finiteelement.FiniteElement` object. This
 requires generation and manipulation of ``FiniteElement`` objects.
 
 Geometrically, an extruded mesh cell is the *product* of a base, "horizontal",
@@ -319,7 +319,7 @@ The ``TensorProductElement`` operator
 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 
 To create an element compatible with an extruded mesh, one should use
-the :py:class:`~ufl.finiteelement.tensorproductelement.TensorProductElement`
+the :py:class:`~finat.ufl.tensorproductelement.TensorProductElement`
 operator. For example,
 
 .. code-block:: python3
@@ -359,7 +359,7 @@ but is piecewise constant horizontally.
 
 A more complicated element, like a Mini horizontal element with linear
 variation in the vertical direction, may be built using the
-:py:class:`~ufl.finiteelement.enrichedelement.EnrichedElement` functionality
+:py:class:`~finat.ufl.enrichedelement.EnrichedElement` functionality
 in either of the following ways:
 
 .. code-block:: python3
@@ -392,7 +392,7 @@ The ``HDivElement`` and ``HCurlElement`` operators
 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 
 For moderately complicated vector-valued elements,
-:py:class:`~ufl.finiteelement.tensorproductelement.TensorProductElement`
+:py:class:`~finat.ufl.tensorproductelement.TensorProductElement`
 does not give enough information to unambiguously produce the desired
 space. As an example, consider the lowest-order *Raviart-Thomas* element on a
 quadrilateral. The degrees of freedom live on the facets, and consist of
@@ -417,7 +417,7 @@ However, this is only scalar-valued. There are two natural vector-valued
 elements that can be generated from this: one of them preserves tangential
 continuity between elements, and the other preserves normal continuity
 between elements. To obtain the Raviart-Thomas element, we must use the
-:py:class:`~ufl.finiteelement.hdivcurl.HDivElement` operator:
+:py:class:`~finat.ufl.hdivcurl.HDivElement` operator:
 
 .. code-block:: python3
 
@@ -433,7 +433,7 @@ between elements. To obtain the Raviart-Thomas element, we must use the
   :align: center
 
   The RT quadrilateral element, requiring the use
-  of :py:class:`~ufl.finiteelement.hdivcurl.HDivElement`
+  of :py:class:`~finat.ufl.hdivcurl.HDivElement`
 
 Another reason to use these operators is when expanding a vector into a
 higher dimensional space. Consider the lowest-order Nedelec element of the
@@ -453,7 +453,7 @@ the product of this with a continuous element on an interval:
 
 This element still only has two components. To expand this into a
 three-dimensional curl-conforming element, we must use the
-:py:class:`~ufl.finiteelement.hdivcurl.HCurlElement` operator; the syntax is:
+:py:class:`~finat.ufl.hdivcurl.HCurlElement` operator; the syntax is:
 
 .. code-block:: python3
 

docs/source/manual.rst

@@ -19,6 +19,7 @@ Manual
    duals
    interpolation
    point-evaluation
+   quadrature
    external_operators
    adjoint
    visualisation

docs/source/quadrature.rst

@@ -0,0 +1,91 @@
+.. only:: html
+
+  .. contents::
+
+Quadrature rules
+================
+
+To numerically compute the integrals in the variational formulation,
+a quadrature rule is required.
+By default, Firedrake obtains a quadrature rule by estimating the polynomial
+degree of the integrands within a :py:class:`ufl.Form`. Sometimes
+this estimate might be quite large, and a warning like this one will be raised:
+
+.. code-block::
+
+   tsfc:WARNING Estimated quadrature degree 13 more than tenfold greater than any argument/coefficient degree (max 1)
+
+Sometimes the estimated quadrature degree might be even in the hundreds or thousands,
+rendering the integration prohibitively expensive. 
+
+Specifying the quadrature rule
+------------------------------
+
+To manually override the default, the
+quadrature degree can be prescribed on each integral :py:class:`~ufl.measure.Measure`,
+
+.. code-block::
+
+   inner(sin(u)**4, v) * dx(degree=4)
+
+Setting ``degree=4`` means that the quadrature rule will be exact only for integrands
+of total polynomial degree up to 4. This, of course, will introduce a greater numerical error than the default.
+
+The default quadrature degree estimation may be overridden by setting the ``"quadrature_degree"`` in
+the ``form_compiler_parameters`` dictionary passed on to
+:py:func:`~.solve`, :py:func:`~.project`, or :py:class:`~.NonlinearVariationalProblem`.
+
+.. code-block::
+
+   F = inner(grad(u), grad(v))*dx(degree=0) + inner(exp(u), v)*dx - inner(1, v)*dx
+
+   solve(F == 0, u, form_compiler_parameters={"quadrature_degree": 4})
+
+In the example above, only the integrals with unspecified quadrature degree
+will be computed on a quadrature rule that exactly integrates polynomials of
+the degree set in ``form_compiler_parameters``.
+
+Another way to specify the quadrature rule is through the ``scheme`` keyword. This could be
+either a :py:class:`~finat.quadrature.QuadratureRule`, or a string. Supported string values
+are ``"default"``, ``"canonical"``, and ``"lump"``. For more details see
+:py:func:`~FIAT.quadrature_schemes.create_quadrature`.  
+
+Lumped quadrature schemes
+-------------------------
+
+Spectral elements, such as Gauss-Legendre-Lobatto and KMV, may be used with
+lumped quadrature schemes to produce a diagonal mass matrix.
+
+.. code-block::
+
+   degree = 3
+   V = FunctionSpace(mesh, "KMV", degree)
+   u = TrialFunction(V)
+   v = TestFunction(V)
+
+   inner(u, v) * dx(scheme="lump", degree=degree)
+
+.. Note::
+
+   To obtain the lumped mass matrix with ``scheme="lump"``,
+   the ``degree`` argument should match the degree of the :py:class:`~finat.ufl.finiteelement.FiniteElement`.
+
+
+The Quadrature Space
+--------------------
+
+It is possible to define a finite element :py:class:`~.Function` on a quadrature rule.
+The ``"Quadrature"`` and ``"Boundary Quadrature"`` spaces are useful to
+interpolate data at quadrature points on cell interiors and cell boundaries,
+respectively.
+
+.. code-block::
+
+   Q = FunctionSpace(mesh, "Quadrature", degree=quad_degree, quad_scheme=quad_scheme)
+
+The ``quad_scheme`` keyword argument again may be either
+:py:class:`~finat.quadrature.QuadratureRule` or a string.
+If a :py:class:`~.Function` in the ``"Quadrature"`` space appears within an
+integral, Firedrake will automatically select the quadrature rule that corresponds
+to ``dx(degree=quad_degree, scheme=quad_scheme)`` to match the one associated
+with the quadrature space.

docs/source/variational-problems.rst

@@ -180,7 +180,7 @@ family such as ``"Raviart-Thomas"``. Secondly, you may use the
 family, which gives a vector-valued space where each component is
 identical to the appropriate scalar-valued
 :py:class:`~.FunctionSpace`.  Thirdly, you can create a
-:py:class:`~ufl.classes.VectorElement` directly (which is itself
+:py:class:`~finat.ufl.mixedelement.VectorElement` directly (which is itself
 *vector-valued* and pass that to the :py:func:`~.FunctionSpace`
 constructor).
 
@@ -275,7 +275,7 @@ Firedrake supports the use of the following finite elements.
     :file: element_list.csv
 
 In addition, the
-:py:class:`~ufl.finiteelement.tensorproductelement.TensorProductElement`
+:py:class:`~finat.ufl.tensorproductelement.TensorProductElement`
 operator can be used to create product elements on extruded meshes.
 
 Element variants
@@ -289,16 +289,17 @@ continuous elements these are the Gauss-Lobatto-Legendre points.
 For CG and DG spaces on simplices, Firedrake offers both equispaced points and
 the better conditioned recursive Legendre points from :cite:`Isaac2020` via the
 `recursivenodes`_ module. These are selected by passing ``variant="equispaced"``
-or ``variant="spectral"`` to the :py:class:`~ufl.classes.FiniteElement` or
+or ``variant="spectral"`` to the :py:class:`~finat.ufl.finiteelement.FiniteElement` or
 :py:func:`~.FunctionSpace` constructors. For example:
 
 .. code-block:: python3
 
     fe = FiniteElement("RTCE", quadrilateral, 2, variant="equispaced")
 
-For CG and DG spaces, and most tensor-product elements, the default is the spectral variant.
+For CG and DG spaces, and most tensor-product elements, the default is ``variant="spectral"``.
 
-Integral-type degrees of freedom are also available through ``variant="integral"``. For instance,
+Integral-type degrees of freedom are also available through ``variant="integral"``. These
+enable orthogonal bases for DG on any cell type:
 
 .. code-block:: python3
 
@@ -310,7 +311,7 @@ can be increased also through the variant option. For instance, ``variant="integ
 will use a quadrature that exactly evaluates the degrees of freedom on polynomials of
 4 degrees higher than that of the space.
 
-The key role of integral-type degrees of freedom is to nearly preserve curl-free or
+Integral-type degrees of freedom are useful to nearly preserve curl-free or
 divergence-free functions when interpolating into the finite element space.
 H(curl) and H(div) elements, such as Nédélec, Brezzi-Douglas-Marini, and Raviart-Thomas,
 have ``variant="integral"`` as default. These elements also support ``variant="point"``,
@@ -710,50 +711,6 @@ we can write:
        c.assign(t)
 
 
-Specifying the quadrature rule
-------------------------------
-
-To numerically compute the integrals in the variational formulation,
-a quadrature rule is required.
-By default, Firedrake obtains a quadrature rule by estimating the polynomial
-degree of the integrands within a :py:class:`ufl.Form`. Sometimes
-this estimate might be quite large, and a warning like this one will be raised:
-
-.. code-block::
-
-   tsfc:WARNING Estimated quadrature degree 13 more than tenfold greater than any argument/coefficient degree (max 1)
-
-This might increase memory and computational requirements. To manually override the default, the
-quadrature degree can be prescribed on each integral,
-
-.. code-block::
-
-   inner(sin(u)**4, v) * dx(degree=4)
-
-This, of course, will introduce a variational crime.
-
-
-Another way to specify the quadrature rule is through the ``scheme`` keyword. This could be
-either a :py:class:`~finat.quadrature.QuadratureRule`, or a string. Supported string values
-are ``"default"``, ``"canonical"``, and ``"lump"``. For more details see
-:py:func:`~FIAT.quadrature_schemes.create_quadrature`.  Spectral elements, such as
-Gauss-Legendre-Lobatto and KMV, may be used with lumped quadrature schemes to
-produce a diagonal mass matrix.
-
-.. code-block::
-
-   degree = 3
-   V = FunctionSpace(mesh, "KMV", degree)
-   u = TrialFunction(V)
-   v = TestFunction(V)
-
-   inner(u, v) * dx(scheme="lump", degree=degree)
-
-.. Note::
-
-   For ``scheme="lump"`` the ``degree`` argument is interpreted as the degree of :py:func:`~.FunctionSpace`.
-
-
 .. _more_complicated_forms:
 
 More complicated forms