Grid construction¶
- dune.mmesh.mmesh(constructor, dimgrid=None, **parameters)[source]
Create an MMesh grid.
- Parameters
constructor – Grid constructor, e.g. (dune.grid.reader.gmsh, ‘grid.msh’) or dune.grid.cartesianDomain([0,0],[1,1],[10,10]).
dimgrid (int, optional) – The world dimension (must be 2 or 3). Might be guessed by constructor.
**parameters – Additional parameters.
- Returns
An MMesh grid instance.
Examples
Creating simple structured grid:
from dune.grid import cartesianDomain from dune.mmesh import mmesh grid = mmesh(cartesianDomain([0,0], [1,1], [10,10]))
Using a .msh file:
from dune.grid import reader from dune.mmesh import mmesh grid = mmesh((reader.gmsh, "grid.msh"), 2) igrid = grid.hierarchicalGrid.interfaceGrid
where grid.msh can be generated by gmsh [GR09]. For instance:
name = "grid.msh" h = 0.1 hf = 0.01 import gmsh gmsh.initialize() gmsh.option.setNumber("Mesh.MshFileVersion", 2.2) gmsh.model.add("name") geo = gmsh.model.geo p1 = geo.addPoint(0, 0, 0, h, 1) p2 = geo.addPoint(1, 0, 0, h, 2) p3 = geo.addPoint(1, 1, 0, h, 3) p4 = geo.addPoint(0, 1, 0, h, 4) l1 = geo.addLine(p1, p2, 1) l2 = geo.addLine(p2, p3, 2) l3 = geo.addLine(p3, p4, 3) l4 = geo.addLine(p4, p1, 4) p5 = geo.addPoint(0.25, 0.5, 0, hf, 5) p6 = geo.addPoint(0.75, 0.5, 0, hf, 6) lf = geo.addLine(p5, p6, 10) geo.addCurveLoop([l1, l2, l3, l4], 1) geo.addPlaneSurface([1], 0) geo.synchronize() gmsh.model.mesh.embed(1, [lf], 2, 0) gmsh.model.mesh.generate(dim=2) gmsh.write(name) gmsh.finalize()
- GR09
Geuzaine and J.-F. Remacle. Gmsh: a three-dimensional finite element mesh generator with built-in pre- and post-processing facilities. International Journal for Numerical Methods in Engineering 79(11), pp. 1309-1331, 2009.