nschloe / meshio

:spider_web: input/output for many mesh formats
MIT License
1.96k stars 402 forks source link

Import of GMSH file with curvelinear boundaries #955

Open awolf-github opened 4 years ago

awolf-github commented 4 years ago

Firts of all thank you for the meshio-project.

I'm using gmsh version 4.6.0 to generate a mesh from a step file (gmsh -3 modell.step) modell.zip

The error accours if I read the mesh, generate a copy of the mesh by using meshio and try to read the copy again.

import meshio
msh = meshio.read("modell.msh")
meshio.write("modell2.msh",msh)
msh2 = meshio.read("modell2.msh")
nschloe commented 4 years ago

The zip contains a step file, not a msh.

awolf-github commented 4 years ago

Yes. But gmsh can deal with it

nschloe commented 4 years ago

The code you pasted doesn't use step files.

awolf-github commented 4 years ago

In firts step I've used gmsh -3 modell.step to generate the mesh file in bash. The output is modell.msh file which is then used in python

nschloe commented 4 years ago

The please attach the mesh, not the step file.

awolf-github commented 4 years ago

See attachment modell.zip

gdmcbain commented 4 years ago

I confirm the issue (with meshio 4.3.2 from conda-forge). The line

meshio.write("modell2.msh", msh)

raises

WARNING:root:Binary Gmsh needs c_size_t (got int64). Converting.
WARNING:root:Binary Gmsh needs c_size_t (got int64). Converting.
WARNING:root:Binary Gmsh needs c_size_t (got int64). Converting.
WARNING:root:Binary Gmsh needs c_size_t (got int64). Converting.
WARNING:root:Binary Gmsh needs c_size_t (got int64). Converting.
WARNING:root:Binary Gmsh needs c_size_t (got int64). Converting.
WARNING:root:Binary Gmsh needs c_size_t (got int64). Converting.
WARNING:root:Binary Gmsh needs c_size_t (got int64). Converting.
WARNING:root:Binary Gmsh needs c_size_t (got int64). Converting.
WARNING:root:Binary Gmsh needs c_size_t (got int64). Converting.
WARNING:root:Binary Gmsh needs c_size_t (got int64). Converting.
WARNING:root:Binary Gmsh needs c_size_t (got int64). Converting.
WARNING:root:Binary Gmsh needs c_size_t (got int64). Converting.
WARNING:root:Binary Gmsh needs c_size_t (got int64). Converting.
Traceback (most recent call last):
  File "/home/gmcbain/projects/meshio/955-conda/msmdir.015615/mwe.py", line 3, in <module>
    meshio.write("modell2.msh", msh)
  File "/home/gmcbain/miniconda3/envs/meshio/lib/python3.9/site-packages/meshio/_helpers.py", line 146, in write
    return writer(filename, mesh, **kwargs)
  File "/home/gmcbain/miniconda3/envs/meshio/lib/python3.9/site-packages/meshio/gmsh/main.py", line 111, in <lambda>
    "gmsh": lambda f, m, **kwargs: write(f, m, "4.1", **kwargs),
  File "/home/gmcbain/miniconda3/envs/meshio/lib/python3.9/site-packages/meshio/gmsh/main.py", line 102, in write
    writer.write(filename, mesh, binary=binary, float_fmt=float_fmt)
  File "/home/gmcbain/miniconda3/envs/meshio/lib/python3.9/site-packages/meshio/gmsh/_gmsh41.py", line 306, in write
    write4_1(filename, mesh, float_fmt=float_fmt, binary=binary)
  File "/home/gmcbain/miniconda3/envs/meshio/lib/python3.9/site-packages/meshio/gmsh/_gmsh41.py", line 345, in write4_1
    _write_entities(
  File "/home/gmcbain/miniconda3/envs/meshio/lib/python3.9/site-packages/meshio/gmsh/_gmsh41.py", line 514, in _write_entities
    physical_tag = cell_data["gmsh:physical"][matching_cell_block[0]][0]
KeyError: 'gmsh:physical'

The last is unfamiliar. Clearly no MSH 4.1 file generated by gmsh -3 *.step is going to have any physical entities, but none were defined and this had not been an issue quite recently. I see that meshio 4.3.1 is the same but back in meshio 4.3.0, it's

    raise WriteError("Can only deal with one cell type for now")
meshio._exceptions.WriteError: Can only deal with one cell type for now

which is (or was) a known issue for writing MSH 4.1, a duplicate of #524, #614, ….

gdmcbain commented 4 years ago

While there might be a regression in meshio 4.3.1 reading Gmsh MSH 4.1 without physical entities #956, there might be an easy remedy in your case: just define one.

I was initially surprised to read the command-line gmsh -3 modell.step, I didn't know that that would have worked at all. I was more used to using Merge, as in the Gmsh GEO code

Merge "modell.step";

extending this with

Merge "modell.step";
Physical Volume(1) = {1};

gives a MSH 4.1 that can be read, rewritten, reread by meshio without a hitch: this mesh has only one cell type (so avoiding #524) but does have at least one physical entity (so avoiding #956).

awolf-github commented 4 years ago

By using the modell.geo file I still get an error message (see below). I'm using gmsh version 4.6.0 to generate the mesh gmsh -3 modell.geo from attached files modell.zip

Python 3.7.9 (default, Aug 31 2020, 12:42:55) 
[GCC 7.3.0] :: Anaconda, Inc. on linux
Type "help", "copyright", "credits" or "license" for more information.
import meshio
meshio.__version__
'4.3.2'

msh = meshio.read('modell.msh')
meshio.write('modell2.msh',msh)
WARNING:root:Binary Gmsh needs c_size_t (got int64). Converting.
>>> msh2 = meshio.read('modell2.msh')
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/home/y1alwo/miniconda3/envs/local_fem/lib/python3.7/site-packages/meshio/_helpers.py", line 68, in read
    return reader_map[file_format](filename)
  File "/home/y1alwo/miniconda3/envs/local_fem/lib/python3.7/site-packages/meshio/gmsh/main.py", line 19, in read
    mesh = read_buffer(f)
  File "/home/y1alwo/miniconda3/envs/local_fem/lib/python3.7/site-packages/meshio/gmsh/main.py", line 48, in read_buffer
    return reader.read_buffer(f, is_ascii, data_size)
  File "/home/y1alwo/miniconda3/envs/local_fem/lib/python3.7/site-packages/meshio/gmsh/_gmsh41.py", line 82, in read_buffer
    field_data,
  File "/home/y1alwo/miniconda3/envs/local_fem/lib/python3.7/site-packages/meshio/gmsh/_gmsh41.py", line 238, in _read_elements
    pt = None if not physical_tags else physical_tags[dim][tag]
KeyError: 1
awolf-github commented 4 years ago

I've tested simple geometries like a cube, a sphere and a wedge and for all of them meshio could generate a copy of the mesh. Since sphere has a curvelinear boundary the title of the issue is not correct. simple_meshes.zip

gdmcbain commented 4 years ago
msh = meshio.read('modell.msh')
meshio.write('modell2.msh',msh)
WARNING:root:Binary Gmsh needs c_size_t (got int64). Converting.
>>> msh2 = meshio.read('modell2.msh')
Traceback (most recent call last):

Right, meshio 4.3.2 successfully reads modell.msh which contains 'gmsh:physical' and 459 tetra and no cells of other type, but then writes it out again with only an irrelevant warning but the output modell2.msh cannot be reread, raising the KeyError: 1, which is different from the original KeyError: 'gmsh:physical' which I've attempted to characterize as #956.

The new modell2.msh is also rejected by Gmsh 4.6.0:

Info    : Running 'gmsh modell2.msh' [Gmsh 4.6.0, 1 node, max. 1 thread]
Info    : Started on Fri Oct 30 12:55:49 2020
Info    : Reading 'modell2.msh'...
Info    : 13 entities
Info    : 170 nodes
Info    : 459 elements
Error   : Unknown entity 1 of dimension 3
Error   : Could not read elements
Error   : Error loading 'modell2.msh'
Info    : Done reading 'modell2.msh'

This is bad; meshio should not write MSH 4.1 files that are rejected by Gmsh. I don't know what's happening here. Maybe something introduced in #922?

This is another regression, trying again with meshio 4.3.0

gmsh -3 modell.geo
conda install -c conda-forge meshio=4.3.0

the above works fine:

import meshio
msh = meshio.read("modell.msh")
meshio.write("modell2.msh", msh)

both modell.msh and modell2.msh are accepted by Gmsh.

gdmcbain commented 4 years ago

I confirm that the cube, sphere, and wedge can be read, rewritten, and reread successfully by meshio 4.3.2; there's something different about the MSH 4.1 that Gmsh generates from modell.step.

awolf-github commented 4 years ago

I've used meshio-convert command this time to generate the mesh. meshio-convert --ascii modell.msh modell2.msh and compared both meshes. It seems that meshio does not recognize the $Entities correctly Here is the original mesh

MeshFormat
4.1 0 8
$EndMeshFormat
$Entities
4 6 4 1
1 13 -3.184081677783e-15 9 0 
2 13 -3.184081677783e-15 4 0 
3 8 -1.959434878636e-15 4 0 
4 6.123233995737e-16 -1.499759782662e-31 8 0 
1 -13.0000001 -13.0000001 8.999999900000001 13.0000001 13.0000001 9.000000099999999 0 2 1 -1 
2 12.9999999 -1.000000031840817e-07 3.9999999 13.0000001 9.999999681591831e-08 9.000000099999999 0 2 2 -1 
3 -13.0000001 -13.0000001 3.9999999 13.0000001 13.0000001 4.0000001 0 2 2 -2 
4 -8.000000099999999 -8.000000099999999 3.9999999 8.000000099999999 8.000000099999999 4.0000001 0 2 3 -3 
5 -9.99999993922529e-08 -1.000000019594349e-07 3.9999999 8.000000099999999 1e-07 8.000000099999999 0 2 3 -4 
6 -1.000041676398853e-07 -1.00004162951209e-07 7.999999899999971 1.000041676398853e-07 1.00004162951209e-07 8.000000100000028 0 2 4 -4 
1 -13.0000001 -13.0000001 3.9999999 13.0000001 13.0000001 9.000000099999999 0 4 1 -2 3 2 
2 -13.0000001 -13.0000001 8.999999900000001 13.0000001 13.0000001 9.000000099999999 0 1 1 
3 -13.0000001 -13.0000001 3.9999999 13.0000001 13.0000001 4.0000001 0 2 3 4 
4 -8.000000100002483 -8.000000100002483 3.999999899996689 8.000000100002483 8.000000100002483 8.000000099999999 0 4 6 -5 4 5 
1 -13.0000001 -13.0000001 3.999999899996689 13.0000001 13.0000001 9.000000099999999 1 1 4 1 2 3 4 
$EndEntities
$Nodes
14 170 1 170

and here the modell2.msh

$MeshFormat
4.1 0 8
$EndMeshFormat
$Entities
4 5 4 0
1 0 0 0 0 
2 0 0 0 0 
3 0 0 0 0 
4 0 0 0 0 
1 0 0 0 0 0 0 0 0
2 0 0 0 0 0 0 0 0
3 0 0 0 0 0 0 0 0
4 0 0 0 0 0 0 0 0
5 0 0 0 0 0 0 0 0
1 0 0 0 0 0 0 0 0
2 0 0 0 0 0 0 0 0
3 0 0 0 0 0 0 0 0
4 0 0 0 0 0 0 0 0
$EndEntities
$Nodes
13 170 1 170

By the way, I'm using now the old msh2.2 format in my application and can therefore avoid the error. From my point of view we can close the issue.

gdmcbain commented 4 years ago

Wow, that doesn't look good. Thank you for persisting with the investigation after your own needs have been met. I don't understand the MSH 4.1 $Entities myself but it is a shame to increasingly often have to recommend the legacy MSH 2.2. I would close this issue and launch a new one but I'm not able to characterize the problem well enough to give it a title. I've made a couple of attempts at a MWE but like your cube, wedge, & sphere, they've worked. As you say, it's not the curvature of the boundaries, and I think it's different from #956 but like it possibly brought in by #922.