From f3f20f3b803412b667c560de084f5582099a0176 Mon Sep 17 00:00:00 2001 From: jorgensd Date: Wed, 1 Oct 2025 23:02:51 +0200 Subject: [PATCH] Start using interlinking properly. Covered: - Changelog - Fundamentals - Complex-mode --- Changelog.md | 28 ++-- chapter1/complex_mode.ipynb | 77 +++++++--- chapter1/complex_mode.py | 87 +++++++---- chapter1/fundamentals.md | 4 +- chapter1/fundamentals_code.ipynb | 240 +++++++++++++++++++------------ chapter1/fundamentals_code.py | 205 +++++++++++++++++--------- docker/Dockerfile | 2 +- index.ipynb | 3 +- pyproject.toml | 6 +- 9 files changed, 434 insertions(+), 218 deletions(-) diff --git a/Changelog.md b/Changelog.md index 78f16da7..6c4bca6e 100644 --- a/Changelog.md +++ b/Changelog.md @@ -11,13 +11,13 @@ ## v0.9.0 -- `scale` in `apply_lifting` has been renamed to `alpha` +- `scale` in {py:func}`apply_lifting` has been renamed to `alpha` - Use `dolfinx.fem.Function.x.petsc_vec` as opposed to `dolfinx.fem.Function.vector` ## v0.8.0 -- Replace all `ufl.FiniteElement` and `ufl.VectorElement` with the appropriate `basix.ufl.element` -- Replace `dolfinx.fem.FunctionSpace` with `dolfinx.fem.functionspace` +- Replace all `ufl.FiniteElement` and `ufl.VectorElement` with the appropriate {py:func}`basix.ufl.element` +- Replace {py:class}`dolfinx.fem.FunctionSpace` with {py:func}`dolfinx.fem.functionspace` ## v0.7.2 @@ -31,12 +31,13 @@ ## v0.7.0 -- Renamed `dolfinx.graph.create_adjacencylist` to `dolfinx.graph.adjacencylist` -- Renamed `dolfinx.plot.create_vtk_mesh` to `dolfinx.plot.vtk_mesh` -- `dolfinx.geometry.BoundingBoxTree` has been changed to `dolfinx.geometry.bb_tree` +- Renamed `dolfinx.graph.create_adjacencylist` to {py:func}`dolfinx.graph.adjacencylist` +- Renamed `dolfinx.plot.create_vtk_mesh` to {py:func}`dolfinx.plot.vtk_mesh` +- Initialization of {py:class}`dolfinx.geometry.BoundingBoxTree` has been changed to {py:func}`dolfinx.geometry.bb_tree` - `create_mesh` with Meshio has been modified. Note that you now need to pass dtype `np.int32` to the cell_data. -- Update dolfinx petsc API. Now one needs to explicitly import `dolfinx.fem.petsc` and `dolfinx.fem.nls`, as PETSc is no longer a strict requirement. Replace `petsc4py.PETSc.ScalarType` with `dolfinx.default_scalar_type` in demos where we do not use `petsc4py` explicitly. - +- Update dolfinx petsc API. Now one needs to explicitly import {py:mod}`dolfinx.fem.petsc` and {py:mod}`dolfinx.fem.nls`, as PETSc is no longer a strict requirement. + Replace `petsc4py.PETSc.ScalarType` with `dolfinx.default_scalar_type` in demos where we do not use {py:mod}`petsc4py` explicitly. + ## v0.6.0 - Remove `ipygany` and `pythreejs` as plotting backends. Using `panel`. @@ -50,7 +51,8 @@ - Using new GMSH interface in DOLFINx (`dolfinx.io.gmshio`) in all demos using GMSH - Added a section on custom Newton-solvers, see [chapter4/newton-solver]. -- Various minor DOLFINx API updates. `dolfinx.mesh.compute_boundary_facets` -> `dolfinx.mesh.exterior_facet_indices` with slightly different functionality. Use `dolfinx.mesh.MeshTagsMetaClass.find` instead of `mt.indices[mt.values==value]`. +- Various minor DOLFINx API updates. `dolfinx.mesh.compute_boundary_facets` -> {py:func}`dolfinx.mesh.exterior_facet_indices` with slightly different functionality. + Use `dolfinx.mesh.MeshTagsMetaClass.find` instead of `mt.indices[mt.values==value]`. - Various numpy updates, use `np.full_like`. - Change all notebooks to use [jupytext](https://jupytext.readthedocs.io/en/latest/install.html) to automatically sync `.ipynb` with `.py` files. - Add example of how to use `DOLFINx` in complex mode, see [chapter1/complex_mode]. @@ -62,11 +64,13 @@ ## 0.4.0 (05.02.2021) - All `pyvista` plotting has been rewritten to use `ipygany` and `pythreejs` as well as using a cleaner interface. -- `dolfinx.plot.create_vtk_topology` has been renamed to `dolfinx.plot.create_vtk_mesh` and can now be directly used as input to `pyvista.UnstructuredGrid`. +- `dolfinx.plot.create_vtk_topology` has been renamed to `dolfinx.plot.create_vtk_mesh` and can now be directly used as input + to {py:class}`pyvista.UnstructuredGrid`. - `dolfinx.fem.Function.compute_point_values` has been deprecated. Interpolation into a CG-1 is now the way of getting vertex values. -- API updates wrt. DOLFINx. `Form`->`form`, `DirichletBC`->`dirichletbc`. +- Instead of initializing class with {py:class}`Form`, use {py:func}`form`. +- Instead of initializing class with {py:class}`DirichletBC` use {py:func}`dirichletbc`. - Updates on error computations in [Error control: Computing convergence rates](chapter4/convergence). -- Added tutorial on interpolation of `ufl.Expression` in [Deflection of a membrane](chapter1/membrane_code). +- Added tutorial on interpolation of {py:class}`ufl.core.expr.Expr` in [Deflection of a membrane](chapter1/membrane_code). - Added tutorial on how to apply constant-valued Dirichlet conditions in [Deflection of a membrane](chapter1/membrane_code). - Various API changes relating to the import structure of DOLFINx diff --git a/chapter1/complex_mode.ipynb b/chapter1/complex_mode.ipynb index 4a8e5011..29dc7d72 100644 --- a/chapter1/complex_mode.ipynb +++ b/chapter1/complex_mode.ipynb @@ -9,32 +9,58 @@ "\n", "Author: Jørgen S. Dokken\n", "\n", - "Many PDEs, such as the [Helmholtz equation](https://docs.fenicsproject.org/dolfinx/v0.4.1/python/demos/demo_helmholtz.html) require complex-valued fields.\n", + "Many PDEs, such as the [Helmholtz equation](https://docs.fenicsproject.org/dolfinx/main/python/demos/demo_helmholtz.html)\n", + "require complex-valued fields.\n", "\n", "For simplicity, let us consider a Poisson equation of the form:\n", "\n", - "$$-\\Delta u = f \\text{ in } \\Omega,$$\n", - "$$ f = -1 - 2j \\text{ in } \\Omega,$$\n", - "$$ u = u_{exact} \\text{ on } \\partial\\Omega,$$\n", - "$$u_{exact}(x, y) = \\frac{1}{2}x^2 + 1j\\cdot y^2,$$\n", + "$$\n", + "\\begin{align}\n", + "-\\Delta u &= f &&\\text{in } \\Omega,\\\\\n", + "f &= -1 - 2j &&\\text{in } \\Omega,\\\\\n", + "u &= u_{exact} &&\\text{on } \\partial\\Omega,\\\\\n", + "u_{exact}(x, y) &= \\frac{1}{2}x^2 + 1j\\cdot y^2,\n", + "\\end{align}\n", + "$$\n", "\n", - "As in [Solving the Poisson equation](./fundamentals) we want to express our partial differential equation as a weak formulation.\n", + "As in [Solving the Poisson equation](./fundamentals) we want to express our partial differential equation\n", + "as a weak formulation.\n", "\n", - "We start by defining our discrete function space $V_h$, such that $u_h\\in V_h$ and $u_h = \\sum_{i=1}^N c_i \\phi_i(x, y)$ where $\\phi_i$ are **real valued** global basis functions of our space $V_h$, and $c_i \\in \\mathcal{C}$ are the **complex valued** degrees of freedom.\n", + "We start by defining our discrete function space $V_h$, such that $u_h\\in V_h$ and\n", + "$u_h = \\sum_{i=1}^N c_i \\phi_i(x, y)$ where $\\phi_i$ are **real valued** global basis\n", + "functions of our space $V_h$, and $c_i \\in \\mathcal{C}$ are the **complex valued** degrees of freedom.\n", "\n", - "Next, we choose a test function $v\\in \\hat V_h$ where $\\hat V_h\\subset V_h$ such that $v\\vert_{\\partial\\Omega}=0$, as done in the first tutorial.\n", - "We now need to define our inner product space. We choose the $L^2$ inner product spaces, which is a _[sesquilinear](https://en.wikipedia.org/wiki/Sesquilinear_form) 2-form_, meaning that $\\langle u, v\\rangle$ is a map from $V_h\\times V_h\\mapsto K$, and $\\langle u, v \\rangle = \\int_\\Omega u \\cdot \\bar v ~\\mathrm{d} x$. As it is sesquilinear, we have the following properties:\n", + "Next, we choose a test function $v\\in \\hat V_h$ where $\\hat V_h\\subset V_h$ such that $v\\vert_{\\partial\\Omega}=0$, as done in the [first tutorial](./fundamentals).\n", + "We now need to define our inner product space.\n", + "We choose the $L^2$ inner product spaces, which is a _[sesquilinear](https://en.wikipedia.org/wiki/Sesquilinear_form) 2-form_,\n", + "meaning that $\\langle u, v\\rangle$ is a map from $V_h\\times V_h\\mapsto K$, and\n", + "$\\langle u, v \\rangle = \\int_\\Omega u \\cdot \\bar v ~\\mathrm{d} x$. As it is sesquilinear, we have the following properties:\n", "\n", - "$$\\langle u , v \\rangle = \\overline{\\langle v, u \\rangle},$$\n", - "$$\\langle u , u \\rangle \\geq 0.$$\n", + "$$\n", + "\\begin{align}\n", + "\\langle u , v \\rangle &= \\overline{\\langle v, u \\rangle},\\\\\n", + "\\langle u , u \\rangle &\\geq 0.\n", + "\\end{align}\n", + "$$\n", "\n", "We can now use this inner product space to do integration by parts\n", "\n", - "$$\\int_\\Omega \\nabla u_h \\cdot \\nabla \\overline{v}~ \\mathrm{dx} = \\int_{\\Omega} f \\cdot \\overline{v} ~\\mathrm{d} s \\qquad \\forall v \\in \\hat{V}_h.$$\n", + "$$\n", + "\\int_\\Omega \\nabla u_h \\cdot \\nabla \\overline{v}~\\mathrm{dx} =\n", + "\\int_{\\Omega} f \\cdot \\overline{v} ~\\mathrm{d} s \\qquad \\forall v \\in \\hat{V}_h.\n", + "$$\n", "\n", "## Installation of FEniCSx with complex number support\n", "\n", - "FEniCSx supports both real and complex numbers, so we can create a function space with real valued or complex valued coefficients.\n" + "FEniCSx supports both real and complex numbers, so we can create a {py:class}`function space `\n", + "with either real valued or complex valued coefficients.\n", + "```{admonition} Function or Coefficient\n", + "In FEniCSx, the term *function* and *coefficient* are used interchangeably.\n", + "A function is a linear combination of basis functions with coefficients, and the coefficients can be real or complex numbers.\n", + "In {py:mod}`ufl`, the term {py:class}`Coefficient `, while in {py:mod}`dolfinx` we use {py:class}`Function`\n", + "to represent the same concept (through inheritance). This is because most people think of finding the **unknown** function that solves a PDE,\n", + "while the coefficients are the set of values that define the function.\n", + "```" ] }, { @@ -63,9 +89,16 @@ "id": "2", "metadata": {}, "source": [ - "However, as we would like to solve linear algebra problems of the form $Ax=b$, we need to be able to use matrices and vectors that support real and complex numbers. As [PETSc](https://petsc.org/release/) is one of the most popular interfaces to linear algebra packages, we need to be able to work with their matrix and vector structures.\n", + "However, as we would like to solve linear algebra problems of the form $Ax=b$, we need to be able to use matrices and vectors that support real and complex numbers.\n", + "As {[PETSc](https://petsc.org/release/)} is the most popular interfaces to linear algebra packages, we need to be able to work with their matrix and vector structures.\n", "\n", - "Unfortunately, PETSc only supports one floating type in their matrices, thus we need to install two versions of PETSc, one that supports `float64` and one that supports `complex128`. In the [docker images](https://hub.docker.com/r/dolfinx/dolfinx) for DOLFINx, both versions are installed, and one can switch between them by calling `source dolfinx-real-mode` or `source dolfinx-complex-mode`. For the `dolfinx/lab` images, one can change the Python kernel to be either the real or complex mode, by going to `Kernel->Change Kernel...` and choosing `Python3 (ipykernel)` (for real mode) or `Python3 (DOLFINx complex)` (for complex mode).\n", + "Unfortunately, PETSc only supports one floating type in their matrices, thus we need to install two versions of PETSc,\n", + "one that supports `float64` and one that supports `complex128`.\n", + "In the [Docker images]https://github.com/orgs/FEniCS/packages/container/package/dolfinx%2Fdolfinx) for DOLFINx, both versions are installed,\n", + "and one can switch between them by calling `source dolfinx-real-mode` or `source dolfinx-complex-mode`.\n", + "For the [dolfinx/lab](https://github.com/FEniCS/dolfinx/pkgs/container/dolfinx%2Flab) images,\n", + "one can change the Python kernel to be either the real or complex mode, by going to\n", + "`Kernel->Change Kernel...` and choosing `Python3 (ipykernel)` (for real mode) or `Python3 (DOLFINx complex)` (for complex mode).\n", "\n", "We check that we are using the correct installation of PETSc by inspecting the scalar type." ] @@ -114,9 +147,14 @@ "id": "6", "metadata": {}, "source": [ - "Note that we have used the `PETSc.ScalarType` to wrap the constant source on the right hand side. This is because we want the integration kernels to assemble into the correct floating type.\n", + "Note that we have used the `PETSc.ScalarType` to wrap the constant source on the right hand side.\n", + "This is because we want the integration kernels to assemble into the correct floating type.\n", "\n", - "Secondly, note that we are using `ufl.inner` to describe multiplication of $f$ and $v$, even if they are scalar values. This is because `ufl.inner` takes the conjugate of the second argument, as decribed by the $L^2$ inner product. One could alternatively write this out explicitly\n", + "Secondly, note that we are using {py:func}`ufl.inner` to describe multiplication of $f$ and $v$,\n", + "even if they are scalar values.\n", + "This is because {py:func}`ufl.inner` takes the conjugate of the second argument,\n", + "as decribed by the $L^2$ inner product.\n", + "One could alternatively write this out explicitly\n", "\n", "### Inner-products and derivatives" ] @@ -138,7 +176,10 @@ "id": "8", "metadata": {}, "source": [ - "Similarly, if we want to use the function `ufl.derivative` to take derivatives of functionals, we need to take some special care. As `ufl.derivative` inserts a `ufl.TestFunction` to represent the variation, we need to take the conjugate of this to be able to use it to assemble vectors.\n" + "Similarly, if we want to use the function {py:func}`ufl.derivative` to take derivatives of functionals,\n", + "we need to take some special care.\n", + "As {py:func}`ufl.derivative` inserts a {py:func}`ufl.TestFunction` to represent the variation,\n", + "we need to take the conjugate of this to be able to use it to assemble vectors." ] }, { diff --git a/chapter1/complex_mode.py b/chapter1/complex_mode.py index 530fd871..63e2a81a 100644 --- a/chapter1/complex_mode.py +++ b/chapter1/complex_mode.py @@ -17,33 +17,58 @@ # # Author: Jørgen S. Dokken # -# Many PDEs, such as the [Helmholtz equation](https://docs.fenicsproject.org/dolfinx/v0.4.1/python/demos/demo_helmholtz.html) require complex-valued fields. +# Many PDEs, such as the [Helmholtz equation](https://docs.fenicsproject.org/dolfinx/main/python/demos/demo_helmholtz.html) +# require complex-valued fields. # # For simplicity, let us consider a Poisson equation of the form: # -# $$-\Delta u = f \text{ in } \Omega,$$ -# $$ f = -1 - 2j \text{ in } \Omega,$$ -# $$ u = u_{exact} \text{ on } \partial\Omega,$$ -# $$u_{exact}(x, y) = \frac{1}{2}x^2 + 1j\cdot y^2,$$ -# -# As in [Solving the Poisson equation](./fundamentals) we want to express our partial differential equation as a weak formulation. -# -# We start by defining our discrete function space $V_h$, such that $u_h\in V_h$ and $u_h = \sum_{i=1}^N c_i \phi_i(x, y)$ where $\phi_i$ are **real valued** global basis functions of our space $V_h$, and $c_i \in \mathcal{C}$ are the **complex valued** degrees of freedom. -# -# Next, we choose a test function $v\in \hat V_h$ where $\hat V_h\subset V_h$ such that $v\vert_{\partial\Omega}=0$, as done in the first tutorial. -# We now need to define our inner product space. We choose the $L^2$ inner product spaces, which is a _[sesquilinear](https://en.wikipedia.org/wiki/Sesquilinear_form) 2-form_, meaning that $\langle u, v\rangle$ is a map from $V_h\times V_h\mapsto K$, and $\langle u, v \rangle = \int_\Omega u \cdot \bar v ~\mathrm{d} x$. As it is sesquilinear, we have the following properties: -# -# $$\langle u , v \rangle = \overline{\langle v, u \rangle},$$ -# $$\langle u , u \rangle \geq 0.$$ +# $$ +# \begin{align} +# -\Delta u &= f &&\text{in } \Omega,\\ +# f &= -1 - 2j &&\text{in } \Omega,\\ +# u &= u_{exact} &&\text{on } \partial\Omega,\\ +# u_{exact}(x, y) &= \frac{1}{2}x^2 + 1j\cdot y^2, +# \end{align} +# $$ +# +# As in [Solving the Poisson equation](./fundamentals) we want to express our partial differential equation +# as a weak formulation. +# +# We start by defining our discrete function space $V_h$, such that $u_h\in V_h$ and +# $u_h = \sum_{i=1}^N c_i \phi_i(x, y)$ where $\phi_i$ are **real valued** global basis +# functions of our space $V_h$, and $c_i \in \mathcal{C}$ are the **complex valued** degrees of freedom. +# +# Next, we choose a test function $v\in \hat V_h$ where $\hat V_h\subset V_h$ such that $v\vert_{\partial\Omega}=0$, as done in the [first tutorial](./fundamentals). +# We now need to define our inner product space. +# We choose the $L^2$ inner product spaces, which is a _[sesquilinear](https://en.wikipedia.org/wiki/Sesquilinear_form) 2-form_, +# meaning that $\langle u, v\rangle$ is a map from $V_h\times V_h\mapsto K$, and +# $\langle u, v \rangle = \int_\Omega u \cdot \bar v ~\mathrm{d} x$. As it is sesquilinear, we have the following properties: +# +# $$ +# \begin{align} +# \langle u , v \rangle &= \overline{\langle v, u \rangle},\\ +# \langle u , u \rangle &\geq 0. +# \end{align} +# $$ # # We can now use this inner product space to do integration by parts # -# $$\int_\Omega \nabla u_h \cdot \nabla \overline{v}~ \mathrm{dx} = \int_{\Omega} f \cdot \overline{v} ~\mathrm{d} s \qquad \forall v \in \hat{V}_h.$$ +# $$ +# \int_\Omega \nabla u_h \cdot \nabla \overline{v}~\mathrm{dx} = +# \int_{\Omega} f \cdot \overline{v} ~\mathrm{d} s \qquad \forall v \in \hat{V}_h. +# $$ # # ## Installation of FEniCSx with complex number support # -# FEniCSx supports both real and complex numbers, so we can create a function space with real valued or complex valued coefficients. -# +# FEniCSx supports both real and complex numbers, so we can create a {py:class}`function space ` +# with either real valued or complex valued coefficients. +# ```{admonition} Function or Coefficient +# In FEniCSx, the term *function* and *coefficient* are used interchangeably. +# A function is a linear combination of basis functions with coefficients, and the coefficients can be real or complex numbers. +# In {py:mod}`ufl`, the term {py:class}`Coefficient `, while in {py:mod}`dolfinx` we use {py:class}`Function` +# to represent the same concept (through inheritance). This is because most people think of finding the **unknown** function that solves a PDE, +# while the coefficients are the set of values that define the function. +# ``` # + from mpi4py import MPI @@ -60,9 +85,16 @@ print(u_c.x.array.dtype) # - -# However, as we would like to solve linear algebra problems of the form $Ax=b$, we need to be able to use matrices and vectors that support real and complex numbers. As [PETSc](https://petsc.org/release/) is one of the most popular interfaces to linear algebra packages, we need to be able to work with their matrix and vector structures. +# However, as we would like to solve linear algebra problems of the form $Ax=b$, we need to be able to use matrices and vectors that support real and complex numbers. +# As {[PETSc](https://petsc.org/release/)} is the most popular interfaces to linear algebra packages, we need to be able to work with their matrix and vector structures. # -# Unfortunately, PETSc only supports one floating type in their matrices, thus we need to install two versions of PETSc, one that supports `float64` and one that supports `complex128`. In the [docker images](https://hub.docker.com/r/dolfinx/dolfinx) for DOLFINx, both versions are installed, and one can switch between them by calling `source dolfinx-real-mode` or `source dolfinx-complex-mode`. For the `dolfinx/lab` images, one can change the Python kernel to be either the real or complex mode, by going to `Kernel->Change Kernel...` and choosing `Python3 (ipykernel)` (for real mode) or `Python3 (DOLFINx complex)` (for complex mode). +# Unfortunately, PETSc only supports one floating type in their matrices, thus we need to install two versions of PETSc, +# one that supports `float64` and one that supports `complex128`. +# In the [Docker images]https://github.com/orgs/FEniCS/packages/container/package/dolfinx%2Fdolfinx) for DOLFINx, both versions are installed, +# and one can switch between them by calling `source dolfinx-real-mode` or `source dolfinx-complex-mode`. +# For the [dolfinx/lab](https://github.com/FEniCS/dolfinx/pkgs/container/dolfinx%2Flab) images, +# one can change the Python kernel to be either the real or complex mode, by going to +# `Kernel->Change Kernel...` and choosing `Python3 (ipykernel)` (for real mode) or `Python3 (DOLFINx complex)` (for complex mode). # # We check that we are using the correct installation of PETSc by inspecting the scalar type. @@ -87,9 +119,14 @@ L = ufl.inner(f, v) * ufl.dx # - -# Note that we have used the `PETSc.ScalarType` to wrap the constant source on the right hand side. This is because we want the integration kernels to assemble into the correct floating type. +# Note that we have used the `PETSc.ScalarType` to wrap the constant source on the right hand side. +# This is because we want the integration kernels to assemble into the correct floating type. # -# Secondly, note that we are using `ufl.inner` to describe multiplication of $f$ and $v$, even if they are scalar values. This is because `ufl.inner` takes the conjugate of the second argument, as decribed by the $L^2$ inner product. One could alternatively write this out explicitly +# Secondly, note that we are using {py:func}`ufl.inner` to describe multiplication of $f$ and $v$, +# even if they are scalar values. +# This is because {py:func}`ufl.inner` takes the conjugate of the second argument, +# as decribed by the $L^2$ inner product. +# One could alternatively write this out explicitly # # ### Inner-products and derivatives @@ -97,8 +134,10 @@ print(L) print(L2) -# Similarly, if we want to use the function `ufl.derivative` to take derivatives of functionals, we need to take some special care. As `ufl.derivative` inserts a `ufl.TestFunction` to represent the variation, we need to take the conjugate of this to be able to use it to assemble vectors. -# +# Similarly, if we want to use the function {py:func}`ufl.derivative` to take derivatives of functionals, +# we need to take some special care. +# As {py:func}`ufl.derivative` inserts a {py:func}`ufl.TestFunction` to represent the variation, +# we need to take the conjugate of this to be able to use it to assemble vectors. J = u_c**2 * ufl.dx F = ufl.derivative(J, u_c, ufl.conj(v)) diff --git a/chapter1/fundamentals.md b/chapter1/fundamentals.md index 093463fe..86fc737b 100644 --- a/chapter1/fundamentals.md +++ b/chapter1/fundamentals.md @@ -5,7 +5,9 @@ Authors: Hans Petter Langtangen, Anders Logg Adapted to FEniCSx by Jørgen S. Dokken The goal of this tutorial is to solve one of the most basic PDEs, the Poisson equation, with a few lines of code in FEniCSx. -We start by introducing some fundamental FEniCSx objects, such as `Function`, `functionspace`, `TrialFunction` and `TestFunction`, and learn how to write a basic PDE solver. +We start by introducing some fundamental FEniCSx objects, such as {py:class}`Function`, +{py:func}`functionspace`, {py:func}`TrialFunction` and {py:func}`TestFunction`, +and learn how to write a basic PDE solver. This will include: - How to formulate a mathematical variational problem diff --git a/chapter1/fundamentals_code.ipynb b/chapter1/fundamentals_code.ipynb index 467d869d..63b72acd 100644 --- a/chapter1/fundamentals_code.ipynb +++ b/chapter1/fundamentals_code.ipynb @@ -21,22 +21,35 @@ "\n", "## Interactive tutorials\n", "```{admonition} Run the tutorial as Jupyter notebook in browser\n", - "As this book has been published as a Jupyter Book, each code can be run in your browser as a Jupyter notebook. To start such a notebook click the rocket symbol in the top right corner of the relevant tutorial.\n", + "As this book has been published as a Jupyter Book, each code can be run in your browser as a Jupyter notebook.\n", + "To start such a notebook click the rocket symbol in the top right corner of the relevant tutorial.\n", "```\n", "\n", - "The Poisson problem has so far featured a general domain $\\Omega$ and general functions $u_D$ for the boundary conditions and $f$ for the right hand side.\n", - "Therefore, we need to make specific choices of $\\Omega, u_D$ and $f$. A wise choice is to construct a problem with a known analytical solution, so that we can check that the computed solution is correct. The primary candidates are lower-order polynomials. The continuous Galerkin finite element spaces of degree $r$ will exactly reproduce polynomials of degree $r$.\n", + "The Poisson problem has so far featured a general domain $\\Omega$ and general functions $u_D$ for\n", + "the boundary conditions and $f$ for the right hand side.\n", + "Therefore, we need to make specific choices of $\\Omega, u_D$ and $f$.\n", + "A wise choice is to construct a problem with a known analytical solution,\n", + "so that we can check that the computed solution is correct.\n", + "The primary candidates are lower-order polynomials.\n", + "The continuous Galerkin finite element spaces of degree $r$ will exactly reproduce polynomials of degree $r$.\n", "\n", " We use this fact to construct a quadratic function in $2D$. In particular we choose\n", + "\n", + "$$\n", "\\begin{align}\n", " u_e(x,y)=1+x^2+2y^2\n", " \\end{align}\n", + "$$\n", "\n", "Inserting $u_e$ in the original boundary problem, we find that\n", + "\n", + "$$\n", "\\begin{align}\n", " f(x,y)= -6,\\qquad u_D(x,y)=u_e(x,y)=1+x^2+2y^2,\n", "\\end{align}\n", + "$$\n", + "\n", "regardless of the shape of the domain as long as we prescribe\n", "$u_e$ on the boundary.\n", "\n", @@ -44,14 +57,19 @@ "\n", "This simple but very powerful method for constructing test problems is called _the method of manufactured solutions_.\n", "First pick a simple expression for the exact solution, plug into\n", - "the equation to obtain the right-hand side (source term $f$). Then solve the equation with this right hand side, and using the exact solution as boundary condition. Finally, we create a program that tries to reproduce the exact solution.\n", + "the equation to obtain the right-hand side (source term $f$).\n", + "Then solve the equation with this right hand side, and using the exact solution as boundary condition.\n", + "Finally, we create a program that tries to reproduce the exact solution.\n", "\n", - "Note that in many cases, it can be hard to determine if the program works if it produces an error of size $10^{-5}$ on a\n", - "$20 \\times 20$ grid. However, since we are using Sobolev spaces,\n", - "we usually know about the numerical errors _asymptotic properties_. For instance that it is proportional to $h^2$ if $h$ is the size of a cell in the mesh. We can then compare the error on meshes with different $h$-values to see if the asymptotic behavior is correct. This technique will be explained in detail in the chapter [Improving your fenics code](./../chapter4/convergence).\n", + "Note that in many cases, it can be hard to determine if the program works if it produces an error of size\n", + "$10^{-5}$ on a $20 \\times 20$ grid.\n", + "However, since we are using Sobolev spaces, we usually know about the numerical errors _asymptotic properties_.\n", + "For instance that it is proportional to $h^2$ if $h$ is the size of a cell in the mesh.\n", + "We can then compare the error on meshes with different $h$-values to see if the asymptotic behavior is correct.\n", + "This technique will be explained in detail in the chapter [Improving your fenics code](./../chapter4/convergence).\n", "\n", - "However, in cases where we have a solution we know that should have no approximation error, we know that the solution should\n", - "be produced to machine precision by the program." + "However, in cases where we have a solution we know that should have no approximation error,\n", + "we know that the solution should be produced to machine precision by the program." ] }, { @@ -59,9 +77,15 @@ "id": "1", "metadata": {}, "source": [ - "A major difference between a traditional FEniCS code and a FEniCSx code, is that one is not advised to use the wildcard import. We will see this throughout this first example.\n", + "A major difference between a traditional FEniCS code and a FEniCSx code,\n", + "is that one is not advised to use the wildcard import.\n", + "We will see this throughout this first example.\n", + "\n", "## Generating simple meshes\n", - "The next step is to define the discrete domain, _the mesh_. We do this by importing one of the built-in mesh generators. We will build a unit square mesh, i.e. a mesh spanning $[0,1]\\times[0,1]$. It can consist of either triangles or quadrilaterals." + "The next step is to define the discrete domain, _the mesh_.\n", + "We do this by importing one of the built-in mesh generators.\n", + "We will build a {py:func}`unit square mesh`, i.e. a mesh spanning $[0,1]\\times[0,1]$.\n", + "It can consist of either triangles or quadrilaterals." ] }, { @@ -86,12 +110,14 @@ "Note that in addition to give how many elements we would like to have in each direction.\n", "We also have to supply the _MPI-communicator_.\n", "This is to specify how we would like the program to behave in parallel.\n", - "If we supply `MPI.COMM_WORLD` we create a single mesh, whose data is distributed over the number of processors we\n", - "would like to use. We can for instance run the program in parallel on two processors by using `mpirun`, as:\n", + "If we supply {py:data}`MPI.COMM_WORLD` we create a single mesh,\n", + "whose data is distributed over the number of processors we would like to use.\n", + "We can for instance run the program in parallel on two processors by using `mpirun`, as:\n", "``` bash\n", " mpirun -n 2 python3 t1.py\n", "```\n", - "However, if we would like to create a separate mesh on each processor, we can use `MPI.COMM_SELF`.\n", + "However, if we would like to create a separate mesh on each processor,\n", + "we can use {py:data}`MPI.COMM_SELF`.\n", "This is for instance useful if we run a small problem, and would like to run it with multiple parameters.\n", "\n", "## Defining the finite element function space\n", @@ -117,6 +143,14 @@ "V = fem.functionspace(domain, (\"Lagrange\", 1))" ] }, + { + "cell_type": "markdown", + "id": "960048ad", + "metadata": {}, + "source": [ + "Further details about specification/customization of this tuple, see {py:class}`dolfinx.fem.ElementMetaData`." + ] + }, { "cell_type": "markdown", "id": "5", @@ -143,24 +177,14 @@ "id": "7", "metadata": {}, "source": [ - "We now have the boundary data (and in this case the solution of the finite element problem) represented in the discrete function space.\n", - "Next we would like to apply the boundary values to all degrees of freedom that are on the boundary of the discrete domain. We start by identifying the facets (line-segments) representing the outer boundary, using `dolfinx.mesh.exterior_facet_indices`." - ] - }, - { - "cell_type": "markdown", - "id": "8", - "metadata": { - "magic_args": "vscode={\"languageId\": \"python\"}" - }, - "source": [] - }, - { - "cell_type": "markdown", - "id": "9", - "metadata": {}, - "source": [ - "Create facet to cell connectivity required to determine boundary facets" + "We now have the boundary data (and in this case the solution of the finite element problem)\n", + "represented in the discrete function space.\n", + "Next we would like to apply the boundary values to all degrees of freedom that are on the\n", + "boundary of the discrete domain.\n", + "We start by identifying the facets (line-segments) representing the outer boundary,\n", + "using {py:func}`dolfinx.mesh.exterior_facet_indices`.\n", + "We start by creating the facet to cell connectivity required to determine boundary facets by\n", + "calling {py:meth}`dolfinx.mesh.Topology.create_connectivity`." ] }, { @@ -183,12 +207,18 @@ "lines_to_next_cell": 2 }, "source": [ - "For the current problem, as we are using the \"Lagrange\" 1 function space, the degrees of freedom are located at the vertices of each cell, thus each facet contains two degrees of freedom.\n", + "For the current problem, as we are using the first order Lagrange function space,\n", + "the degrees of freedom are located at the vertices of each cell, thus each facet contains two degrees of freedom.\n", "\n", - "To find the local indices of these degrees of freedom, we use `dolfinx.fem.locate_dofs_topological`, which takes in the function space, the dimension of entities in the mesh we would like to identify and the local entities.\n", + "To find the local indices of these degrees of freedom, we use {py:func}`dolfinx.fem.locate_dofs_topological`\n", + "which takes in the function space, the dimension of entities in the mesh we would like to identify and the local entities.\n", "```{admonition} Local ordering of degrees of freedom and mesh vertices\n", "Many people expect there to be a 1-1 correspondence between the mesh coordinates and the coordinates of the degrees of freedom.\n", - "However, this is only true in the case of `Lagrange` 1 elements on a first order mesh. Therefore, in DOLFINx we use separate local numbering for the mesh coordinates and the dof coordinates. To obtain the local dof coordinates we can use `V.tabulate_dof_coordinates()`, while the ordering of the local vertices can be obtained by `mesh.geometry.x`.\n", + "However, this is only true in the case of `Lagrange` 1 elements on a first order mesh.\n", + "Therefore, in DOLFINx we use separate local numbering for the mesh coordinates and the dof coordinates.\n", + "To obtain the local dof coordinates we can use\n", + "{py:meth}`V.tabulate_dof_coordinates()`,\n", + "while the ordering of the local vertices can be obtained by {py:attr}`mesh.geometry.x`.\n", "```\n", "With this data at hand, we can create the Dirichlet boundary condition" ] @@ -213,19 +243,20 @@ "source": [ "## Defining the trial and test function\n", "\n", - "In mathematics, we distinguish between trial and test spaces $V$ and $\\hat{V}$. The only difference in the present problem is the boundary conditions.\n", - "In FEniCSx, we do not specify boundary conditions as part of the function space, so it is sufficient to use a common space for the trial and test function.\n", + "In mathematics, we distinguish between trial and test spaces $V$ and $\\hat{V}$.\n", + "The only difference in the present problem is the boundary conditions.\n", + "In FEniCSx, we do not specify boundary conditions as part of the function space,\n", + "so it is sufficient to use a common space for the trial and test function.\n", "\n", - "We use the [Unified Form Language](https://github.com/FEniCS/ufl/) (UFL) to specify the varitional formulations. See {cite}`ufl2014` for more details." + "We use the {py:mod}`Unified Form Language` (UFL) to specify the varitional formulations.\n", + "See {cite}`ufl2014` for more details." ] }, { "cell_type": "code", "execution_count": null, "id": "14", - "metadata": { - "lines_to_next_cell": 0 - }, + "metadata": {}, "outputs": [], "source": [ "import ufl\n", @@ -240,7 +271,7 @@ "metadata": {}, "source": [ "## Defining the source term\n", - "As the source term is constant over the domain, we use `dolfinx.fem.Constant`" + "As the source term is constant over the domain, we use {py:class}`dolfinx.fem.Constant`" ] }, { @@ -261,11 +292,14 @@ "metadata": {}, "source": [ "```{admonition} Compilation speed-up\n", - "Instead of wrapping $-6$ in a `dolfinx.fem.Constant`, we could simply define $f$ as `f=-6`.\n", - "However, if we would like to change this parameter later in the simulation, we would have to redefine our variational formulation.\n", - "The `dolfinx.fem.Constant` allows us to update the value in $f$ by using `f.value=5`.\n", - "Additionally, by indicating that $f$ is a constant, we speed of compilation of the variational formulations required for the created linear system.\n", + "Instead of wrapping $-6$ in a {py:class}`dolfinx.fem.Constant`, we could simply define $f$ as `f=-6`.\n", + "However, if we would like to change this parameter later in the simulation,\n", + "we would have to redefine our variational formulation.\n", + "The {py:attr}`dolfinx.fem.Constant.value` allows us to update the value in $f$ by using `f.value=5`.\n", + "Additionally, by indicating that $f$ is a constant, we speed of compilation of the variational\n", + "formulations required for the created linear system.\n", "```\n", + "\n", "## Defining the variational problem\n", "As we now have defined all variables used to describe our variational problem, we can create the weak formulation" ] @@ -288,32 +322,42 @@ "source": [ "Note that there is a very close correspondence between the Python syntax and the mathematical syntax\n", "$\\int_{\\Omega} \\nabla u \\cdot \\nabla v ~\\mathrm{d} x$ and $\\int_{\\Omega}fv~\\mathrm{d} x$.\n", - "The integration over the domain $\\Omega$ is defined by using `ufl.dx`, an integration measure over all cells of the mesh.\n", + "The integration over the domain $\\Omega$ is defined by using {py:func}`ufl.dx`, an integration\n", + "{py:class}`measure` over all cells of the mesh.\n", "\n", - "This is the key strength of FEniCSx: the formulas in the variational formulation translate directly to very similar Python code, a feature that makes it easy to specify and solve complicated PDE problems.\n", + "This is the key strength of FEniCSx:\n", + "the formulas in the variational formulation translate directly to very similar Python code,\n", + "a feature that makes it easy to specify and solve complicated PDE problems.\n", "\n", "## Expressing inner products\n", - "The inner product $\\int_\\Omega \\nabla u \\cdot \\nabla v ~\\mathrm{d} x$ can be expressed in various ways in UFL. We have used the notation `ufl.dot(ufl.grad(u), ufl.grad(v))*ufl.dx`.\n", - "The dot product in UFL computes the sum (contraction) over the last index of the first factor and first index of the second factor.\n", - "In this case, both factors are tensors of rank one (vectors) and so the sum is just over the single index of both $\\nabla u$ and $\\nabla v$.\n", - "To compute an inner product of matrices (with two indices), on must use use the function `ufl.inner` instead of `ufl.dot`.\n", - "For vectors, `ufl.dot` and `ufl.inner` are equivalent.\n", + "The inner product $\\int_\\Omega \\nabla u \\cdot \\nabla v ~\\mathrm{d} x$ can be expressed in various ways in UFL.\n", + "We have used the notation `ufl.dot(ufl.grad(u), ufl.grad(v))*ufl.dx`.\n", + "The {py:func}`dot` product in UFL computes the sum (contraction) over the last index\n", + "of the first factor and first index of the second factor.\n", + "In this case, both factors are tensors of rank one (vectors) and so the sum is just over\n", + "the single index of both $\\nabla u$ and $\\nabla v$.\n", + "To compute an inner product of matrices (with two indices),\n", + "one must use the function {py:func}`ufl.inner` instead of {py:func}`ufl.dot`.\n", + "For real-valued vectors, {py:func}`ufl.dot` and {py:func}`ufl.inner` are equivalent.\n", "\n", "```{admonition} Complex numbers\n", "In DOLFINx, one can solve complex number problems by using an installation of PETSc using complex numbers.\n", - "For variational formulations with complex numbers, one cannot use `ufl.dot` to compute inner products.\n", - "One has to use `ufl.inner`, with the test-function as the second input argument for `ufl.inner`.\n", + "For variational formulations with complex numbers, one cannot use {py:func}`ufl.dot` to compute inner products.\n", + "One has to use {py:func}`ufl.inner`, with the test-function as the second input argument for {py:func}`ufl.inner`.\n", "See [Running DOLFINx in complex mode](./complex_mode) for more information.\n", "```\n", "\n", "\n", "## Forming and solving the linear system\n", "\n", - "Having defined the finite element variational problem and boundary condition, we can create our `dolfinx.fem.petsc.LinearProblem`, as class for solving\n", - "the variational problem: Find $u_h\\in V$ such that $a(u_h, v)==L(v) \\quad \\forall v \\in \\hat{V}$. We will use PETSc as our linear algebra backend, using a direct solver (LU-factorization).\n", + "Having defined the finite element variational problem and boundary condition,\n", + "we can create our {py:class}`LinearProblem` to solve the variational problem:\n", + "Find $u_h\\in V$ such that $a(u_h, v)==L(v) \\quad \\forall v \\in \\hat{V}$.\n", + "We will use {py:mod}`PETSc` as our linear algebra backend, using a direct solver (LU-factorization).\n", "See the [PETSc-documentation](https://petsc.org/main/docs/manual/ksp/?highlight=ksp#ksp-linear-system-solvers) of the method for more information.\n", "PETSc is not a required dependency of DOLFINx, and therefore we explicitly import the DOLFINx wrapper for interfacing with PETSc.\n", - "To ensure that the options passed to the LinearProblem is only used for the given KSP solver, we pass a **unique** option prefix as well." + "To ensure that the options passed to the {py:class}`LinearProblem`\n", + "is only used for the given KSP solver, we pass a **unique** option prefix as well." ] }, { @@ -323,16 +367,8 @@ "metadata": {}, "outputs": [], "source": [ - "from dolfinx.fem.petsc import LinearProblem" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "id": "21", - "metadata": {}, - "outputs": [], - "source": [ + "from dolfinx.fem.petsc import LinearProblem\n", + "\n", "problem = LinearProblem(\n", " a,\n", " L,\n", @@ -345,36 +381,43 @@ }, { "cell_type": "markdown", - "id": "22", + "id": "c82e1ec7", "metadata": {}, "source": [ - "Using `problem.solve()` we solve the linear system of equations and return a `dolfinx.fem.Function` containing the solution.\n", + "Using {py:meth}`problem.solve()` we solve the linear system of equations and\n", + "return a {py:class}`Function` containing the solution.\n", + "\n", "## Computing the error\n", - "Finally, we want to compute the error to check the accuracy of the solution. We do this by comparing the finite element solution `u` with the exact solution.\n", - "We do this by interpolating the exact solution into the the $P_2$-function space." + "Finally, we want to compute the error to check the accuracy of the solution.\n", + "We do this by comparing the finite element solution `u` with the exact solution.\n", + "First we interpolate the exact solution into a function space that contains it" ] }, { "cell_type": "code", "execution_count": null, - "id": "23", + "id": "21", "metadata": {}, "outputs": [], "source": [ "V2 = fem.functionspace(domain, (\"Lagrange\", 2))\n", - "uex = fem.Function(V2)\n", + "uex = fem.Function(V2, name=\"u_exact\")\n", "uex.interpolate(lambda x: 1 + x[0] ** 2 + 2 * x[1] ** 2)" ] }, { "cell_type": "markdown", - "id": "24", + "id": "22", "metadata": {}, "source": [ - "We compute the error in two different ways. First, we compute the $L^2$-norm of the error, defined by $E=\\sqrt{\\int_\\Omega (u_D-u_h)^2\\mathrm{d} x}$.\n", - "We use UFL to express the $L^2$-error, and use `dolfinx.fem.assemble_scalar` to compute the scalar value.\n", - "In DOLFINx, `assemble_scalar` only assembles over the cells on the local process. This means that if we use 2 processes to solve our problem, we need to gather the solution to one (or all the processes.\n", - "We can do this with the `MPI.allreduce` function." + "We compute the error in two different ways.\n", + "First, we compute the $L^2$-norm of the error, defined by $E=\\sqrt{\\int_\\Omega (u_D-u_h)^2\\mathrm{d} x}$.\n", + "We use UFL to express the $L^2$-error, and use {py:func}`dolfinx.fem.assemble_scalar` to compute the scalar value.\n", + "In DOLFINx, {py:func}`assemble_scalar`\n", + "only assembles over the cells on the local process.\n", + "This means that if we use 2 processes to solve our problem,\n", + "we need to accumulate the local contributions to get the global error (on one or all processes).\n", + "We can do this with the {func}`Comm.allreduce` function." ] }, { @@ -397,9 +440,13 @@ "Secondly, we compute the maximum error at any degree of freedom.\n", "As the finite element function $u$ can be expressed as a linear combination of basis functions $\\phi_j$, spanning the space $V$:\n", "$ u = \\sum_{j=1}^N U_j\\phi_j.$\n", - "By writing `problem.solve()` we compute all the coefficients $U_1,\\dots, U_N$. These values are known as the _degrees of freedom_ (dofs). We can access the degrees of freedom by accessing the underlying vector in `uh`.\n", - "However, as a second order function space has more dofs than a linear function space, we cannot compare these arrays directly.\n", - "As we already have interpolated the exact solution into the first order space when creating the boundary condition, we can compare the maximum values at any degree of freedom of the approximation space." + "By writing `problem.solve()` we compute all the coefficients $U_1,\\dots, U_N$.\n", + "These values are known as the _degrees of freedom_ (dofs).\n", + "We can access the degrees of freedom by accessing the underlying vector in `uh`.\n", + "However, as a second order function space has more dofs than a linear function space,\n", + "we cannot compare these arrays directly.\n", + "As we already have interpolated the exact solution into the first order space when creating the boundary condition,\n", + "we can compare the maximum values at any degree of freedom of the approximation space." ] }, { @@ -423,9 +470,11 @@ "source": [ "## Plotting the mesh using pyvista\n", "We will visualizing the mesh using [pyvista](https://docs.pyvista.org/), an interface to the VTK toolkit.\n", - "We start by converting the mesh to a format that can be used with `pyvista`.\n", - "To do this we use the function `dolfinx.plot.vtk_mesh`. The first step is to create an unstructured grid that can be used by `pyvista`.\n", - "We need to start a virtual framebuffer for plotting through docker containers. You can print the current backend and change it with `pyvista.set_jupyter_backend(backend)`" + "We start by converting the mesh to a format that can be used with {py:mod}`pyvista`.\n", + "To do this we use the function {py:func}`dolfinx.plot.vtk_mesh`.\n", + "It creates the data required to create a {py:class}`pyvista.UnstructuredGrid`.\n", + "On Linux (in docker) we need to start a virtual framebuffer for plotting through docker containers.\n", + "You can print the current backend and change it with {py:func}`pyvista.set_jupyter_backend`." ] }, { @@ -461,8 +510,8 @@ "metadata": {}, "source": [ "There are several backends that can be used with pyvista, and they have different benefits and drawbacks.\n", - "See the [pyvista documentation](https://docs.pyvista.org/user-guide/jupyter/index.html#state-of-3d-interactive-jupyterlab-plotting) for more information and installation details.\n", - "n this example and the rest of the tutorial we will use [panel](https://github.com/holoviz/panel)." + "See the [pyvista documentation](https://docs.pyvista.org/user-guide/jupyter/index.html#state-of-3d-interactive-jupyterlab-plotting)\n", + "for more information and installation details." ] }, { @@ -470,8 +519,9 @@ "id": "32", "metadata": {}, "source": [ - "We can now use the `pyvista.Plotter` to visualize the mesh. We visualize it by showing it in 2D and warped in 3D.\n", - "In the jupyter notebook environment, we use the default setting of `pyvista.OFF_SCREEN=False`, which will render plots directly in the notebook." + "We can now use the {py:class}`pyvista.Plotter` to visualize the mesh. We visualize it by showing it in 2D and warped in 3D.\n", + "In the jupyter notebook environment, we use the default setting of `pyvista.OFF_SCREEN=False`,\n", + "which will render plots directly in the notebook." ] }, { @@ -496,7 +546,11 @@ "metadata": {}, "source": [ "## Plotting a function using pyvista\n", - "We want to plot the solution `uh`. As the function space used to defined the mesh is disconnected from the function space defining the mesh, we create a mesh based on the dof coordinates for the function space `V`. We use `dolfinx.plot.vtk_mesh` with the function space as input to create a mesh with mesh geometry based on the dof coordinates." + "We want to plot the solution `uh`.\n", + "As the function space used to defined the mesh is decoupled from the representation of the mesh,\n", + "we create a mesh based on the dof coordinates for the function space `V`.\n", + "We use {py:func}`dolfinx.plot.vtk_mesh` with the function space as input to create a mesh with\n", + "mesh geometry based on the dof coordinates." ] }, { @@ -514,7 +568,7 @@ "id": "36", "metadata": {}, "source": [ - "Next, we create the `pyvista.UnstructuredGrid` and add the dof-values to the mesh." + "Next, we create the {py:class}`pyvista.UnstructuredGrid` and add the dof-values to the mesh." ] }, { @@ -562,7 +616,9 @@ "metadata": {}, "source": [ "## External post-processing\n", - "For post-processing outside the python code, it is suggested to save the solution to file using either `dolfinx.io.VTXWriter` or `dolfinx.io.XDMFFile` and using [Paraview](https://www.paraview.org/). This is especially suggested for 3D visualization." + "For post-processing outside the python code, it is suggested to save the solution to file using either\n", + "{py:class}`dolfinx.io.VTXWriter` or {py:class}`dolfinx.io.XDMFFile` and using [Paraview](https://www.paraview.org/).\n", + "This is especially suggested for 3D visualization." ] }, { diff --git a/chapter1/fundamentals_code.py b/chapter1/fundamentals_code.py index 975647ad..94a29c56 100644 --- a/chapter1/fundamentals_code.py +++ b/chapter1/fundamentals_code.py @@ -29,22 +29,35 @@ # # ## Interactive tutorials # ```{admonition} Run the tutorial as Jupyter notebook in browser -# As this book has been published as a Jupyter Book, each code can be run in your browser as a Jupyter notebook. To start such a notebook click the rocket symbol in the top right corner of the relevant tutorial. +# As this book has been published as a Jupyter Book, each code can be run in your browser as a Jupyter notebook. +# To start such a notebook click the rocket symbol in the top right corner of the relevant tutorial. # ``` # -# The Poisson problem has so far featured a general domain $\Omega$ and general functions $u_D$ for the boundary conditions and $f$ for the right hand side. -# Therefore, we need to make specific choices of $\Omega, u_D$ and $f$. A wise choice is to construct a problem with a known analytical solution, so that we can check that the computed solution is correct. The primary candidates are lower-order polynomials. The continuous Galerkin finite element spaces of degree $r$ will exactly reproduce polynomials of degree $r$. +# The Poisson problem has so far featured a general domain $\Omega$ and general functions $u_D$ for +# the boundary conditions and $f$ for the right hand side. +# Therefore, we need to make specific choices of $\Omega, u_D$ and $f$. +# A wise choice is to construct a problem with a known analytical solution, +# so that we can check that the computed solution is correct. +# The primary candidates are lower-order polynomials. +# The continuous Galerkin finite element spaces of degree $r$ will exactly reproduce polynomials of degree $r$. # # We use this fact to construct a quadratic function in $2D$. In particular we choose +# +# $$ # \begin{align} # u_e(x,y)=1+x^2+2y^2 # \end{align} +# $$ # # Inserting $u_e$ in the original boundary problem, we find that +# +# $$ # \begin{align} # f(x,y)= -6,\qquad u_D(x,y)=u_e(x,y)=1+x^2+2y^2, # \end{align} +# $$ +# # regardless of the shape of the domain as long as we prescribe # $u_e$ on the boundary. # @@ -52,18 +65,29 @@ # # This simple but very powerful method for constructing test problems is called _the method of manufactured solutions_. # First pick a simple expression for the exact solution, plug into -# the equation to obtain the right-hand side (source term $f$). Then solve the equation with this right hand side, and using the exact solution as boundary condition. Finally, we create a program that tries to reproduce the exact solution. +# the equation to obtain the right-hand side (source term $f$). +# Then solve the equation with this right hand side, and using the exact solution as boundary condition. +# Finally, we create a program that tries to reproduce the exact solution. # -# Note that in many cases, it can be hard to determine if the program works if it produces an error of size $10^{-5}$ on a -# $20 \times 20$ grid. However, since we are using Sobolev spaces, -# we usually know about the numerical errors _asymptotic properties_. For instance that it is proportional to $h^2$ if $h$ is the size of a cell in the mesh. We can then compare the error on meshes with different $h$-values to see if the asymptotic behavior is correct. This technique will be explained in detail in the chapter [Improving your fenics code](./../chapter4/convergence). +# Note that in many cases, it can be hard to determine if the program works if it produces an error of size +# $10^{-5}$ on a $20 \times 20$ grid. +# However, since we are using Sobolev spaces, we usually know about the numerical errors _asymptotic properties_. +# For instance that it is proportional to $h^2$ if $h$ is the size of a cell in the mesh. +# We can then compare the error on meshes with different $h$-values to see if the asymptotic behavior is correct. +# This technique will be explained in detail in the chapter [Improving your fenics code](./../chapter4/convergence). # -# However, in cases where we have a solution we know that should have no approximation error, we know that the solution should -# be produced to machine precision by the program. +# However, in cases where we have a solution we know that should have no approximation error, +# we know that the solution should be produced to machine precision by the program. -# A major difference between a traditional FEniCS code and a FEniCSx code, is that one is not advised to use the wildcard import. We will see this throughout this first example. +# A major difference between a traditional FEniCS code and a FEniCSx code, +# is that one is not advised to use the wildcard import. +# We will see this throughout this first example. +# # ## Generating simple meshes -# The next step is to define the discrete domain, _the mesh_. We do this by importing one of the built-in mesh generators. We will build a unit square mesh, i.e. a mesh spanning $[0,1]\times[0,1]$. It can consist of either triangles or quadrilaterals. +# The next step is to define the discrete domain, _the mesh_. +# We do this by importing one of the built-in mesh generators. +# We will build a {py:func}`unit square mesh`, i.e. a mesh spanning $[0,1]\times[0,1]$. +# It can consist of either triangles or quadrilaterals. # + from mpi4py import MPI @@ -76,12 +100,14 @@ # Note that in addition to give how many elements we would like to have in each direction. # We also have to supply the _MPI-communicator_. # This is to specify how we would like the program to behave in parallel. -# If we supply `MPI.COMM_WORLD` we create a single mesh, whose data is distributed over the number of processors we -# would like to use. We can for instance run the program in parallel on two processors by using `mpirun`, as: +# If we supply {py:data}`MPI.COMM_WORLD` we create a single mesh, +# whose data is distributed over the number of processors we would like to use. +# We can for instance run the program in parallel on two processors by using `mpirun`, as: # ``` bash # mpirun -n 2 python3 t1.py # ``` -# However, if we would like to create a separate mesh on each processor, we can use `MPI.COMM_SELF`. +# However, if we would like to create a separate mesh on each processor, +# we can use {py:data}`MPI.COMM_SELF`. # This is for instance useful if we run a small problem, and would like to run it with multiple parameters. # # ## Defining the finite element function space @@ -100,6 +126,8 @@ V = fem.functionspace(domain, ("Lagrange", 1)) # - +# Further details about specification/customization of this tuple, see {py:class}`dolfinx.fem.ElementMetaData`. + # ## Dirichlet boundary conditions # Next, we create a function that will hold the Dirichlet boundary data, and use interpolation to # fill it with the appropriate data. @@ -107,26 +135,32 @@ uD = fem.Function(V) uD.interpolate(lambda x: 1 + x[0] ** 2 + 2 * x[1] ** 2) -# We now have the boundary data (and in this case the solution of the finite element problem) represented in the discrete function space. -# Next we would like to apply the boundary values to all degrees of freedom that are on the boundary of the discrete domain. We start by identifying the facets (line-segments) representing the outer boundary, using `dolfinx.mesh.exterior_facet_indices`. - -# + [markdown] magic_args="vscode={\"languageId\": \"python\"}" -# -# - - -# Create facet to cell connectivity required to determine boundary facets +# We now have the boundary data (and in this case the solution of the finite element problem) +# represented in the discrete function space. +# Next we would like to apply the boundary values to all degrees of freedom that are on the +# boundary of the discrete domain. +# We start by identifying the facets (line-segments) representing the outer boundary, +# using {py:func}`dolfinx.mesh.exterior_facet_indices`. +# We start by creating the facet to cell connectivity required to determine boundary facets by +# calling {py:meth}`dolfinx.mesh.Topology.create_connectivity`. tdim = domain.topology.dim fdim = tdim - 1 domain.topology.create_connectivity(fdim, tdim) boundary_facets = mesh.exterior_facet_indices(domain.topology) -# For the current problem, as we are using the "Lagrange" 1 function space, the degrees of freedom are located at the vertices of each cell, thus each facet contains two degrees of freedom. +# For the current problem, as we are using the first order Lagrange function space, +# the degrees of freedom are located at the vertices of each cell, thus each facet contains two degrees of freedom. # -# To find the local indices of these degrees of freedom, we use `dolfinx.fem.locate_dofs_topological`, which takes in the function space, the dimension of entities in the mesh we would like to identify and the local entities. +# To find the local indices of these degrees of freedom, we use {py:func}`dolfinx.fem.locate_dofs_topological` +# which takes in the function space, the dimension of entities in the mesh we would like to identify and the local entities. # ```{admonition} Local ordering of degrees of freedom and mesh vertices # Many people expect there to be a 1-1 correspondence between the mesh coordinates and the coordinates of the degrees of freedom. -# However, this is only true in the case of `Lagrange` 1 elements on a first order mesh. Therefore, in DOLFINx we use separate local numbering for the mesh coordinates and the dof coordinates. To obtain the local dof coordinates we can use `V.tabulate_dof_coordinates()`, while the ordering of the local vertices can be obtained by `mesh.geometry.x`. +# However, this is only true in the case of `Lagrange` 1 elements on a first order mesh. +# Therefore, in DOLFINx we use separate local numbering for the mesh coordinates and the dof coordinates. +# To obtain the local dof coordinates we can use +# {py:meth}`V.tabulate_dof_coordinates()`, +# while the ordering of the local vertices can be obtained by {py:attr}`mesh.geometry.x`. # ``` # With this data at hand, we can create the Dirichlet boundary condition @@ -136,10 +170,13 @@ # ## Defining the trial and test function # -# In mathematics, we distinguish between trial and test spaces $V$ and $\hat{V}$. The only difference in the present problem is the boundary conditions. -# In FEniCSx, we do not specify boundary conditions as part of the function space, so it is sufficient to use a common space for the trial and test function. +# In mathematics, we distinguish between trial and test spaces $V$ and $\hat{V}$. +# The only difference in the present problem is the boundary conditions. +# In FEniCSx, we do not specify boundary conditions as part of the function space, +# so it is sufficient to use a common space for the trial and test function. # -# We use the [Unified Form Language](https://github.com/FEniCS/ufl/) (UFL) to specify the varitional formulations. See {cite}`ufl2014` for more details. +# We use the {py:mod}`Unified Form Language` (UFL) to specify the varitional formulations. +# See {cite}`ufl2014` for more details. # + @@ -148,8 +185,9 @@ u = ufl.TrialFunction(V) v = ufl.TestFunction(V) # - + # ## Defining the source term -# As the source term is constant over the domain, we use `dolfinx.fem.Constant` +# As the source term is constant over the domain, we use {py:class}`dolfinx.fem.Constant` # + from dolfinx import default_scalar_type @@ -158,11 +196,14 @@ # - # ```{admonition} Compilation speed-up -# Instead of wrapping $-6$ in a `dolfinx.fem.Constant`, we could simply define $f$ as `f=-6`. -# However, if we would like to change this parameter later in the simulation, we would have to redefine our variational formulation. -# The `dolfinx.fem.Constant` allows us to update the value in $f$ by using `f.value=5`. -# Additionally, by indicating that $f$ is a constant, we speed of compilation of the variational formulations required for the created linear system. +# Instead of wrapping $-6$ in a {py:class}`dolfinx.fem.Constant`, we could simply define $f$ as `f=-6`. +# However, if we would like to change this parameter later in the simulation, +# we would have to redefine our variational formulation. +# The {py:attr}`dolfinx.fem.Constant.value` allows us to update the value in $f$ by using `f.value=5`. +# Additionally, by indicating that $f$ is a constant, we speed of compilation of the variational +# formulations required for the created linear system. # ``` +# # ## Defining the variational problem # As we now have defined all variables used to describe our variational problem, we can create the weak formulation @@ -171,33 +212,44 @@ # Note that there is a very close correspondence between the Python syntax and the mathematical syntax # $\int_{\Omega} \nabla u \cdot \nabla v ~\mathrm{d} x$ and $\int_{\Omega}fv~\mathrm{d} x$. -# The integration over the domain $\Omega$ is defined by using `ufl.dx`, an integration measure over all cells of the mesh. +# The integration over the domain $\Omega$ is defined by using {py:func}`ufl.dx`, an integration +# {py:class}`measure` over all cells of the mesh. # -# This is the key strength of FEniCSx: the formulas in the variational formulation translate directly to very similar Python code, a feature that makes it easy to specify and solve complicated PDE problems. +# This is the key strength of FEniCSx: +# the formulas in the variational formulation translate directly to very similar Python code, +# a feature that makes it easy to specify and solve complicated PDE problems. # # ## Expressing inner products -# The inner product $\int_\Omega \nabla u \cdot \nabla v ~\mathrm{d} x$ can be expressed in various ways in UFL. We have used the notation `ufl.dot(ufl.grad(u), ufl.grad(v))*ufl.dx`. -# The dot product in UFL computes the sum (contraction) over the last index of the first factor and first index of the second factor. -# In this case, both factors are tensors of rank one (vectors) and so the sum is just over the single index of both $\nabla u$ and $\nabla v$. -# To compute an inner product of matrices (with two indices), on must use use the function `ufl.inner` instead of `ufl.dot`. -# For vectors, `ufl.dot` and `ufl.inner` are equivalent. +# The inner product $\int_\Omega \nabla u \cdot \nabla v ~\mathrm{d} x$ can be expressed in various ways in UFL. +# We have used the notation `ufl.dot(ufl.grad(u), ufl.grad(v))*ufl.dx`. +# The {py:func}`dot` product in UFL computes the sum (contraction) over the last index +# of the first factor and first index of the second factor. +# In this case, both factors are tensors of rank one (vectors) and so the sum is just over +# the single index of both $\nabla u$ and $\nabla v$. +# To compute an inner product of matrices (with two indices), +# one must use the function {py:func}`ufl.inner` instead of {py:func}`ufl.dot`. +# For real-valued vectors, {py:func}`ufl.dot` and {py:func}`ufl.inner` are equivalent. # # ```{admonition} Complex numbers # In DOLFINx, one can solve complex number problems by using an installation of PETSc using complex numbers. -# For variational formulations with complex numbers, one cannot use `ufl.dot` to compute inner products. -# One has to use `ufl.inner`, with the test-function as the second input argument for `ufl.inner`. +# For variational formulations with complex numbers, one cannot use {py:func}`ufl.dot` to compute inner products. +# One has to use {py:func}`ufl.inner`, with the test-function as the second input argument for {py:func}`ufl.inner`. # See [Running DOLFINx in complex mode](./complex_mode) for more information. # ``` # # # ## Forming and solving the linear system # -# Having defined the finite element variational problem and boundary condition, we can create our `dolfinx.fem.petsc.LinearProblem`, as class for solving -# the variational problem: Find $u_h\in V$ such that $a(u_h, v)==L(v) \quad \forall v \in \hat{V}$. We will use PETSc as our linear algebra backend, using a direct solver (LU-factorization). +# Having defined the finite element variational problem and boundary condition, +# we can create our {py:class}`LinearProblem` to solve the variational problem: +# Find $u_h\in V$ such that $a(u_h, v)==L(v) \quad \forall v \in \hat{V}$. +# We will use {py:mod}`PETSc` as our linear algebra backend, using a direct solver (LU-factorization). # See the [PETSc-documentation](https://petsc.org/main/docs/manual/ksp/?highlight=ksp#ksp-linear-system-solvers) of the method for more information. # PETSc is not a required dependency of DOLFINx, and therefore we explicitly import the DOLFINx wrapper for interfacing with PETSc. -# To ensure that the options passed to the LinearProblem is only used for the given KSP solver, we pass a **unique** option prefix as well. +# To ensure that the options passed to the {py:class}`LinearProblem` +# is only used for the given KSP solver, we pass a **unique** option prefix as well. +# + from dolfinx.fem.petsc import LinearProblem problem = LinearProblem( @@ -208,20 +260,28 @@ petsc_options_prefix="Poisson", ) uh = problem.solve() +# - -# Using `problem.solve()` we solve the linear system of equations and return a `dolfinx.fem.Function` containing the solution. +# Using {py:meth}`problem.solve()` we solve the linear system of equations and +# return a {py:class}`Function` containing the solution. +# # ## Computing the error -# Finally, we want to compute the error to check the accuracy of the solution. We do this by comparing the finite element solution `u` with the exact solution. -# We do this by interpolating the exact solution into the the $P_2$-function space. +# Finally, we want to compute the error to check the accuracy of the solution. +# We do this by comparing the finite element solution `u` with the exact solution. +# First we interpolate the exact solution into a function space that contains it V2 = fem.functionspace(domain, ("Lagrange", 2)) -uex = fem.Function(V2) +uex = fem.Function(V2, name="u_exact") uex.interpolate(lambda x: 1 + x[0] ** 2 + 2 * x[1] ** 2) -# We compute the error in two different ways. First, we compute the $L^2$-norm of the error, defined by $E=\sqrt{\int_\Omega (u_D-u_h)^2\mathrm{d} x}$. -# We use UFL to express the $L^2$-error, and use `dolfinx.fem.assemble_scalar` to compute the scalar value. -# In DOLFINx, `assemble_scalar` only assembles over the cells on the local process. This means that if we use 2 processes to solve our problem, we need to gather the solution to one (or all the processes. -# We can do this with the `MPI.allreduce` function. +# We compute the error in two different ways. +# First, we compute the $L^2$-norm of the error, defined by $E=\sqrt{\int_\Omega (u_D-u_h)^2\mathrm{d} x}$. +# We use UFL to express the $L^2$-error, and use {py:func}`dolfinx.fem.assemble_scalar` to compute the scalar value. +# In DOLFINx, {py:func}`assemble_scalar` +# only assembles over the cells on the local process. +# This means that if we use 2 processes to solve our problem, +# we need to accumulate the local contributions to get the global error (on one or all processes). +# We can do this with the {func}`Comm.allreduce` function. L2_error = fem.form(ufl.inner(uh - uex, uh - uex) * ufl.dx) error_local = fem.assemble_scalar(L2_error) @@ -230,9 +290,13 @@ # Secondly, we compute the maximum error at any degree of freedom. # As the finite element function $u$ can be expressed as a linear combination of basis functions $\phi_j$, spanning the space $V$: # $ u = \sum_{j=1}^N U_j\phi_j.$ -# By writing `problem.solve()` we compute all the coefficients $U_1,\dots, U_N$. These values are known as the _degrees of freedom_ (dofs). We can access the degrees of freedom by accessing the underlying vector in `uh`. -# However, as a second order function space has more dofs than a linear function space, we cannot compare these arrays directly. -# As we already have interpolated the exact solution into the first order space when creating the boundary condition, we can compare the maximum values at any degree of freedom of the approximation space. +# By writing `problem.solve()` we compute all the coefficients $U_1,\dots, U_N$. +# These values are known as the _degrees of freedom_ (dofs). +# We can access the degrees of freedom by accessing the underlying vector in `uh`. +# However, as a second order function space has more dofs than a linear function space, +# we cannot compare these arrays directly. +# As we already have interpolated the exact solution into the first order space when creating the boundary condition, +# we can compare the maximum values at any degree of freedom of the approximation space. error_max = numpy.max(numpy.abs(uD.x.array - uh.x.array)) # Only print the error on one process @@ -242,9 +306,11 @@ # ## Plotting the mesh using pyvista # We will visualizing the mesh using [pyvista](https://docs.pyvista.org/), an interface to the VTK toolkit. -# We start by converting the mesh to a format that can be used with `pyvista`. -# To do this we use the function `dolfinx.plot.vtk_mesh`. The first step is to create an unstructured grid that can be used by `pyvista`. -# We need to start a virtual framebuffer for plotting through docker containers. You can print the current backend and change it with `pyvista.set_jupyter_backend(backend)` +# We start by converting the mesh to a format that can be used with {py:mod}`pyvista`. +# To do this we use the function {py:func}`dolfinx.plot.vtk_mesh`. +# It creates the data required to create a {py:class}`pyvista.UnstructuredGrid`. +# On Linux (in docker) we need to start a virtual framebuffer for plotting through docker containers. +# You can print the current backend and change it with {py:func}`pyvista.set_jupyter_backend`. # + import pyvista @@ -261,11 +327,12 @@ # - # There are several backends that can be used with pyvista, and they have different benefits and drawbacks. -# See the [pyvista documentation](https://docs.pyvista.org/user-guide/jupyter/index.html#state-of-3d-interactive-jupyterlab-plotting) for more information and installation details. -# n this example and the rest of the tutorial we will use [panel](https://github.com/holoviz/panel). +# See the [pyvista documentation](https://docs.pyvista.org/user-guide/jupyter/index.html#state-of-3d-interactive-jupyterlab-plotting) +# for more information and installation details. -# We can now use the `pyvista.Plotter` to visualize the mesh. We visualize it by showing it in 2D and warped in 3D. -# In the jupyter notebook environment, we use the default setting of `pyvista.OFF_SCREEN=False`, which will render plots directly in the notebook. +# We can now use the {py:class}`pyvista.Plotter` to visualize the mesh. We visualize it by showing it in 2D and warped in 3D. +# In the jupyter notebook environment, we use the default setting of `pyvista.OFF_SCREEN=False`, +# which will render plots directly in the notebook. plotter = pyvista.Plotter() plotter.add_mesh(grid, show_edges=True) @@ -276,11 +343,15 @@ figure = plotter.screenshot("fundamentals_mesh.png") # ## Plotting a function using pyvista -# We want to plot the solution `uh`. As the function space used to defined the mesh is disconnected from the function space defining the mesh, we create a mesh based on the dof coordinates for the function space `V`. We use `dolfinx.plot.vtk_mesh` with the function space as input to create a mesh with mesh geometry based on the dof coordinates. +# We want to plot the solution `uh`. +# As the function space used to defined the mesh is decoupled from the representation of the mesh, +# we create a mesh based on the dof coordinates for the function space `V`. +# We use {py:func}`dolfinx.plot.vtk_mesh` with the function space as input to create a mesh with +# mesh geometry based on the dof coordinates. u_topology, u_cell_types, u_geometry = plot.vtk_mesh(V) -# Next, we create the `pyvista.UnstructuredGrid` and add the dof-values to the mesh. +# Next, we create the {py:class}`pyvista.UnstructuredGrid` and add the dof-values to the mesh. u_grid = pyvista.UnstructuredGrid(u_topology, u_cell_types, u_geometry) u_grid.point_data["u"] = uh.x.array.real @@ -300,7 +371,9 @@ plotter2.show() # ## External post-processing -# For post-processing outside the python code, it is suggested to save the solution to file using either `dolfinx.io.VTXWriter` or `dolfinx.io.XDMFFile` and using [Paraview](https://www.paraview.org/). This is especially suggested for 3D visualization. +# For post-processing outside the python code, it is suggested to save the solution to file using either +# {py:class}`dolfinx.io.VTXWriter` or {py:class}`dolfinx.io.XDMFFile` and using [Paraview](https://www.paraview.org/). +# This is especially suggested for 3D visualization. # + from dolfinx import io diff --git a/docker/Dockerfile b/docker/Dockerfile index eb728854..5d039524 100644 --- a/docker/Dockerfile +++ b/docker/Dockerfile @@ -1,6 +1,6 @@ # Execute from root of repo as: docker buildx build --platform=linux/arm64,linux/amd64 -f docker/Dockerfile ./docker/ --progress=plain -FROM ghcr.io/fenics/dolfinx/lab:v0.9.0 +FROM ghcr.io/fenics/dolfinx/lab:nightly ARG TARGETPLATFORM diff --git a/index.ipynb b/index.ipynb index e4681a36..1cce4752 100644 --- a/index.ipynb +++ b/index.ipynb @@ -8,7 +8,8 @@ "# The FEniCSx tutorial\n", "Author: Jørgen S. Dokken\n", "\n", - "These webpages give a concise overview of the functionality of [DOLFINx](https://github.com/FEniCS/dolfinx/), including a gentle introduction to the finite element method. This webpage is an adaptation of the FEniCS tutorial {cite}`FenicsTutorial`.\n", + "These webpages give a concise overview of the functionality of [DOLFINx](https://github.com/FEniCS/dolfinx/), including a gentle introduction to the finite element method.\n", + "This webpage started as an adaptation of the FEniCS tutorial {cite}`FenicsTutorial`, but has evolved into a larger subset of tutorials.\n", "\n", "DOLFINx can be used as either C++ or Python software, but this tutorial will focus on Python programming, as it is the simplest and most effective approach for beginners. After having gone through this tutorial, the reader should familiarize themselves with the DOLFINx [documentation](https://docs.fenicsproject.org/dolfinx/main/python/), which includes the API and numerous demos.\n", "\n", diff --git a/pyproject.toml b/pyproject.toml index 04004c8f..ef80ddcd 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -4,7 +4,7 @@ build-backend = "setuptools.build_meta" [project] name = "DOLFINx_Tutorial" -version = "0.9.0" +version = "0.10.0.dev0" dependencies = [ "jupyter-book", "meshio", @@ -13,14 +13,14 @@ dependencies = [ "pandas", "tqdm", "pyvista[all]>=0.43.0", - "fenics-dolfinx>=0.9.0", + "fenics-dolfinx>0.9.0", "sphinx-codeautolink", ] [project.optional-dependencies] dev = ["pdbpp", "ipython", "jupytext", "ruff", "pre-commit"] netgen = [ - "ngsPETSc@git+https://github.com/NGSolve/ngsPETSc.git@dokken/gmshio-dolfinx", + "ngsPETSc@git+https://github.com/NGSolve/ngsPETSc.git", ] # Has to be optional as we cant get netgen on linux/arm64 [tool.setuptools]