The Dune Grid Format (DGF)
[I/O]
General
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.Usage
There are two ways of constructing Dune grids using DGF files:
- By including the file gridtype.hh from the dune/grid/io/file/dgfparser directory:
by defining one of the symbolsALBERTAGRID
,ALUGRID_CUBE
,ALUGRID_SIMPLEX
,ALUGRID_CONFORM
,SGRID
,UGGRID
, orYASPGRID
and the integerGRIDDIM
one obtains a definition of the typeGridType
and the required header files for the desired grid and the macrogrid parser are included. - By directly including one of the files dgf*.hh from the dune/grid/io/file/dgfparser directory:
in this case only the required header files for the desired grid and the macrogrid parser are included but no typedef is made.
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() );
Dune::GridPtr<GridType> gridptr(filename);
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;
Remarks:
- The last argument
should
be of the typeDune::MPIHelper::MPICommunicator
which defaults toMPI_COMM_WORLD
for parallel runs or some default value for serial runs. - If the file given through the first argument is not a dgf file a suitable constructure on the
GridType
class is called - if one is available.
Format Description
We assume in the following that the typeGridType
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).
First line
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 constructure in the GridType class.Blocks
In the following all blocks are briefly described in alphabetical order. The construction of the grid is detailed below. In general lines are parser sequentially, with the exception of lines starting with a special keyword which are treated separately.
- Boundarydomain
Each line consists of an integer greater than zero and two vectors describing an interval indimworld
space. The first entry of the line is a boundary identifier.- A special keyword is default followed by a positive integer which is to be used as a default boundary identifier.
- Boundarysegments
Each line consists of a positive integer denoting a boundary id and a set of vertex indices (see Vertex block) describing one boundary patch of the macrogrid. - Cube
Each line consists ofdimworld2
vertex indices (see Vertex block) describing one cube of the macrogrid.- If the ordering of the local vertices of the cube elements does not follow the Dune reference element then the mapping between the local numbering used and the dune reference cube has to be prescribed by a line starting with the keyword map followed by
dimworld2
values describing the mapping between the local vertex numbers. - By means of the special keyword parameters followed by a positive integer
n
it is possible to addn
parameters to each element of the grid. These double valued parameters then have to be added to the definition lines for the elementsi behind the vertex numbers.
- If the ordering of the local vertices of the cube elements does not follow the Dune reference element then the mapping between the local numbering used and the dune reference cube has to be prescribed by a line starting with the keyword map followed by
- Interval
Each interval is described by three lines in this block: the first two lines give the lower and upper vector of an interval, the third line consists ofdimworld
integers. used to determine the initial partition of the interval into equal sized cubes. More than one interval can be described in this fashion. - Simplex
Each line consists ofdimworld+1
vertex indices (see Vertex block) describing one simplex of the macrogrid.- Note that no ordering of local vertices is required.
- Parameters can be added to each element using the parameters keyword as described for cube elements.
- Simplexgenerator
Using this block a simplex grid can be automatically generated using one of the freely available grid generation tools Tetgen (http://tetgen.berlios.de) fordimworld=3
or Triangle (http://www.cs.cmu.edu/~quake/triangle.html) fordimworld=2
. For more detail see Using Tetgen/Triangle. - Vertex
Each line consists of a vector representing a vertex of the macrogrid. The vertices are consecutively numbered. leading to vertex indices which can be referenced in the Simplex, Cube, and Boundarysegment blocks.- By default the numbering starts from zero; to change this behaviour the keyword firstindex followed by a positive integer denoting the number of the first vertex can be used.
- Using the parameters keyword it is possible to add a set of parameters to each vertex of the grid.
- GridParameter
Using this block a grid specific parameters can be passed to certain grids. The following options are implemented at the moment:
- For YaspGrid two options can be choosen:
1. overlap defining the overlap of the grid (default value is zero)
2. periodic defining which dimension should have periodic boundaries, i.e. passing periodic 0 1 will set periodic boundaries for x and y direction. - For UGGrid one option can be choosen:
1. closure (valid values are none or green, which is the default value) will set the closure type of the returned UGGrid.
See the examplegrid5.dgf file for an example.
- For YaspGrid two options can be choosen:
The Grid Construction Process
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 grid given by GridType
can handle.
Cartesian grids
(Dune::SGrid , or Dune::YaspGrid )The grid is constructed using only the information from the first three lines of the Interval block.
Simplex grids
(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:
- The file is parser for an Interval block; if present a Cartesian grid is build for each interval defined in this block and each element is partitioned either into two triangles or into six tetrahedron. The user has to make sure that this process leads to a conforming grid.
- If no Interval block is found, the grid is generated using the information from the Vertex and the Simplex or Cube blocks. If a non-empty Simplex block is found the element information is only taken from there; in the case where no Simplex block is found or in the case where it is empty, the Cube block is read and each cube is partitioned into simplex elements - in 3d this process will not always lead to a correct macro triangulation! Note that no specific ordering of the local numbering of each simplex is required but the cubes must either conform with the Dune reference element or the mapping must be prescribed using the map keyword.
- If the macrogrid was constructed in the third step of the generation process detailed above (i.e. using the vertex block), then the file is first parsed for Boundarysegment block. Here boundary ids can be individually assigned to each boundary segment of the macrogrid by specifying the vertex ids. Boundary segments can either be described as simplex elements or as cube elements.
- All Boundary segments which have not yet been assigned an identifier and lie inside the interval defined by the first line of the Boundarydomain block are assigned the corresponding id. This process is then repeated with the remaining boundary segments using the following intervals defined in the Boundarydomain block. If after this process boundary segments without id remain, then the default id is used if one was specified in the Boundarydomain block - if no default was given, the behavior of the parser is not well defined.
Bisection:
the refinement edge is always chosen to be the longest edge.
Cube grids
(Dune::ALUCubeGrid<3,3>)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.
Mixed grids
(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:
- The file is first parser for an Interval block; if present a Cartesian grid is build; a simplex grid is generated if a
Simplex
block is present otherwise a cube grid is constructed. - If no Interval block is found, the grid is generated using the information from the Vertex block; cube elements are constructed if the Cube block is present otherwise the Simplex blocks is read. If both a Cube and a Simplex block is found, then only the element information from the Cube block is used and each cube is split into simplex elements so that a simplex grid is constructed.
Accessing parameters
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 ofen
the method returns either the element or vertex parameters of that entity in the DGF file which has minimal distance from the midpoint of the entity en
. Note that in this implementation the search procedure requires is O(N)
where N is the number of entities defined in the DGF file - in future releases this will be replaced by a nearest neighbor search. The number of parameters for a given codimension can be retrived using the method Dune::GridPtr<GridType>::nofParameters.Using Tetgen/Triangle
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 th DGF file:- For the first approach use the token file to give a filemane and the type of the file (e.g. node, mesh, ply...). Example: If the filetype is givn, it will be appended to name and this will be passed to Tetgen/Triangle. Additional parameters for the grid generators can be given by the the parameter token. If no file type is given it is assumed that in a previous run of Tetgen/Triangle files name.node and name.ele were generated and these will be used to described the vertices and elements of the Dune grid.
file name mesh
- In the second approach the vertex and the interval blocks (if present) are used to generate the vertices of the grid; the Cube and Simplex blocks are used for element information in the interior and the boundarysegment block is used for boundary information:
- If only a vertex and/or interval block is found the resulting .node file is passed directly to Tetgen/Triangle and a tessellation of the convex hull of the points is generated.
- A more detailed description of the domain is possible by using the Boundarysegment block together with the vertex block. Planar polyhedral boundary faces of the domain can be prescribed in this way. Note that in this case the whole domain boundary has to be defined (see the description of the .poly files in the Tetgen documentation). In 3d each polyhedral face (p0,..,pn) is automatically closed by adding the segment between p0 and pn. To define the whole boundary of a 2d domain using only one polyhedron the closing segment has to be added, i.e., (p0,..,pn,p0).
- If a cube or simplex block is found the element information is also passed to tetgen/triangle together with the parameters - if given. Note that triangle can only handle one region atribute in its .poly files so that only the first parameter is the simplex or cube block can be retrived.
- max-area followed by a positive real number used as an upper bound for the area of all simplicies of the mesh.
- min-angle followed by a positive number. In 2d this limits the angles in resulting mesh from below; in 3d this bounds the radius-edge ratio from above.
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
- The identifier path (followed by a path name) can be used to give search path for Triangle/Tetgen.
- The identifier display followed by 1 can be used to get a first impression of the resulting mesh using the visualization tools distributed with Triangle/Tetgen.
Download
Work in progress
- There should be a mechanism to fix the desired refinement edge for simplex grids. An automatic selection is performed using the longest edge, but it should be possible to turn off this automatic selection.
- A callback to the user between parsing the DGF file and the construction of the dune grid; here e.g. a post-processing of the vertex coordinates could be included.
Examples
In two space dimensions:In three space dimensions:
Manual Grid Construction
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

The resulting grid
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

The resulting grid
Simplex
block Simplex
#

The resulting grid
Automated Grid Construction
Automatic tessellation using Triangle, with vertices defined as in the example dgfexample1:
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 resulting grid
min-angle 30
Simplexgenerator
block

The resulting grid
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

The resulting grid
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

The resulting grid
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

The resulting grid
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

The resulting grid
Using Parameters
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:

The resulting grid with element and vertex parameters
Interval Domain
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 periodic boundaries in x direction periodic 0 % set closure to none for UGGrid closure none # BOUNDARYDOMAIN default 1 % default boundary id #BOUNDARYDOMAIN # examplegrid5.dgf

The resulting grid using SGrid<2,2>

The resulting grid using AlbertaGrid<2,2>
Simplex
Block Simplex
#
Grid Generation in 3d
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

The resulting grid using ALUCubeGrid<3,3>

The resulting grid using ALUSimplexGrid<3,3>
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

The resulting grid using ALUSimplexGrid<3,3>
Now we can also include some quality enhancement:
First: min-angle = 1.2 (remove the first % in the Simplexgenerator block)

The resulting grid using ALUSimplexGrid<3,3>

The resulting grid using ALUSimplexGrid<3,3>
This examples show different ways to define a grid for the following domain and boundary ids:

An example domain
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

The resulting cube grid

The resulting simplex grid
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

The resulting cube grid

The resulting simplex grid
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:

The resulting grid

The resulting grid
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 resulting grid with element parameters
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

The resulting grid with vertex parameters
Importing Grids written in a Tetgen/Triangle format
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

The resulting grid using ALUSimplexGrid<3,3>

The resulting grid using ALUSimplexGrid<3,3>

The resulting grid using ALUSimplexGrid<3,3>

The resulting grid using ALUSimplexGrid<3,3>