dune-grid
2.2.0
|
Classes for reading a macrogrid file in the dune macrogrid format (dgf) More...
Classes | |
class | Dune::DGFWriter< GV > |
write a GridView to a DGF file More... | |
Modules | |
DGF grid parameter for different grids |
Classes for reading a macrogrid file in the dune macrogrid format (dgf)
The DGF format allows the simple description of macrogrids, which can be used to construct Dune grids independent of the underlying implementation. Due to the generality of the approach only a subset of the description language for each grid implementation is available through this interface.
There are two ways of constructing Dune grids using DGF files:
ALBERTAGRID
, ALUGRID_CUBE
, ALUGRID_SIMPLEX
, ALUGRID_CONFORM
, SGRID
, UGGRID
, or YASPGRID
and the integer GRIDDIM
one obtains a definition of the type GridType
and the required header files for the desired grid and the macrogrid parser are included.After a grid type (denoted with GridType
in the following) is selected in some way, the grid can be constructed either by calling
Dune::GridPtr<GridType> gridptr(filename, mpiHelper.getCommunicator() );
or
Dune::GridPtr<GridType> gridptr(filename);
or
Dune::GridPtr<GridType> gridptr; ... gridptr=Dune::GridPtr<GridType>(filename);
where in the second and third example MPIHelper::getCommunicator()
is selected as default value; filename
is the name of the dgf file. This creates an auto pointer like object Dune::GridPtr<GridType>
holding a pointer to a the grid instance described through the dgf file. Access to the grid is gained by calling the operator * of GridPtr
.
GridType & grid = *gridptr;
Like in all auto_ptr like objects the grid instance is destroyed if the GridPtr goes out of scope. To be able to destruct the GridPtr object without losing the grid call the release method:
GridType *grid; { GridPtr< GridType > gridPtr( filename.c_str(), mpiHelper.getCommunicator() ); grid = gridPtr.release(); }
Remarks:
Dune::MPIHelper::MPICommunicator
which defaults to MPI_COMM_WORLD
for parallel runs or some default value for serial runs.GridType
class is called - if one is available.We assume in the following that the type GridType
is suitably defined, denoting a dune grid type in dimworld
space dimension. A point in dimworld space is called a vector.
In general dgf files consists of one or more blocks, each starting with a keyword and ending with a line starting with a # symbol. In a block each full line is parsed from the beginning up to the first occupance of a % symbol, which can be used to include comments. Trailing whitespaces are ignored during line parsing. Also the parsing is not case sensitive.
Some example files are given below (Examples).
DGF files must start with the keyword DGF. Files passed to the grid parser not starting with this keyword are directly passed to a suitable constructor in the GridType class.
In the following all blocks are briefly described in alphabetical order. The construction of the grid is detailed below. In general lines are parse sequentially, with the exception of lines starting with a special keyword which are treated separately.
dimworld
space. The first entry of the line is a boundary identifier.dimworld2
vertex indices (see Vertex block) describing one cube of the macrogrid.0
...2dimworld describing the mapping from the local vertex numbers in the reference elements used in the dgf file to the local vertex numbers in dune reference element.n
it is possible to add n
parameters to each element of the grid. These double valued parameters then have to be added to the definition lines for the elements behind the vertex numbers.dimworld
integers. used to determine the initial partition of the interval into equal sized cubes. More than one interval can be described in this fashion.dimworld+1
vertex indices (see Vertex block) describing one simplex of the macrogrid.dimworld=3
or Triangle (http://www.cs.cmu.edu/~quake/triangle.html) for dimworld=2
. For more detail see Using Tetgen/Triangle.1 0, 0 1 + 1 0
For simplicity we first describe how the grid is manually constructed, i.e., if no Simplexgenerator block is found. Details on how Tetgen/Triangle can be used is detailed below (see Using Tetgen/Triangle). How to access the element and vertex parameter is detailed in Section Accessing parameters. The details of the construction are logged in the file dgfparser.log.
The construction of the grid depends on the type of elements the given by GridType
can handle.
(Dune::SGrid , or Dune::YaspGrid )
The grid is constructed using only the information from the first three lines of the Interval block.
(Dune::AlbertaGrid, or Dune::ALUSimplexGrid<3,3>, and Dune::ALUSimplexGrid<2,2>)
The vertices and elements of the grid are constructed in the following three steps:
If the simplex generation process was successful, boundary ids are assigned to all boundary faces of the macrogrid.
Remark:
Bisection:
The grid is constructed using the information from the Interval block, if present; otherwise the Vertex and Cube blocks are used.
The boundary ids are assigned in the same manner as for simplex grids described above.
(Dune::UGGrid
)
Note that in this version only grids consisting of one element type are constructed even if the implemented grid allows for mixed elements.
The vertices and elements of the grid are constructed in one of the following stages:
Simplex
block is present otherwise a cube grid is constructed.Boundary ids are assigned in the same manner as described for simplex grids.
In addition to the element and vertex information it is possible to add a set of parameters through the DGF file. Using the construction mechanism via the Dune::GridPtr<GridType> class these parameters can be accessed using the Dune::GridPtr<GridType>::parameters(const Entity& en) method. Depending on the codimentsion of en
the method returns either the element or vertex parameters of that entity in the DGF file. The number of parameters for a given codimension can be retrived using the method Dune::GridPtr<GridType>::nofParameters. In the case of parallel runs only the process with rank zero reads the dgf file and thus the parameters are only available on that processor. Do not call loadBalance (or for that matter globalRefine, adapt) on the grid instance before retrieving the parameters. To simplify the handling of parameters in a parallel run the Dune::GridPtr<GridType>::loadBalance method can be used to distribute the macro grid and the parameters. Afterwards the parameters can be extracted on each processor using the method GridPtr<GridType>::parameters.
The freely available simplex grid generators are direcltly called via system call through the dgfparser. Therefore one should either add the path containing the executables of Triangle and/or Tetgen to the environment variable PATH or use the path option described below. One can either use a file in Tetgen/Triangle format to directly generate a macro grid or one can prescribe vertices and boundary faces which will be used to generate the grid directly in the DGF file:
file name mesh
Some identifiers can be used to influence the quality of the generated mesh:
In this case the grid is constructed using two calls to Tetgen/Triangle. The details of the call are logged in the dgfparser.log file.
Note: vertex parameters are interpolated in triangle but tetgen assigns a zero to all newly inserted vertices; therfore quality enhancement and vertex parameters should not be combined in 3d. On the other hand element parameters and boundary ids can be used together with quality enhancement.
The remaining identifiers are
Note that parameters can be attached to the vertices and elements of the grid as described in the Triangle/Tetgen manual. Then can be retrieved as described in Section Accessing parameters.
Download
In two space dimensions:
In three space dimensions:
A tessellation of the unit square into six simplex entities. Some boundary segments on the lower and right boundary are given their own id the remaining are given a default id.
DGF Vertex % the verticies of the grid -1 -1 % vertex 0 -0.2 -1 % vertex 1 1 -1 % vertex 2 1 -0.6 % vertex 3 1 1 % vertex 4 -1 1 % vertex 5 0.2 -0.2 % vertex 6 # SIMPLEX % a simplex grid 0 1 5 % triangle 0, verticies 0,1,5 1 3 6 % triangle 1 1 2 3 % triangle 2 6 3 4 % triangle 3 6 4 5 % triangle 4 6 5 1 % triangle 5 # BOUNDARYSEGMENTS 2 1 2 % between vertex 1,2 use id 2 2 2 3 % between vertex 2,3 use id 2 4 3 4 % between vertex 3,4 use id 4 # BOUNDARYDOMAIN default 1 % all other boundary segments have id 1 # # examplegrid1s.dgf
A tessellation into cubes using the same vertices as before,
DGF Vertex % the verticies of the grid -1 -1 % vertex 0 -0.2 -1 % vertex 1 1 -1 % vertex 2 1 -0.6 % vertex 3 1 1 % vertex 4 -1 1 % vertex 5 0.2 -0.2 % vertex 6 # CUBE % a cube grid 0 1 5 6 % cube 0, verticies 0,1,5,6 1 2 6 3 % cube 1 6 3 5 4 % cube 2 # BOUNDARYSEGMENTS 2 1 2 % between vertex 1,2 use id 2 2 2 3 % between vertex 2,3 use id 2 4 3 4 % between vertex 3,4 use id 4 # BOUNDARYDOMAIN default 1 % all other boundary segments have id 1 # # examplegrid1c.dgf
Using the last input file with a simplex grid or by adding an empty Simplex
block
Simplex
#
leads to the following macro triangulation
Automatic tessellation using Triangle, with vertices defined as in the example Manual Grid Construction :
DGF Simplexgenerator min-angle 30 display 0 % to use Triangle from a certain path, uncomment this path line %path $HOME/bin % path to Triangle # Vertex % the verticies of the grid -1 -1 % vertex 0 -0.2 -1 % vertex 1 1 -1 % vertex 2 1 -0.6 % vertex 3 1 1 % vertex 4 -1 1 % vertex 5 0.2 -0.2 % vertex 6 # BOUNDARYDOMAIN default 1 % all other boundary segments have id 1 # # examplegrid1gen.dgf
The quality of the grid can be enhanced by adding the line
min-angle 30
in the Simplexgenerator
block
Automatic tessellation using Triangle, with vertices are defined on a Cartesian grid with two additional vertices in the top right corner and one vertex outside the unit square.
All boundary are given a default id.
DGF Interval -1 -1 % first corner 1 1 % second corner 4 4 % 4 cells in x and 4 in y direction # Vertex % some additional points 0.75 0.75 % inside the interval 0.9 0.9 % also inside the interval 1.25 1.25 % and one outside of the interval # Simplexgenerator display 0 % show result using Triangle viewer % to use Triangle from a certain path, uncomment path line %path $HOME/bin % path to Triangle # BOUNDARYDOMAIN default 1 % all boundaries have id 1 #BOUNDARYDOMAIN # examplegrid2a.dgf
Adding some quality enhancement. The boundaries are numbered counterclockwise starting at the left boundary from one to four.
DGF Interval -1 -1 % first corner 1 1 % second corner 4 4 % 4 cells in x and 4 in y direction # Vertex % some additional points 0.75 0.75 % inside the interval 0.9 0.9 % also inside the interval 1.25 1.25 % and one outside of the interval # Simplexgenerator min-angle 30 % quality enhancment display 0 % show result using Triangle viewer % area restriction % to use Triangle from a certain path, uncomment this path line %path $HOME/bin % path to Triangle # BOUNDARYDOMAIN 1 1 -1 1.25 1.25 % right boundary has id 1 2 -1 1 1.25 1.25 % upper boundary has id 2 3 -1 -1 -1 1 % left boundary has id 3 4 -1 -1 1 -1 % lower boundary has id 4 #BOUNDARYDOMAIN # examplegrid2b.dgf
Using both quality enhancement and a maximum area restriction. The bottom boundary is given the id 1 all other boundaries have id 2; here we do not use a default value.
DGF Interval -1 -1 % first corner 1 1 % second corner 4 4 % 4 cells in x and 4 in y direction # Vertex % some additional points 0.75 0.75 % inside the interval 0.9 0.9 % also inside the interval 1.25 1.25 % and one outside of the interval # Simplexgenerator min-angle 30 % quality enhancment max-area 0.01 % area restriction display 0 % show result using Triangle viewer % to use Triangle from a certain path, uncomment this path line %path $HOME/bin % path to Triangle # BOUNDARYDOMAIN 1 -1 -1 1 -1 % lower boundary has id 1 2 -10 -10 10 10 % all others have id 2 (could use default keyword...) #BOUNDARYDOMAIN # examplegrid2c.dgf
A similar grid is generated by prescribing the boundary of the domain:
DGF Vertex % some additional points -1 -1 1 -1 1 1 -1 1 0.75 0.75 % inside the interval 0.9 0.9 % also inside the interval 1.25 1.25 % and one outside of the interval # Simplexgenerator min-angle 30 % quality enhancment max-area 0.01 % area restriction display 0 % show result using Triangle viewer % to use Triangle from a certain path, uncomment this path line %path $HOME/bin % path to Triangle # BOUNDARYSEGMENTS 1 0 1 2 1 6 3 0 # # examplegrid2d.dgf
We use the same domain as in the previous example but include vertex and element parameters:
DGF vertex parameters 2 -1 -1 0 -1 1 -1 1 1 1 1 2 1 -1 1 1 -1 0.75 0.75 2 0.75 0.90000000000000002 0.90000000000000002 2 0.90000000000000002 1.25 1.25 -1 1.25 # simplex parameters 1 0 4 3 1 4 0 1 -1 3 4 5 1 5 1 2 1 1 5 4 -1 2 6 3 1 6 2 1 -1 5 2 3 -1 # boundarydomain default 1 # Simplexgenerator min-angle 30 % quality enhancment max-area 0.01 % area restriction display 0 % show result using Triangle viewer % to use Triangle from a certain path, uncomment this path line % path $HOME/bin % path to Triangle # # examplegrid2e.dgf
The results in piecewise constant data on the elements and piecewise linear date using the vertex parameters:
A automatic tessellation of the unit square using a Cartesian Grid. All boundaries have id 1.
DGF Interval 0 0 % first corner 1 1 % second corner 4 8 % 4 cells in x and 8 in y direction # Simplex # GridParameter % set overlap to 1 overlap 1 % set closure to none for UGGrid closure none % generate copies in UGGrid copies yes % set heap size for UGGrid heapsize 1000 # PERIODICFACETRANSFORMATION % set periodic boundaries in x direction % 1 0, 0 1 + 1 0 # BOUNDARYDOMAIN default 1 % default boundary id #BOUNDARYDOMAIN # examplegrid5.dgf
If UGGrid<2,2> is used the result would be the same as for SGrid<2,2>. If an empty Simplex
Block
Simplex
#
is added than the same simplex grid as for AlbertaGrid<2,2> would be constructed.
The following example shows a DGF file that projects 3 sides of a quadrilateral onto the surrounding circle:
DGF VERTEX -1 -1 1 -1 -1 1 1 1 # CUBE 0 1 2 3 # PROJECTION function p(x) = sqrt(2) * x / |x| function id(x) = x segment 0 1 id default p # # example-projection.dgf
An automatic tessellation of the unit square using a Cartesian Grid is shown. All boundaries have id 1 except boundary segment on the lower boundary which have value 2.
First we use only the interval block:
DGF Interval 0 0 0 % first corner 10 10 10 % second corner 4 4 3 % 5 cells in x, 2 in y, and 10 in z direction # simplex # BOUNDARYDOMAIN default 1 % default boundary id % all boundary segments in the interval [(-1,-1,-1),(11,11,0)] % have boundary id 2, i.e., the bottom boundary 2 -1 -1 -1 11 11 0 # # examplegrid6.dgf
Now the vertices are still defined through the interval block; the simplicies are constructed using Tetgen (note the comment symbol % in the Simplexgenerator block):
DGF Simplexgenerator %min-angle 1.1 % quality enhancment %max-area 0.1 % maximum element area restriction display 0 % show result using the Tetgen viewer % to use TetGen from a certain path, uncomment this path line %path $HOME/bin % path to Tetgen # Interval 0 0 0 % first corner 10 10 10 % second corner 2 5 10 % 5 cells in x, 2 in y, and 10 in z direction # BOUNDARYDOMAIN default 1 % default boundary id % all boundary segments in the interval [(-1,-1,-1),(11,11,0)] % have boundary id 2, i.e., the bottom boundary 2 -1 -1 -1 11 11 0 # # examplegrid7.dgf
Note that in the grid would be the same as in the above example if ALUCubeGrid<3,3> where used.
Now we can also include some quality enhancement:
First: min-angle = 1.2 (remove the first % in the Simplexgenerator block)
Second: min-angle = 1.2 and max-area = 0.1 (remove both % in the Simplexgenerator block)
This examples show different ways to define a grid for the following domain and boundary ids:
Manual grid construction: Note that the reference element used to define the cube is not the one used in dune so that the mapping has to be given; furthermore the vertex index are numbered starting with 1. If a simplex grid is to be constructed some care must be taken in the numbering of the nodes so that a conforming grid is constructed.
DGF vertex firstindex 1 0 0 5 # 1 0 0 0 # 2 0 -10 0 # 3 0 -10 10 # 4 0 10 10 # 5 0 10 5 # 6 30 10 0 # 7 30 -10 0 # 8 30 -10 10 # 9 30 10 10 # 10 25 10 0 # 11 25 0 0 # 12 25 0 5 # 13 25 10 5 # 14 0 0 10 # 15 30 0 0 # 16 30 0 5 # 17 30 0 10 # 18 30 10 5 # 19 30 -10 5 # 20 25 10 10 # 21 25 0 10 # 22 25 -10 10 # 23 25 -10 5 # 24 25 -10 0 # 25 0 -10 5 # 26 # Cube map 0 1 3 2 4 5 7 6 21 14 13 22 5 6 1 15 23 22 13 24 4 15 1 26 25 24 13 12 3 26 1 2 17 19 7 16 13 14 11 12 17 16 8 20 13 12 25 24 17 20 9 18 13 24 23 22 17 18 10 19 13 22 21 14 # boundarysegments 2 1 2 3 26 2 26 4 15 1 2 1 15 5 6 1 7 16 17 19 1 19 17 18 10 1 8 16 17 20 1 20 17 18 9 4 21 10 19 14 4 11 7 19 14 4 5 6 14 21 4 15 5 22 21 4 4 15 22 23 4 23 9 18 22 4 22 18 10 21 # boundarydomain default 3 # # examplegrid10.dgf
Using multiple intervals: Here the boundary ids are not set correctly; this could be done using different boundarydomains.
DGF Interval % first interval 0 0 5 % first corner 25 10 10 % second corner 5 1 2 % 8 cells in x, 2 in y, and 2 in z direction % second interval 25 -10 0 % first corner 30 10 10 % second corner 3 2 4 % 2 cells in x, 2 in y, and 4 in z direction % third interval 0 -0 0 % first corner 25 -10 10 % second corner 5 1 4 % 2 cells in x, 2 in y, and 4 in z direction # BOUNDARYDOMAIN default 3 % default boundary id % ... further domain ... # # examplegrid11.dgf
By only defining the boundary faces and using Tetgen:
DGF Simplexgenerator min-angle 1.1 % quality enhancment max-area 0.5 % maximum element area restriction display 0 % show result using the Tetgen viewer % to use TetGen from a certain path, uncomment this path line % path $HOME/bin % path to Tetgen # vertex 0 0 5 0 0 0 0 -10 0 0 -10 10 0 10 10 0 10 5 30 10 0 30 -10 0 30 -10 10 30 10 10 25 10 0 25 0 0 25 0 5 25 10 5 # boundarysegments % 2 out, 1 in, 3 slip, 4 reflect 2 0 1 2 3 4 5 % ausgang 1 6 7 8 9 % eingang 3 10 11 12 13 % balken 3 11 12 0 1 % balken 3 13 12 0 5 % balken 3 3 4 9 8 % seite 3 2 3 8 7 % unten 4 6 10 13 5 4 9 % sym. oben 4 6 10 11 1 2 7 % sym. seite # # examplegrid12.dgf
Without the max-area and min-angle keywords the following grid is constructed:
With the quality enhancement active:
Finally we demonstrate the use of parameters. With a simple modification of the above example
DGF vertex firstindex 1 0 0 5 # 1 0 0 0 # 2 0 -10 0 # 3 0 -10 10 # 4 0 10 10 # 5 0 10 5 # 6 30 10 0 # 7 30 -10 0 # 8 30 -10 10 # 9 30 10 10 # 10 25 10 0 # 11 25 0 0 # 12 25 0 5 # 13 25 10 5 # 14 0 0 10 # 15 30 0 0 # 16 30 0 5 # 17 30 0 10 # 18 30 10 5 # 19 30 -10 5 # 20 25 10 10 # 21 25 0 10 # 22 25 -10 10 # 23 25 -10 5 # 24 25 -10 0 # 25 0 -10 5 # 26 # Cube map 0 1 3 2 4 5 7 6 parameters 2 21 14 13 22 5 6 1 15 0 1 23 22 13 24 4 15 1 26 1 -1 25 24 13 12 3 26 1 2 2 1 17 19 7 16 13 14 11 12 3 -1 17 16 8 20 13 12 25 24 4 1 17 20 9 18 13 24 23 22 5 -1 17 18 10 19 13 22 21 14 6 1 # boundarysegments 2 1 2 3 26 2 26 4 15 1 2 1 15 5 6 1 7 16 17 19 1 19 17 18 10 1 8 16 17 20 1 20 17 18 9 4 21 10 19 14 4 11 7 19 14 4 5 6 14 21 4 15 5 22 21 4 4 15 22 23 4 23 9 18 22 4 22 18 10 21 # boundarydomain default 3 # # examplegrid10a.dgf
we can define parameters on each element:
The following example
DGF vertex parameters 1 firstindex 1 0 0 5 5 # 1 0 0 0 0 # 2 0 -10 0 -10 # 3 0 -10 10 0 # 4 0 10 10 20 # 5 0 10 5 15 # 6 30 10 0 40 # 7 30 -10 0 20 # 8 30 -10 10 30 # 9 30 10 10 50 # 10 25 10 0 35 # 11 25 0 0 25 # 12 25 0 5 30 # 13 25 10 5 40 # 14 0 0 10 10 # 15 30 0 0 30 # 16 30 0 5 35 # 17 30 0 10 40 # 18 30 10 5 45 # 19 30 -10 5 25 # 20 25 10 10 45 # 21 25 0 10 35 # 22 25 -10 10 25 # 23 25 -10 5 20 # 24 25 -10 0 15 # 25 0 -10 5 -5 # 26 # Cube map 0 1 3 2 4 5 7 6 21 14 13 22 5 6 1 15 23 22 13 24 4 15 1 26 25 24 13 12 3 26 1 2 17 19 7 16 13 14 11 12 17 16 8 20 13 12 25 24 17 20 9 18 13 24 23 22 17 18 10 19 13 22 21 14 # boundarysegments 2 1 2 3 26 2 26 4 15 1 2 1 15 5 6 1 7 16 17 19 1 19 17 18 10 1 8 16 17 20 1 20 17 18 9 4 21 10 19 14 4 11 7 19 14 4 5 6 14 21 4 15 5 22 21 4 4 15 22 23 4 23 9 18 22 4 22 18 10 21 # boundarydomain default 3 # # examplegrid10b.dgf
defines parameters for each vertex of the grid, leading to a piecewise linear function
Here a .mesh file used to generate a 3d simplex grid through Tetgen. The mesh file is taken from http://www-c.inria.fr/Eric.Saltel/download/download.php
DGF Simplexgenerator % to use TetGen from a certain path, uncomment this path line % path $HOME/dune/modules/tetgen1.4.0 % path to Tetgen % read information from file (uncomment each line seperatly) % file BBEETH1M.d_cut mesh % file Orb_cut mesh % file bunny.p65.param_skin mesh % file pmdc.1 % all grids are for dimension 3 dimension 3 # BOUNDARYDOMAIN default 1 # # examplegrid9.dgf