Run this tutorial
Click here to run this tutorial on mybinder.org:
Example: working with gmsh grids#
1: creating a gmsh file#
We use pyMORs PolygonalDomaindescription and discretize_gmsh to obtain a grid file that gmsh can read. Note that dune-grid can only handle gmsh version 2 files, we have thus installed gmsh version 2.16 in this virtualenv. For newer versions of gmsh, you need to follow these instructions.
# wurlitzer: display dune's output in the notebook
%load_ext wurlitzer
%matplotlib notebook
import numpy as np
np.warnings.filterwarnings('ignore') # silence numpys warnings
from pymor.analyticalproblems.domaindescriptions import PolygonalDomain
L_shaped_domain = PolygonalDomain(points=[
[0, 0],
[0, 1],
[-1, 1],
[-1, -1],
[1, -1],
[1, 0],
], boundary_types={'dirichlet': [1, 2, 3, 4, 5, 6]})
from pymor.discretizers.builtin.domaindiscretizers.gmsh import discretize_gmsh
grid, bi = discretize_gmsh(L_shaped_domain, msh_file_path='L_shaped_domain.msh')
Info : Running 'gmsh /tmp/tmphrc9f6p5.geo -2 -algo del2d -clscale 1.0 -o L_shaped_domain.msh' [Gmsh 2.16.0, 1 node, max. 1 thread]
Info : Started on Wed May 4 12:02:41 2022
Info : Reading '/tmp/tmphrc9f6p5.geo'...
Info : Done reading '/tmp/tmphrc9f6p5.geo'
Info : Finalized high order topology of periodic connections
Info : Meshing 1D...
Info : Meshing curve 1 (Line)
Info : Meshing curve 2 (Line)
Info : Meshing curve 3 (Line)
Info : Meshing curve 4 (Line)
Info : Meshing curve 5 (Line)
Info : Meshing curve 6 (Line)
Info : Done meshing 1D (0.000184 s)
Info : Meshing 2D...
Info : Meshing surface 8 (Plane, Delaunay)
Info : Done meshing 2D (0.00416803 s)
Info : 81 vertices 166 elements
Info : Writing 'L_shaped_domain.msh'...
Info : Done writing 'L_shaped_domain.msh'
Info : Stopped on Wed May 4 12:02:41 2022
2: discretization using pyMOR#
We may use pyMOR to solve an elliptic PDE on this grid.
from pymor.analyticalproblems.functions import ConstantFunction
from pymor.analyticalproblems.elliptic import StationaryProblem
problem = StationaryProblem(
domain=L_shaped_domain,
rhs=ConstantFunction(1, dim_domain=2),
diffusion=ConstantFunction(1, dim_domain=2),
)
from pymor.discretizers.builtin import discretize_stationary_cg
fom, fom_data = discretize_stationary_cg(problem, grid=grid, boundary_info=bi)
fom.visualize(fom.solve())
3: using the gmsh grid in dune#
dune-grid only supports gmsh version 2 files, and only a subset of the specification.
This virtualenv includes the gmsh version 2.16 (as visible in the output of the discretize_gmsh command above), but we still need to clean up the mesh file for dune-grid to correctly parse it.
In particular, we need to remove the boundary type definition (which we do not require, we have our own boundary info), which is achieved by the following bash code (Note that you have to provide the same filename here as in the call to discretize_gmsh):
# remove all lines between $PhysicalNames and $EndPhysicalNames ...
!sed '/^\$PhysicalNames/,/^\$EndPhysicalNames/{//!d;};' -i L_shaped_domain.msh
# ... and remove those two lines as well:
!sed '/^\$PhysicalNames/d' -i L_shaped_domain.msh
!sed '/^\$EndPhysicalNames/d' -i L_shaped_domain.msh
from dune.xt.grid import make_gmsh_grid, Dim, Simplex, visualize_grid
grid = make_gmsh_grid('L_shaped_domain.msh', Dim(2), Simplex())
Reading 2d Gmsh grid...
version 2.2 Gmsh file detected
file contains 81 nodes
file contains 160 elements
number of real vertices = 81
number of boundary elements = 32
number of elements = 128
This grid can now be used as any other grid, e.g. for visualization …
_ = visualize_grid(grid)
… or discretization:
from dune.gdt import visualize_function
from discretize_elliptic_cg import discretize_elliptic_cg_dirichlet_zero
u_h = discretize_elliptic_cg_dirichlet_zero(grid, diffusion=1, source=1)
_ = visualize_function(u_h)
Download the code:
example__gmsh_grid.md
example__gmsh_grid.ipynb