-
Notifications
You must be signed in to change notification settings - Fork 170
Added demo for performing adaptive multigrid #4535
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
base: main
Are you sure you want to change the base?
Conversation
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Looks really good! Small comments.
demos/adaptive_multigrid/demo.py.rst
Outdated
================================================= | ||
|
||
The purpose of this demo is to show how to use Firedrake's multigrid solver on a hierarchy of adaptively refined Netgen meshes. | ||
We will first have a look at how to use the AdaptiveMeshHierarchy to construct the mesh hierarchy with Netgen meshes, then we will consider a solution to the Poisson problem on an L shaped domain. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Be consistent about "L shaped" vs "L-shaped" (I prefer the latter)
demos/adaptive_multigrid/demo.py.rst
Outdated
ngmsh = geo.GenerateMesh(maxh=0.5) | ||
mesh = Mesh(ngmsh) | ||
|
||
It is important to convert the initial Netgen mesh into a Firedrake mesh before constructing the AdaptiveMeshHierarchy. To call the constructor to the hierarchy, we must pass the initial mesh as a list. Our initial mesh looks like this: |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Why must it be a list?
demos/adaptive_multigrid/demo.py.rst
Outdated
|
||
- \nabla^2 u = f \text{ in } \Omega, \quad u = 0 \text{ on } \partial \Omega | ||
|
||
Our approach strongly follows the similar problem in this `lecture course <https://github.com/pefarrell/icerm2024>`_. We define the function solve_poisson. The first lines correspond to finding a solution in the CG1 space. The variational problem is formulated with F, where f is the constant function equal to 1. Since we want Dirichlet boundary conditions, we construct the DirichletBC object and apply it to the entire boundary: :: |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Backticks on solve_poisson
demos/adaptive_multigrid/demo.py.rst
Outdated
}, | ||
} | ||
|
||
For more information about patch relaxation, see `Using patch relaxation for multigrid <https://www.firedrakeproject.org/demos/poisson_mg_patches.py.html>`_. The initial solution is shown below |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Full stop
demos/adaptive_multigrid/demo.py.rst
Outdated
.. math:: | ||
\eta_K^2 = h_K^2 \int_K | f + \nabla^2 u_h |^2 \mathrm{d}x + \frac{h_K}{2} \int_{\partial K \setminus \partial \Omega} \llbracket \nabla u_h \cdot n \rrbracket^2 \mathrm{d}s, | ||
|
||
where :math:`K` is the element, :math:`h_K` is the diameter of the element, :math:`n` is the normal, and :math:`\llbracket \cdot \rrbracket` is the jump operator. The a-posteriori estimator is computed using the solution at the current level :math:`h`. We can use a trick to compute the estimator on each element. We transform the above estimator into the variational problem |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Say something like "Integrating over the domain and using the fact that the functions are piecewise constant on each cell, we transform ..."
demos/adaptive_multigrid/demo.py.rst
Outdated
error_est = sqrt(eta_.dot(eta_)) | ||
return (eta, error_est) | ||
|
||
The next step is to choose which elements to refine. For this we Dörfler marking, developed by Professor Willy Dörfler: |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Remove ", developed by Professor Willy Dörfler" and instead cite the paper
demos/adaptive_multigrid/demo.py.rst
Outdated
if i != max_iterations - 1: | ||
amh.adapt(eta, theta) | ||
|
||
To perform Dörfler marking, refine the current mesh, and add the mesh to the AdaptiveMeshHierarchy, we us the amh.adapt(eta, theta) method. In this method the input is the recently computed error estimator :math:`eta` and the Dörfler marking parameter :math:`theta`. The method always performs this on the current fine mesh in the hierarchy. There is another method for adding a mesh to the hierarchy. This is the amh.add_mesh(mesh). In this method, refinement on the mesh is performed externally by some custom procedure and the resulting mesh directly gets added to the hierarchy. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Backticks on amh.add_mesh(mesh)
demos/adaptive_multigrid/demo.py.rst
Outdated
+------------------------------------+------------------------------------+ | ||
|
||
|
||
The convergence follows the expected behavior: |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
expected (optimal) behaviour
demos/adaptive_multigrid/demo.py.rst
Outdated
error_est = sqrt(eta_.dot(eta_)) | ||
return (eta, error_est) | ||
|
||
The next step is to choose which elements to refine. For this we Dörfler marking [Dörfler1996]: |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
missing 'use'
Does this PR need to do anything to link it on https://www.firedrakeproject.org/documentation.html ? |
Yes. You need to add it to: You also need to enable testing by adding it to: |
Description
Only added demo for performing adaptive multigrid using the AdaptiveMeshHierarchy and AdaptiveTransferManager.
Depends on #4489, #4511, and NGSolve/ngsPETSc#87