What is to date the best way to create a 3D mesh for use in FEniCS?

mesh
1

Dear all,
I’ve been using mshr to create (with success) 2D meshes. However to create a 3D version of the same mesh I’m struggling to join a sphere with a box (i.e.: there is no flow of heat from one portion of the mesh to the other). It was pointed out to me that mshr is currently not maintained so I’m looking into alternatives. What is to date the best workflow to create complex 3D meshes for use in FEniCS?

Thank you for your help!

2

I would suggest looking at Gmsh and/or Netgen/NGSolve. There are several other quality mesh generators too. The above two have a very nice python API. I would have suggested using pygmsh (pre 7.0) but it uses gmsh python API as of the 7.0 release, so not much of a difference, but feel free to check it out too.

1 Like
3

To add to @bhaveshshrimali’s answer:
I have made tutorials on how to use pygmsh (Pre version 7) and the Gmsh Python api, which can be found here.

2 Likes
4

Thank you both very much! I’ll look into these resources.

5

I forgot to mention I’ve been using an anaconda environment for FEniCS. Would similar steps work in there or it’s better to switch to Docker altogether?

6

Simulator steps would work with anaconda. You just have to install the dependencies (Gmsh, meshio, pygmsh)

8

Thanks! I was able to follow your tutorial successfully in a conda environment. I’m only having some issues with getting the ball to join the box. I assume this is coming from geom.add_physical(union, 1) since my sphere is only meant to join the box partially. Is there any alternative or do I have to go with gmsh GUI?

9

Please make a minimal example that illustrates your issue.

10

I would also recommend Salome, which is very similar to CAD programs such as FreeCad. If you use this however, you will need to use gmsh to convert the .med file to a .msh

1 Like
11

Probably you can use meshio to do the conversion.

1 Like
12 
import pygmsh.opencascade as opencascade
from pygmsh import generate_mesh
import meshio
import numpy as np


def create_mesh(mesh, cell_type):
    cells = np.vstack([cell.data for cell in mesh.cells if cell.type==cell_type])
    data = np.hstack([mesh.cell_data_dict["gmsh:physical"][key]
                       for key in mesh.cell_data_dict["gmsh:physical"].keys() if key==cell_type])
    mesh = meshio.Mesh(points=mesh.points, cells={cell_type: cells},
                           cell_data={"name_to_read":[data]})
    return mesh

rbottom = [0.0, 0.0, 0.0]                                          
rtop = [50.0, 50.0, 2.0]
sphereRadius = 3.0
sphereCenter = [25, 25, 4.45]

char_length = 0.1

geom = opencascade.Geometry()

box =  geom.add_box(rbottom, rtop, char_length=char_length)

ball = geom.add_ball(sphereCenter, sphereRadius, char_length=char_length)

union = geom.boolean_union([box, ball])

geom.add_physical(union, 1)

mesh3D = generate_mesh(geom, geo_filename="mesh3D.geo", msh_filename="mesh3D.msh", verbose=False)

tetra_mesh = create_mesh(mesh3D, "tetra")
meshio.write("mesh3D.xdmf", tetra_mesh)

This is pretty much a copy of your tutorial, it’s the closest I got to make it work.

13

Thanks, I’ll have a look into that!

14

The behavior of pygmsh seems really weird to me (The resolution of the ball seems to be ignored). Have you tried doing the same operations in the GMSH API?

15

No, I haven’t. I’ll try with the API or maybe the GUI.