From 476c61c3b7c743803ea25086ac1dad0ae9b8d7a7 Mon Sep 17 00:00:00 2001 From: Abdullah Mujahid Date: Thu, 22 Feb 2024 16:44:46 +0100 Subject: [PATCH] C++ demos document generation (#3052) * Fix typo and add doc for demo (#3012) * Rename README to README.md in demo * Fix typos in doc of custom_kernel * Add basic documentation for hyperelasticity * Build documentation (#3012) * create a CMakeLists for documentation * supply INPUT and OUTPUT DIRS in Doxygen from Cmake * add FindSphinx * restructure source/index.rst to include demo.rst * api.rst was previously index.rst * Add rst files for demo (#3012) * include flags for code blocks in demo cpp file(s) * create demo explanation as rst in /doc dir * refer code snippets in demo rst file * Setup jupyter_process and biharmonic demo * add jupyter_process for cpp demo * convert biharmonic demo doc to myst * index demo.rst * copy python conf.py to cpp docs * Switch off Doxygen warn as error * Add problematic files back * Copy demo into docs subdirs * C++ demos are organized in directory * Add all demos * add jupytext specifications * update jupytext_process to include all demos * Clean jupytext header, clean jupytext_process * Switch back to ccpp.yml execution * delete CMakeLists.txt * revert Doxyfile * revert conf.py * Remove extra file * Cleaning up * Pure conversion and some clean up * Ruff format jupytext_process.py * Clean up * Switch WARN_AS_ERROR to ON in Doxyfile * conf.py cleanup * Ruff format conf.py --------- Co-authored-by: Jack S. Hale --- cpp/demo/{README => README.md} | 0 cpp/demo/biharmonic/biharmonic.py | 8 +- cpp/demo/biharmonic/main.cpp | 241 +++++++++--------- cpp/demo/custom_kernel/main.cpp | 15 +- cpp/demo/hyperelasticity/hyperelasticity.py | 8 +- cpp/demo/hyperelasticity/main.cpp | 15 +- cpp/demo/interpolation-io/main.cpp | 58 +++-- .../interpolation_different_meshes/main.cpp | 2 + cpp/demo/mixed_topology/main.cpp | 3 + cpp/demo/poisson/main.cpp | 134 +++++----- cpp/demo/poisson/poisson.py | 10 +- cpp/demo/poisson_matrix_free/main.cpp | 31 ++- cpp/demo/poisson_matrix_free/poisson.py | 2 +- cpp/doc/Doxyfile | 3 +- cpp/doc/source/api.rst | 32 +++ cpp/doc/source/conf.py | 21 +- cpp/doc/source/demo.rst | 26 ++ cpp/doc/source/index.rst | 16 +- cpp/doc/source/jupytext_process.py | 63 +++++ 19 files changed, 421 insertions(+), 267 deletions(-) rename cpp/demo/{README => README.md} (100%) create mode 100644 cpp/doc/source/api.rst create mode 100644 cpp/doc/source/demo.rst create mode 100644 cpp/doc/source/jupytext_process.py diff --git a/cpp/demo/README b/cpp/demo/README.md similarity index 100% rename from cpp/demo/README rename to cpp/demo/README.md diff --git a/cpp/demo/biharmonic/biharmonic.py b/cpp/demo/biharmonic/biharmonic.py index 09c2b4e5a3f..f47c263983a 100644 --- a/cpp/demo/biharmonic/biharmonic.py +++ b/cpp/demo/biharmonic/biharmonic.py @@ -1,9 +1,7 @@ -# UFL input for the Biharmonic equation -# ===================================== -# # The first step is to define the variational problem at hand. We define -# the variational problem in UFL terms in a separate form file -# :download:`biharmonic.py`. We begin by defining the finite element:: +# the variational problem in UFL terms in a separate form file. +# We begin by defining the finite element:: + from basix.ufl import element from ufl import ( CellDiameter, diff --git a/cpp/demo/biharmonic/main.cpp b/cpp/demo/biharmonic/main.cpp index 46364eb467f..7a8bf8458a3 100644 --- a/cpp/demo/biharmonic/main.cpp +++ b/cpp/demo/biharmonic/main.cpp @@ -1,5 +1,4 @@ -// Biharmonic equation (C++) -// ========================= +// # Biharmonic equation (C++) // // This demo illustrates how to: // @@ -7,118 +6,133 @@ // * Use a discontinuous Galerkin method // * Solve a fourth-order differential equation // -// The solution for :math:`u` in this demo will look as follows: +// The solution for $u$ in this demo will look as follows: // -// .. image:: ../biharmonic_u.png -// :scale: 75 % +// ```{figure} ../biharmonic_u.png +// :scale: 75 % +// :alt: biharmonic +// solution +// ``` // -// Equation and problem definition -// ------------------------------- +// ## Equation and problem definition // -// The biharmonic equation is a fourth-order elliptic equation. On the -// domain :math:`\Omega \subset \mathbb{R}^{d}`, :math:`1 \le d \le 3`, -// it reads +// ### Strong formulation // -// .. math:: -// \nabla^{4} u = f \quad {\rm in} \ \Omega, +// The biharmonic equation is a fourth-order elliptic equation. +// On the domain $\Omega \subset \mathbb{R}^{d}$, $1 \le d \le 3$, it reads // -// where :math:`\nabla^{4} \equiv \nabla^{2} \nabla^{2}` is the -// biharmonic operator and :math:`f` is a prescribed source term. To -// formulate a complete boundary value problem, the biharmonic equation +// $$ +// \nabla^{4} u = f \quad {\rm in} \ \Omega, +// $$ +// +// where $\nabla^{4} \equiv \nabla^{2} \nabla^{2}$ is the biharmonic operator +// and $f$ is a prescribed source term. +// To formulate a complete boundary value problem, the biharmonic equation // must be complemented by suitable boundary conditions. // +// ### Weak formulation +// // Multiplying the biharmonic equation by a test function and integrating -// by parts twice leads to a problem second-order derivatives, which -// would requires :math:`H^{2}` conforming (roughly :math:`C^{1}` -// continuous) basis functions. To solve the biharmonic equation using -// Lagrange finite element basis functions, the biharmonic equation can -// be split into two second-order equations (see the Mixed Poisson demo -// for a mixed method for the Poisson equation), or a variational -// formulation can be constructed that imposes weak continuity of normal -// derivatives between finite element cells. The demo uses a -// discontinuous Galerkin approach to impose continuity of the normal -// derivative weakly. -// -// Consider a triangulation :math:`\mathcal{T}` of the domain -// :math:`\Omega`, where the set of interior facets is denoted by -// :math:`\mathcal{E}_h^{\rm int}`. Functions evaluated on opposite -// sides of a facet are indicated by the subscripts ':math:`+`' and -// ':math:`-`'. Using the standard continuous Lagrange finite element -// space -// -// .. math:: -// V = \left\{v \in H^{1}_{0}(\Omega)\,:\, v \in P_{k}(K) \ \forall \ K \in -// \mathcal{T} \right\} +// by parts twice leads to a problem of second-order derivatives, which would +// require $H^{2}$ conforming (roughly $C^{1}$ continuous) basis functions. +// To solve the biharmonic equation using Lagrange finite element basis +// functions, the biharmonic equation can be split into two second-order +// equations (see the Mixed Poisson demo for a mixed method for the Poisson +// equation), or a variational formulation can be constructed that imposes +// weak continuity of normal derivatives between finite element cells. +// This demo uses a discontinuous Galerkin approach to impose continuity +// of the normal derivative weakly. +// +// Consider a triangulation $\mathcal{T}$ of the domain $\Omega$, where +// the set of interior facets is denoted by $\mathcal{E}_h^{\rm int}$. +// Functions evaluated on opposite sides of a facet are indicated by the +// subscripts $+$ and $-$. +// Using the standard continuous Lagrange finite element space +// +// $$ +// V = \left\{v \in H^{1}_{0}(\Omega)\,:\, v \in P_{k}(K) \ +// \forall \ K \in \mathcal{T} \right\} +// $$ // // and considering the boundary conditions // -// .. math:: -// u &= 0 \quad {\rm on} \ \partial\Omega \\ -// \nabla^{2} u &= 0 \quad {\rm on} \ \partial\Omega +// \begin{align} +// u &= 0 \quad {\rm on} \ \partial\Omega, \\ +// \nabla^{2} u &= 0 \quad {\rm on} \ \partial\Omega, +// \end{align} // -// a weak formulation of the biharmonic problem reads: find :math:`u \in -// V` such that +// a weak formulation of the biharmonic problem reads: find $u \in V$ such that // -// .. math:: -// a(u,v)=L(v) \quad \forall \ v \in V, +// $$ +// a(u,v)=L(v) \quad \forall \ v \in V, +// $$ // // where the bilinear form is // -// .. math:: -// a(u, v) = \sum_{K \in \mathcal{T}} \int_{K} \nabla^{2} u \nabla^{2} v \, -// {\rm d}x \ -// +\sum_{E \in \mathcal{E}_h^{\rm int}}\left(\int_{E} \frac{\alpha}{h_E} -// [\!\![ \nabla u ]\!\!] [\!\![ \nabla v ]\!\!] \, {\rm d}s -// - \int_{E} \left<\nabla^{2} u \right>[\!\![ \nabla v ]\!\!] \, {\rm d}s -// - \int_{E} [\!\![ \nabla u ]\!\!] \left<\nabla^{2} v \right> \, {\rm -// d}s\right) +// \begin{align*} +// a(u, v) &= +// \sum_{K \in \mathcal{T}} \int_{K} \nabla^{2} u \nabla^{2} v \, {\rm d}x \\ +// &\qquad+\sum_{E \in \mathcal{E}_h^{\rm int}}\left(\int_{E} \frac{\alpha}{h_E} +// [\!\![ \nabla u ]\!\!] [\!\![ \nabla v ]\!\!] \, {\rm d}s +// - \int_{E} \left<\nabla^{2} u \right>[\!\![ \nabla v ]\!\!] \, {\rm d}s +// - \int_{E} [\!\![ \nabla u ]\!\!] \left<\nabla^{2} v \right> \, +// {\rm d}s\right) +// \end{align*} // // and the linear form is // -// .. math:: -// L(v) = \int_{\Omega} fv \, {\rm d}x +// $$ +// L(v) = \int_{\Omega} fv \, {\rm d}x. +// $$ // -// Furthermore, :math:`\left< u \right> = \frac{1}{2} (u_{+} + u_{-})`, -// :math:`[\!\![ w ]\!\!] = w_{+} \cdot n_{+} + w_{-} \cdot n_{-}`, -// :math:`\alpha \ge 0` is a penalty parameter and :math:`h_E` is a measure of -// the cell size. +// Furthermore, $\left< u \right> = \frac{1}{2} (u_{+} + u_{-})$, +// $[\!\![ w ]\!\!] = w_{+} \cdot n_{+} + w_{-} \cdot n_{-}$, +// $\alpha \ge 0$ is a penalty parameter and +// $h_E$ is a measure of the cell size. // // The input parameters for this demo are defined as follows: // -// * :math:`\Omega = [0,1] \times [0,1]` (a unit square) -// * :math:`\alpha = 8.0` (penalty parameter) -// * :math:`f = 4.0 \pi^4\sin(\pi x)\sin(\pi y)` (source term) +// - $\Omega = [0,1] \times [0,1]$ (a unit square) +// - $\alpha = 8.0$ (penalty parameter) +// - $f = 4.0 \pi^4\sin(\pi x)\sin(\pi y)$ (source term) +// // // -// Implementation -// -------------- +// ## Implementation // // The implementation is split in two files: a form file containing the // definition of the variational forms expressed in UFL and a C++ file // containing the actual solver. // -// Running this demo requires the files: :download:`main.cpp`, -// :download:`biharmonic.py` and :download:`CMakeLists.txt`. +// Running this demo requires the files: {download}`demo_biharmonic/main.cpp`, +// {download}`demo_biharmonic/biharmonic.py` and {download}`demo_biharmonic/CMakeLists.txt`. // -// UFL form file -// ^^^^^^^^^^^^^ +// ### UFL form file // -// The UFL file is implemented in :download:`biharmonic.py`, and the -// explanation of the UFL file can be found at :doc:`here -// `. +// The UFL file is implemented in {download}`demo_biharmonic/biharmonic.py`. +// ````{admonition} UFL form implemented in python +// :class: dropdown +// ![ufl-code] +// ```` // +// ````{note} +// TODO: explanation on how to run cmake and/or shell commands for ffcx +// To compile biharmonic.py using FFCx with an option +// for PETSc scalar type `float64` one woud execute the command +// ```bash +// ffcx biharmonic.py --scalar_type=float64 +// ``` +// ```` // -// C++ program -// ^^^^^^^^^^^ +// ### C++ program // -// The main solver is implemented in the :download:`main.cpp` file. +// The main solver is implemented in the {download}`demo_biharmonic/main.cpp` file. // // At the top we include the DOLFINx header file and the generated // header file "biharmonic.h" containing the variational forms for the // Biharmonic equation, which are defined in the UFL form file. For // convenience we also include the DOLFINx namespace. // -// .. code-block:: cpp #include "biharmonic.h" #include @@ -136,12 +150,10 @@ using U = typename dolfinx::scalar_value_type_t; // Inside the ``main`` function, we begin by defining a mesh of the // domain. As the unit square is a very standard domain, we can use a -// built-in mesh provided by the :cpp:class:`UnitSquareMesh` factory. In +// built-in mesh provided by the {cpp:class}`UnitSquareMesh` factory. In // order to create a mesh consisting of 32 x 32 squares with each square // divided into two triangles, and the finite element space (specified // in the form file) defined relative to this mesh, we do as follows -// -// .. code-block:: cpp int main(int argc, char* argv[]) { @@ -149,25 +161,22 @@ int main(int argc, char* argv[]) PetscInitialize(&argc, &argv, nullptr, nullptr); { - // Create mesh +// Create mesh auto part = mesh::create_cell_partitioner(mesh::GhostMode::shared_facet); auto mesh = std::make_shared>( mesh::create_rectangle(MPI_COMM_WORLD, {{{0.0, 0.0}, {1.0, 1.0}}}, {32, 32}, mesh::CellType::triangle, part)); - // A function space object, which is defined in the generated code, - // is created: - // - // .. code-block:: cpp +// A function space object, which is defined in the generated code, +// is created: +// - // Create function space +// Create function space auto V = std::make_shared>( fem::create_functionspace(functionspace_form_biharmonic_a, "u", mesh)); - // The source function ::math:`f` and the penalty term - // ::math:`\alpha` are declared: - // - // .. code-block:: cpp +// The source function $f$ and the penalty term +// $\alpha$ are declared: auto f = std::make_shared>(V); f->interpolate( @@ -184,26 +193,26 @@ int main(int argc, char* argv[]) }); auto alpha = std::make_shared>(8.0); - // Define variational forms +// Define variational forms +// + auto a = std::make_shared>(fem::create_form( *form_biharmonic_a, {V, V}, {}, {{"alpha", alpha}}, {})); auto L = std::make_shared>( fem::create_form(*form_biharmonic_L, {V}, {{"f", f}}, {}, {})); - // Now, the Dirichlet boundary condition (:math:`u = 0`) can be - // created using the class :cpp:class:`DirichletBC`. A - // :cpp:class:`DirichletBC` takes two arguments: the value of the - // boundary condition, and the part of the boundary on which the - // condition applies. In our example, the value of the boundary - // condition (0.0) can represented using a :cpp:class:`Function`, - // and the Dirichlet boundary is defined by the indices of degrees - // of freedom to which the boundary condition applies. The - // definition of the Dirichlet boundary condition then looks as - // follows: - // - // .. code-block:: cpp +// Now, the Dirichlet boundary condition ($u = 0$) can be +// created using the class {cpp:class}`DirichletBC`. A +// {cpp:class}`DirichletBC` takes two arguments: the value of the +// boundary condition, and the part of the boundary on which the +// condition applies. In our example, the value of the boundary +// condition (0.0) can represented using a {cpp:class}`Function`, +// and the Dirichlet boundary is defined by the indices of degrees +// of freedom to which the boundary condition applies. The +// definition of the Dirichlet boundary condition then looks as +// follows: - // Define boundary condition +// Define boundary condition auto facets = mesh::locate_entities_boundary( *mesh, 1, [](auto x) @@ -225,16 +234,14 @@ int main(int argc, char* argv[]) *V->mesh()->topology_mutable(), *V->dofmap(), 1, facets); auto bc = std::make_shared>(0.0, bdofs, V); - // Now, we have specified the variational forms and can consider the - // solution of the variational problem. First, we need to define a - // :cpp:class:`Function` ``u`` to store the solution. (Upon - // initialization, it is simply set to the zero function.) Next, we - // can call the ``solve`` function with the arguments ``a == L``, - // ``u`` and ``bc`` as follows: - // - // .. code-block:: cpp +// Now, we have specified the variational forms and can consider the +// solution of the variational problem. First, we need to define a +// {cpp:class}`Function` ``u`` to store the solution. (Upon +// initialization, it is simply set to the zero function.) Next, we +// can call the ``solve`` function with the arguments ``a == L``, +// ``u`` and ``bc`` as follows: - // Compute solution +// Compute solution fem::Function u(V); auto A = la::petsc::Matrix(fem::petsc::create_matrix(*a), false); la::Vector b(L->function_spaces()[0]->dofmap()->index_map, @@ -266,17 +273,17 @@ int main(int argc, char* argv[]) la::petsc::Vector _b(la::petsc::create_vector_wrap(b), false); lu.solve(_u.vec(), _b.vec()); - // Update ghost values before output +// Update ghost values before output +// u.x()->scatter_fwd(); - // The function ``u`` will be modified during the call to solve. A - // :cpp:class:`Function` can be saved to a file. Here, we output the - // solution to a ``VTK`` file (specified using the suffix ``.pvd``) - // for visualisation in an external program such as Paraview. - // - // .. code-block:: cpp +// The function ``u`` will be modified during the call to solve. A +// {cpp:class}`Function` can be saved to a file. Here, we output the +// solution to a ``VTK`` file (specified using the suffix ``.pvd``) +// for visualisation in an external program such as Paraview. - // Save solution in VTK format +// Save solution in VTK format +// io::VTKFile file(MPI_COMM_WORLD, "u.pvd", "w"); file.write({u}, 0.0); } diff --git a/cpp/demo/custom_kernel/main.cpp b/cpp/demo/custom_kernel/main.cpp index 70822ad086e..77c8ebda697 100644 --- a/cpp/demo/custom_kernel/main.cpp +++ b/cpp/demo/custom_kernel/main.cpp @@ -1,9 +1,8 @@ -// Custom cell kernel assembly (C++) +// # Custom cell kernel assembly (C++) // // This demo shows various methods to define custom cell kernels in C++ and // have them assembled into DOLFINx linear algebra data structures. // -// .. code-block:: cpp #include #include @@ -28,10 +27,8 @@ using mdspan2_t = MDSPAN_IMPL_STANDARD_NAMESPACE::mdspan>; -// .. code-block:: cpp - /// @brief Compute the P1 element mass matrix on the reference cell. -/// @tparam V Scalar type. +/// @tparam T Scalar type. /// @param phi Basis functions. /// @param w Integration weights. /// @return Element reference matrix (row-major storage). @@ -48,7 +45,7 @@ std::array A_ref(mdspand_t phi, std::span w) } /// @brief Compute the P1 RHS vector for f=1 on the reference cell. -/// @tparam V Scalar type. +/// @tparam T Scalar type. /// @param phi Basis functions. /// @param w Integration weights. /// @return RHS reference vector. @@ -121,7 +118,8 @@ double assemble_vector0(std::shared_ptr> V, auto kernel, /// be important for performance for lightweight kernels. /// /// @tparam T Scalar type. -/// @param V Function space. +/// @param g mesh geometry. +/// @param dofmap dofmap. /// @param kernel Element kernel to execute. /// @param cells Cells to execute the kernel over. /// @return Frobenius norm squared of the matrix. @@ -150,7 +148,8 @@ double assemble_matrix1(const mesh::Geometry& g, const fem::DofMap& dofmap, /// be important for performance for lightweight kernels. /// /// @tparam T Scalar type. -/// @param V Function space. +/// @param g mesh geometry. +/// @param dofmap dofmap. /// @param kernel Element kernel to execute. /// @param cells Cells to execute the kernel over. /// @return l2 norm squared of the vector. diff --git a/cpp/demo/hyperelasticity/hyperelasticity.py b/cpp/demo/hyperelasticity/hyperelasticity.py index 9a622900ec8..6e3dd7ff69c 100644 --- a/cpp/demo/hyperelasticity/hyperelasticity.py +++ b/cpp/demo/hyperelasticity/hyperelasticity.py @@ -1,13 +1,9 @@ -# UFL input for hyperleasticity -# ============================= -# -# The first step is to define the variational problem at hand. We define -# the variational problem in UFL terms in a separate form file -# :download:`hyperElasticity.py`. +# The first step is to define the variational problem at hand. # # We are interested in solving for a discrete vector field in three # dimensions, so first we need the appropriate finite element space and # trial and test functions on this space:: + from basix.ufl import element from ufl import ( Coefficient, diff --git a/cpp/demo/hyperelasticity/main.cpp b/cpp/demo/hyperelasticity/main.cpp index dc006fb4d9a..c84cddd4d4e 100644 --- a/cpp/demo/hyperelasticity/main.cpp +++ b/cpp/demo/hyperelasticity/main.cpp @@ -1,3 +1,17 @@ +// # Hyperelasticity (C++) +// + +// ### UFL form file +// +// The UFL file is implemented in {download}`demo_hyperelasticity/hyperelasticity.py`. +// ````{admonition} UFL form implemented in python +// :class: dropdown +// ![ufl-code] +// ```` +// +// ### C++ program +// + #include "hyperelasticity.h" #include #include @@ -125,7 +139,6 @@ int main(int argc, char* argv[]) // mesh, we initialize the (finite element) function space defined by the // generated code. // - // .. code-block:: cpp // Create mesh and define function space auto mesh = std::make_shared>(mesh::create_box( diff --git a/cpp/demo/interpolation-io/main.cpp b/cpp/demo/interpolation-io/main.cpp index 7342de4a483..a4ef748ee41 100644 --- a/cpp/demo/interpolation-io/main.cpp +++ b/cpp/demo/interpolation-io/main.cpp @@ -1,3 +1,5 @@ +// # Interpolation-io (C++) +// // Copyright (C) 2022-2023 Garth N. Wells // // This file is part of DOLFINx (https://www.fenicsproject.org) @@ -88,14 +90,14 @@ void interpolate_nedelec(std::shared_ptr> mesh, // Create a Nedelec finite element Function auto u = std::make_shared>(V); - // Interpolate the vector field - // u = [x[0], x[1]] if x[0] < 0.5 - // [x[0] + 1, x[1]] if x[0] >= 0.5 - // in the Nedelec space. - // - // Note that the x1 component of this field is continuous, and the x0 - // component is discontinuous across x0 = 0.5. This function lies in - // the Nedelec space when there are cell edges aligned to x0 = 0.5. +// Interpolate the vector field +// u = [x[0], x[1]] if x[0] < 0.5 +// [x[0] + 1, x[1]] if x[0] >= 0.5 +// in the Nedelec space. +// +// Note that the x1 component of this field is continuous, and the x0 +// component is discontinuous across x0 = 0.5. This function lies in +// the Nedelec space when there are cell edges aligned to x0 = 0.5. // Find cells with all vertices satisfying (0) x0 <= 0.5 and (1) x0 >= 0.5 auto cells0 @@ -137,11 +139,11 @@ void interpolate_nedelec(std::shared_ptr> mesh, }, cells1); - // Nedelec spaces are not generally supported by visualisation tools. - // Simply evaluating a Nedelec function at cell vertices can - // mis-represent the function. However, we can represented a Nedelec - // function exactly in a discontinuous Lagrange space which we can - // then visualise. We do this here. +// Nedelec spaces are not generally supported by visualisation tools. +// Simply evaluating a Nedelec function at cell vertices can +// mis-represent the function. However, we can represented a Nedelec +// function exactly in a discontinuous Lagrange space which we can +// then visualise. We do this here. // First create a degree 2 vector-valued discontinuous Lagrange space // (which contains the N2 space): @@ -161,11 +163,11 @@ void interpolate_nedelec(std::shared_ptr> mesh, // space: u_l->interpolate(*u); - // Output the discontinuous Lagrange space in VTX format. When - // plotting the x0 component the field will appear discontinuous at x0 - // = 0.5 (jump in the normal component between cells) and the x1 - // component will appear continuous (continuous tangent component - // between cells). +// Output the discontinuous Lagrange space in VTX format. When +// plotting the x0 component the field will appear discontinuous at x0 +// = 0.5 (jump in the normal component between cells) and the x1 +// component will appear continuous (continuous tangent component +// between cells). #ifdef HAS_ADIOS2 io::VTXWriter outfile(mesh->comm(), filename.replace_extension("bp"), {u_l}, "BP4"); @@ -181,12 +183,12 @@ int main(int argc, char* argv[]) dolfinx::init_logging(argc, argv); MPI_Init(&argc, &argv); - // The main body of the function is scoped to ensure that all objects - // that depend on an MPI communicator are destroyed before MPI is - // finalised at the end of this function. +// The main body of the function is scoped to ensure that all objects +// that depend on an MPI communicator are destroyed before MPI is +// finalised at the end of this function. { - // Create meshes. For what comes later in this demo we need to - // ensure that a boundary between cells is located at x0=0.5 +// Create meshes. For what comes later in this demo we need to +// ensure that a boundary between cells is located at x0=0.5 // Create mesh using float for geometry coordinates auto mesh0 @@ -203,16 +205,16 @@ int main(int argc, char* argv[]) mesh::CellType::triangle, mesh::create_cell_partitioner(mesh::GhostMode::none))); - // Interpolate a function in a scalar Lagrange space and output the - // result to file for visualisation using different types +// Interpolate a function in a scalar Lagrange space and output the +// result to file for visualisation using different types interpolate_scalar(mesh0, "u32"); interpolate_scalar(mesh1, "u64"); interpolate_scalar>(mesh0, "u_complex64"); interpolate_scalar>(mesh1, "u_complex128"); - // Interpolate a function in a H(curl) finite element space, and - // then interpolate the H(curl) function in a discontinuous Lagrange - // space for visualisation using different types +// Interpolate a function in a H(curl) finite element space, and +// then interpolate the H(curl) function in a discontinuous Lagrange +// space for visualisation using different types interpolate_nedelec(mesh0, "u_nedelec32"); interpolate_nedelec(mesh1, "u_nedelec64"); interpolate_nedelec>(mesh0, "u_nedelec_complex64"); diff --git a/cpp/demo/interpolation_different_meshes/main.cpp b/cpp/demo/interpolation_different_meshes/main.cpp index 8262b96f1e0..6135e493fab 100644 --- a/cpp/demo/interpolation_different_meshes/main.cpp +++ b/cpp/demo/interpolation_different_meshes/main.cpp @@ -1,3 +1,5 @@ +// # Interpolation different meshes (C++) +// // Copyright (C) 2022 Igor A. Baratta and Massimiliano Leoni // // This file is part of DOLFINx (https://www.fenicsproject.org) diff --git a/cpp/demo/mixed_topology/main.cpp b/cpp/demo/mixed_topology/main.cpp index 87e1445a4b7..709445e76bd 100644 --- a/cpp/demo/mixed_topology/main.cpp +++ b/cpp/demo/mixed_topology/main.cpp @@ -1,3 +1,6 @@ +// # Mixed topology (C++) +// + #include #include #include diff --git a/cpp/demo/poisson/main.cpp b/cpp/demo/poisson/main.cpp index e6be92e8334..e2a13e18128 100644 --- a/cpp/demo/poisson/main.cpp +++ b/cpp/demo/poisson/main.cpp @@ -1,5 +1,4 @@ -// Poisson equation (C++) -// ====================== +// # Poisson equation (C++) // // This demo illustrates how to: // @@ -8,7 +7,7 @@ // * Define Expressions // * Define a FunctionSpace // -// The solution for :math:`u` in this demo will look as follows: +// The solution for $u$ in this demo will look as follows: // // .. image:: ../poisson_u.png // :scale: 75 % @@ -18,44 +17,46 @@ // ------------------------------- // // The Poisson equation is the canonical elliptic partial differential -// equation. For a domain :math:`\Omega \subset \mathbb{R}^n` with -// boundary :math:`\partial \Omega = \Gamma_{D} \cup \Gamma_{N}`, the +// equation. For a domain $\Omega \subset \mathbb{R}^n$ with +// boundary $\partial \Omega = \Gamma_{D} \cup \Gamma_{N}$, the // Poisson equation with particular boundary conditions reads: // -// .. math:: +// \begin{align*} // - \nabla^{2} u &= f \quad {\rm in} \ \Omega, \\ // u &= 0 \quad {\rm on} \ \Gamma_{D}, \\ // \nabla u \cdot n &= g \quad {\rm on} \ \Gamma_{N}. \\ +// \end{align*} // -// Here, :math:`f` and :math:`g` are input data and :math:`n` denotes the +// Here, $f$ and $g$ are input data and $n$ denotes the // outward directed boundary normal. The most standard variational form -// of Poisson equation reads: find :math:`u \in V` such that +// of Poisson equation reads: find $u \in V$ such that // -// .. math:: +// $$ // a(u, v) = L(v) \quad \forall \ v \in V, +// $$ +// where $V$ is a suitable function space and // -// where :math:`V` is a suitable function space and -// -// .. math:: +// \begin{align*} // a(u, v) &= \int_{\Omega} \nabla u \cdot \nabla v \, {\rm d} x, \\ // L(v) &= \int_{\Omega} f v \, {\rm d} x // + \int_{\Gamma_{N}} g v \, {\rm d} s. +// \end{align*} // -// The expression :math:`a(u, v)` is the bilinear form and :math:`L(v)` -// is the linear form. It is assumed that all functions in :math:`V` -// satisfy the Dirichlet boundary conditions (:math:`u = 0 \ {\rm on} \ -// \Gamma_{D}`). +// The expression $a(u, v)$ is the bilinear form and $L(v)$ +// is the linear form. It is assumed that all functions in $V$ +// satisfy the Dirichlet boundary conditions ($u = 0 \ {\rm on} \ +// \Gamma_{D}$). // // In this demo, we shall consider the following definitions of the input // functions, the domain, and the boundaries: // -// * :math:`\Omega = [0,1] \times [0,1]` (a unit square) -// * :math:`\Gamma_{D} = \{(0, y) \cup (1, y) \subset \partial \Omega\}` +// * $\Omega = [0,1] \times [0,1]$ (a unit square) +// * $\Gamma_{D} = \{(0, y) \cup (1, y) \subset \partial \Omega\}$ // (Dirichlet boundary) -// * :math:`\Gamma_{N} = \{(x, 0) \cup (x, 1) \subset \partial \Omega\}` +// * $\Gamma_{N} = \{(x, 0) \cup (x, 1) \subset \partial \Omega\}$ // (Neumann boundary) -// * :math:`g = \sin(5x)` (normal derivative) -// * :math:`f = 10\exp(-((x - 0.5)^2 + (y - 0.5)^2) / 0.02)` (source term) +// * $g = \sin(5x)$ (normal derivative) +// * $f = 10\exp(-((x - 0.5)^2 + (y - 0.5)^2) / 0.02)$ (source term) // // // Implementation @@ -65,27 +66,26 @@ // definition of the variational forms expressed in UFL and a C++ file // containing the actual solver. // -// Running this demo requires the files: :download:`main.cpp`, -// :download:`poisson.py` and :download:`CMakeLists.txt`. -// +// Running this demo requires the files: {download}`demo_poisson/main.cpp`, +// {download}`demo_poisson/poisson.py` and {download}`demo_poisson/CMakeLists.txt`. // -// UFL form file -// ^^^^^^^^^^^^^ // -// The UFL file is implemented in :download:`poisson.py`, and the -// explanation of the UFL file can be found at :doc:`here `. +// ### UFL form file // +// The UFL file is implemented in {download}`demo_poisson/poisson.py`. +// ````{admonition} UFL form implemented in python +// :class: dropdown +// ![ufl-code] +// ```` // -// C++ program -// ^^^^^^^^^^^ +// ### C++ program // -// The main solver is implemented in the :download:`main.cpp` file. +// The main solver is implemented in the {download}`demo_poisson/main.cpp` file. // // At the top we include the DOLFINx header file and the generated header // file "Poisson.h" containing the variational forms for the Poisson // equation. For convenience we also include the DOLFINx namespace. // -// .. code-block:: cpp #include "poisson.h" #include @@ -100,19 +100,17 @@ using T = PetscScalar; using U = typename dolfinx::scalar_value_type_t; // Then follows the definition of the coefficient functions (for -// :math:`f` and :math:`g`), which are derived from the -// :cpp:class:`Expression` class in DOLFINx +// $f$ and $g$), which are derived from the +// {cpp:class}`Expression` class in DOLFINx // -// .. code-block:: cpp // Inside the ``main`` function, we begin by defining a mesh of the // domain. As the unit square is a very standard domain, we can use a -// built-in mesh provided by the :cpp:class:`UnitSquareMesh` factory. In +// built-in mesh provided by the {cpp:class}`UnitSquareMesh` factory. In // order to create a mesh consisting of 32 x 32 squares with each square // divided into two triangles, and the finite element space (specified in // the form file) defined relative to this mesh, we do as follows // -// .. code-block:: cpp int main(int argc, char* argv[]) { @@ -129,13 +127,12 @@ int main(int argc, char* argv[]) auto V = std::make_shared>( fem::create_functionspace(functionspace_form_poisson_a, "u", mesh)); - // Next, we define the variational formulation by initializing the - // bilinear and linear forms (:math:`a`, :math:`L`) using the previously - // defined :cpp:class:`FunctionSpace` ``V``. Then we can create the - // source and boundary flux term (:math:`f`, :math:`g`) and attach these - // to the linear form. - // - // .. code-block:: cpp +// Next, we define the variational formulation by initializing the +// bilinear and linear forms ($a$, $L$) using the previously +// defined {cpp:class}`FunctionSpace` ``V``. Then we can create the +// source and boundary flux term ($f$, $g$) and attach these +// to the linear form. +// // Prepare and set Constants for the bilinear form auto kappa = std::make_shared>(2.0); @@ -148,18 +145,17 @@ int main(int argc, char* argv[]) auto L = std::make_shared>(fem::create_form( *form_poisson_L, {V}, {{"f", f}, {"g", g}}, {}, {})); - // Now, the Dirichlet boundary condition (:math:`u = 0`) can be created - // using the class :cpp:class:`DirichletBC`. A :cpp:class:`DirichletBC` - // takes two arguments: the value of the boundary condition, - // and the part of the boundary on which the condition applies. - // In our example, the value of the boundary condition (0.0) can - // represented using a :cpp:class:`Function`, and the Dirichlet boundary - // is defined by the indices of degrees of freedom to which the boundary - // condition applies. - // The definition of the Dirichlet boundary condition then looks - // as follows: - // - // .. code-block:: cpp +// Now, the Dirichlet boundary condition ($u = 0$) can be created +// using the class {cpp:class}`DirichletBC`. A {cpp:class}`DirichletBC` +// takes two arguments: the value of the boundary condition, +// and the part of the boundary on which the condition applies. +// In our example, the value of the boundary condition (0.0) can +// represented using a {cpp:class}`Function`, and the Dirichlet boundary +// is defined by the indices of degrees of freedom to which the boundary +// condition applies. +// The definition of the Dirichlet boundary condition then looks +// as follows: +// // Define boundary condition @@ -205,14 +201,13 @@ int main(int argc, char* argv[]) return {f, {f.size()}}; }); - // Now, we have specified the variational forms and can consider the - // solution of the variational problem. First, we need to define a - // :cpp:class:`Function` ``u`` to store the solution. (Upon - // initialization, it is simply set to the zero function.) Next, we can - // call the ``solve`` function with the arguments ``a == L``, ``u`` and - // ``bc`` as follows: - // - // .. code-block:: cpp +// Now, we have specified the variational forms and can consider the +// solution of the variational problem. First, we need to define a +// {cpp:class}`Function` ``u`` to store the solution. (Upon +// initialization, it is simply set to the zero function.) Next, we can +// call the ``solve`` function with the arguments ``a == L``, ``u`` and +// ``bc`` as follows: +// auto u = std::make_shared>(V); auto A = la::petsc::Matrix(fem::petsc::create_matrix(*a), false); @@ -248,12 +243,11 @@ int main(int argc, char* argv[]) // Update ghost values before output u->x()->scatter_fwd(); - // The function ``u`` will be modified during the call to solve. A - // :cpp:class:`Function` can be saved to a file. Here, we output the - // solution to a ``VTK`` file (specified using the suffix ``.pvd``) for - // visualisation in an external program such as Paraview. - // - // .. code-block:: cpp +// The function ``u`` will be modified during the call to solve. A +// {cpp:class}`Function` can be saved to a file. Here, we output the +// solution to a ``VTK`` file (specified using the suffix ``.pvd``) for +// visualisation in an external program such as Paraview. +// // Save solution in VTK format io::VTKFile file(MPI_COMM_WORLD, "u.pvd", "w"); diff --git a/cpp/demo/poisson/poisson.py b/cpp/demo/poisson/poisson.py index 141191f793d..f35680950ad 100644 --- a/cpp/demo/poisson/poisson.py +++ b/cpp/demo/poisson/poisson.py @@ -1,9 +1,7 @@ -# UFL input for the Poisson equation -# ================================== -# # The first step is to define the variational problem at hand. We define # the variational problem in UFL terms in a separate form file -# :download:`poisson.py`. We begin by defining the finite element:: +# {download}`demo_poisson/poisson.py`. We begin by defining the finite element:: + from basix.ufl import element from ufl import ( Coefficient, @@ -28,8 +26,8 @@ # polynomials on triangles). # # Next, we use this element to initialize the trial and test functions -# (:math:`u` and :math:`v`) and the coefficient functions (:math:`f` and -# :math:`g`):: +# ($u$ and $v$) and the coefficient functions ($f$ and +# $g$):: coord_element = element("Lagrange", "triangle", 1, shape=(2,)) mesh = Mesh(coord_element) diff --git a/cpp/demo/poisson_matrix_free/main.cpp b/cpp/demo/poisson_matrix_free/main.cpp index c996f701d5a..0f5d2b7c77f 100644 --- a/cpp/demo/poisson_matrix_free/main.cpp +++ b/cpp/demo/poisson_matrix_free/main.cpp @@ -4,27 +4,42 @@ // // SPDX-License-Identifier: LGPL-3.0-or-later -// ===================================== -// Matrix-free Conjugate Gradient solver -// ===================================== +// # Matrix-free Conjugate Gradient solver // // This demo illustrates how to: // * Solve a linear partial differential equation using a matrix-free CG solver // * Create and apply Dirichlet boundary conditions // * Compute errors // -// .. math:: +// \begin{align*} // - \nabla^{2} u &= f \quad {\rm in} \ \Omega, \\ // u &= u_D \quad {\rm on} \ \Gamma_{D} +// \end{align*} // // where -// .. math:: +// \begin{align*} // u_D &= 1 + x^2 + 2y^2, \\ // f = -6 +// \end{align*} // -// .. note:: This demo illustrates the use of a matrix-free Conjugate -// Gradient solver. Many practical problems will also require -// a preconditioner to create an efficient solver. +// ```{note} +// This demo illustrates the use of a matrix-free Conjugate +// Gradient solver. Many practical problems will also require +// a preconditioner to create an efficient solver. +// ``` +// + +// ### UFL form file +// +// The UFL file is implemented in {download}`demo_poisson_matrix_free/poisson.py`. +// ````{admonition} UFL form implemented in python +// :class: dropdown +// ![ufl-code] +// ```` +// +// ### C++ program +// + #include "poisson.h" #include diff --git a/cpp/demo/poisson_matrix_free/poisson.py b/cpp/demo/poisson_matrix_free/poisson.py index 76b10faf375..d4da0537c70 100644 --- a/cpp/demo/poisson_matrix_free/poisson.py +++ b/cpp/demo/poisson_matrix_free/poisson.py @@ -1,5 +1,5 @@ # UFL input for the Matrix-free Poisson Demo -# ================================== + from basix.ufl import element from ufl import ( Coefficient, diff --git a/cpp/doc/Doxyfile b/cpp/doc/Doxyfile index 3e4cf5e5bbf..f092a8d9be9 100644 --- a/cpp/doc/Doxyfile +++ b/cpp/doc/Doxyfile @@ -909,7 +909,6 @@ WARN_LOGFILE = # Note: If this tag is empty the current directory is searched. INPUT = ../dolfinx - # This tag can be used to specify the character encoding of the source files # that doxygen parses. Internally doxygen uses the UTF-8 encoding. Doxygen uses # libiconv (or the iconv built into libc) for the transcoding. See the libiconv @@ -1010,7 +1009,7 @@ EXCLUDE = ../dolfinx/io/pugixml.hpp \ ../dolfinx/io/pugixml.cpp \ ../dolfinx/io/pugiconfig.hpp \ ../dolfinx/common/loguru.hpp \ - ../dolfinx/common/loguru.cpp + ../dolfinx/common/loguru.cpp \ # The EXCLUDE_SYMLINKS tag can be used to select whether or not files or # directories that are symbolic links (a Unix file system feature) are excluded diff --git a/cpp/doc/source/api.rst b/cpp/doc/source/api.rst new file mode 100644 index 00000000000..358fe05b27f --- /dev/null +++ b/cpp/doc/source/api.rst @@ -0,0 +1,32 @@ +============================= +DOLFINx C++ API documentation +============================= + +The is experimental documentation for the C++ API. The full Doxygen +generated documentation is `here `_. + +Installation +============ + +See https://docs.fenicsproject.org/dolfinx/main/python/ for installation +instructions. + +API documentation +================= + +.. toctree:: + :maxdepth: 1 + + common + geometry + graph + fem + io + la + mesh + refinement + +* :ref:`genindex` + + +.. See https://github.com/sphinx-doc/sphinx/issues/10152 diff --git a/cpp/doc/source/conf.py b/cpp/doc/source/conf.py index c75b0cd2e21..187a74457dd 100644 --- a/cpp/doc/source/conf.py +++ b/cpp/doc/source/conf.py @@ -9,18 +9,26 @@ # If extensions (or modules to document with autodoc) are in another directory, # add these directories to sys.path here. If the directory is relative to the # documentation root, use os.path.abspath to make it absolute, like shown here. -# -# import os -# import sys -# sys.path.insert(0, os.path.abspath('.')) + +import os +import sys + +sys.path.insert(0, os.path.abspath(".")) + +import jupytext_process # isort:skip +myst_heading_anchors = 3 + +jupytext_process.process() + # -- Project information ----------------------------------------------------- project = "DOLFINx" copyright = "2022, FEniCS Project" author = "FEniCS Project" +# TODO: automate version tagging? # The full version, including alpha/beta/rc tags release = "0.3.1" @@ -32,6 +40,7 @@ # ones. extensions = [ "sphinx.ext.mathjax", + "myst_parser", "breathe", ] @@ -59,6 +68,10 @@ # html_static_path = ['_static'] html_static_path = [] +myst_enable_extensions = [ + "dollarmath", + "amsmath", +] breathe_projects = {"DOLFINx": "../xml/"} breathe_default_project = "DOLFINx" diff --git a/cpp/doc/source/demo.rst b/cpp/doc/source/demo.rst new file mode 100644 index 00000000000..173cdb8053b --- /dev/null +++ b/cpp/doc/source/demo.rst @@ -0,0 +1,26 @@ +.. DOLFINx C++ demos + +Demos +===== + +These demos illustrate DOLFINx usage. +Starting with :doc:`Biharmonic equation ` is recommended. + + +Introductory +------------------- + +.. note:: + Not yet categorized! + +.. toctree:: + :maxdepth: 1 + + demos/demo_poisson.md + demos/demo_hyperelasticity.md + demos/demo_custom_kernel.md + demos/demo_poisson_matrix_free.md + demos/demo_biharmonic.md + demos/demo_interpolation_different_meshes.md + demos/demo_interpolation-io.md + demos/demo_mixed_topology.md diff --git a/cpp/doc/source/index.rst b/cpp/doc/source/index.rst index 358fe05b27f..c713c3e796c 100644 --- a/cpp/doc/source/index.rst +++ b/cpp/doc/source/index.rst @@ -1,5 +1,5 @@ ============================= -DOLFINx C++ API documentation +DOLFINx C++ documentation ============================= The is experimental documentation for the C++ API. The full Doxygen @@ -15,16 +15,10 @@ API documentation ================= .. toctree:: - :maxdepth: 1 - - common - geometry - graph - fem - io - la - mesh - refinement + :maxdepth: 2 + + demo + api * :ref:`genindex` diff --git a/cpp/doc/source/jupytext_process.py b/cpp/doc/source/jupytext_process.py new file mode 100644 index 00000000000..44a21de9f91 --- /dev/null +++ b/cpp/doc/source/jupytext_process.py @@ -0,0 +1,63 @@ +# Copyright (C) 2017-2023 Garth N. Wells, Jack S. Hale +# +# This file is part of DOLFINx (https://www.fenicsproject.org) +# +# SPDX-License-Identifier: LGPL-3.0-or-later + +import pathlib +import shutil + +import jupytext + + +def process(): + """Convert C++ demos in the Jupytext 'light' format into MyST + flavoured markdown (and ipynb?) using Jupytext. These files can then be + included in Sphinx documentation. + + """ + # Make demo doc directory + demo_doc_dir = pathlib.Path("./demos") + demo_doc_dir.mkdir(parents=True, exist_ok=True) + + # Directories to scan demo code files + demo_dirs = pathlib.Path("../../demo") + + # Iterate over subdirectories containing demos + for demo_subdir in demo_dirs.iterdir(): + if demo_subdir.is_dir(): + fname = pathlib.Path("/demo_" + demo_subdir.name) + demo_doc_subdir = demo_doc_dir / fname.name + demo_doc_subdir.mkdir(parents=True, exist_ok=True) + # Process each demo using jupytext/myst + for demo_file in demo_subdir.glob("main.cpp"): + # Copy demo files into documentation demo directory + shutil.copy(demo_file, demo_doc_subdir) + cpp_demo = jupytext.read(demo_file) + cpp_myst_text = jupytext.writes(cpp_demo, fmt="myst") + + # myst-parser does not process blocks with {code-cell} + cpp_myst_text = cpp_myst_text.replace("{code-cell}", "cpp") + + # Similarly for python file, dump python myst text in cpp myst + for pydemo_file in demo_subdir.glob("*.py"): + shutil.copy(pydemo_file, demo_doc_subdir) + python_demo = jupytext.read(pydemo_file) + python_myst_text = jupytext.writes(python_demo, fmt="myst") + python_myst_text = python_myst_text.replace("{code-cell}", "python") + cpp_myst_text = cpp_myst_text.replace("![ufl-code]", python_myst_text) + + for cmake_file in demo_subdir.glob("CMakeLists.txt"): + shutil.copy(cmake_file, demo_doc_subdir) + + myst_file = (demo_doc_dir / fname.name).with_suffix(".md") + with open(myst_file, "w") as fw: + fw.write(cpp_myst_text) + + # There is a posibility to use jupyter-notebooks with C++/C kernels + # ipynb_file = (demo_doc_dir / fname.name).with_suffix(".ipynb") + # jupytext.write(cpp_demo, ipynb_file, fmt="ipynb") + + +if __name__ == "__main__": + process()