[Gmsh] Problems with gmsh 2D quad grids and deal.ii

[Johannes Reinhardt]

Hello everyone,

I encounter problems when trying to load 2D quad meshes generated by
gmsh into a deal.ii based code. I could roughly find out, what the
causes of the problems are, but I don't know whether it is a problem
with gmsh, my use of it, or a problem in deal.ii code. Therefore I have
sent this mail to the mailing lists of both gmsh and deal.ii. The rest
of this rather lengthy mail explains my findings in detail.

I need to calculate a certain quantity as a function of a change in the
geometry. For this reason I want to generate a large number of meshes
using gmsh, so I defined some functions in a common file, and then use
a python script to generate many short .geo files that just include the
common file and call the function with different parameters.

The meshes are 2D quad meshes, because my calculation, that is based on 
deal.ii, which can only deal with quad meshes. The meshes are
relatively detailed already, because I do not use adaptive refinement
and need a good approximation of the curved boundaries of the various
surfaces.

However when I try to load the resulting .msh files in my calculation
(which uses deal.ii), some of the meshes fail to load. Of the
294 meshes I generated, 39 fail.

Of these 31 fail with the error message:

The violated condition was: 
n_negative_cells==0 || n_negative_cells==cells.size()

and 8 fail with the error message:

The violated condition was: 
    needed_lines.find(std::make_pair(line_vertices.second,
    line_vertices.first)) == needed_lines.end()

The first error seems to be connected to by a cell orientation
issue that leads to a invalid cell:
http://www.dealii.org/7.3.0/doxygen/deal.II/classGridReordering.html#ae3a479aa86fabf1f761ab05eb1e42838

The second seems to be due to a vertex ordering problem, I quote from
deal.II/source/grid/tria.cc around line 1640, where the error is thrown:

// assert that the line was
// not already inserted in
// reverse order. This
// happens in spite of the
// vertex rotation above,
// if the sense of the cell
// was incorrect.
//
// Here is what usually
// happened when this
// exception is thrown:
// consider these two cells
// and the vertices
//  3---4---5
//  |   |   |
//  0---1---2
// If in the input vector
// the two cells are given
// with vertices <0 1 4 3>
// and <4 1 2 5>, in the
// first cell the middle
// line would have
// direction 1->4, while in
// the second it would be
// 4->1.  This will cause
// the exception.

The first error seems to be rather common, when looking it up, I found
a tool to postprocess .msh files to resolve it:
https://code.google.com/p/tethex/

When using this tool with my meshes, the first error is indeed avoided,
and of my 294 meshes, but now 9 meshes fail with the second error, so
for one file the first error masked the second one.

I use deal.ii 7.3 and the Linux 64 bit release of gmsh 2.7.1 on Xubuntu
13.04 64bit.

I have attached some files from my investigations:

common_params.geo contains functions that are used to build the
geometry and common parameters.
generate_geo.py is used to generate individual .geos that include
common_params
make_meshes.sh is a shell script that sets up directories, runs the
python scripts, runs gmsh and runs tethex
mass_grid.cc is a small snippet that uses the deal.ii grid classes to
load meshes and logs what is going wrong
failure.log is a log about which .msh files can not be loaded with and
without using tethex and why.

I uploaded a selection of the msh files that should cover all the cases
here (~3MB):

http://ist-dein-freund.de/geometries.tar.gz

Thank you for your efforts in advance

Johannes
-------------- next part --------------
A non-text attachment was scrubbed...
Name: common_params.geo
Type: application/octet-stream
Size: 2886 bytes
Desc: not available
URL: <http://www.geuz.org/pipermail/gmsh/attachments/20130529/93c897a9/attachment.geo>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: failure.log
Type: text/x-log
Size: 2996 bytes
Desc: not available
URL: <http://www.geuz.org/pipermail/gmsh/attachments/20130529/93c897a9/attachment.log>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: generate_geo.py
Type: text/x-python
Size: 853 bytes
Desc: not available
URL: <http://www.geuz.org/pipermail/gmsh/attachments/20130529/93c897a9/attachment.py>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: make_meshes.sh
Type: application/x-shellscript
Size: 181 bytes
Desc: not available
URL: <http://www.geuz.org/pipermail/gmsh/attachments/20130529/93c897a9/attachment.sh>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: mass_grid.cc
Type: text/x-c++src
Size: 3434 bytes
Desc: not available
URL: <http://www.geuz.org/pipermail/gmsh/attachments/20130529/93c897a9/attachment.cc>