Partition creation of 3D mesh

General
1

Hi, I’m currently creating a 3D grid with three materials, top, middle and bottom. In my code, if it seems impossible to separate three materials cleanly, because it requires quite a lot of grids, is there a better solution?For example, create a grid in different areas.

from fenics import *
from mshr import *
import matplotlib.pyplot as plt

length = 1  
width = 0.5   
height = 1  
radius = 0.2 

domain_1 = Box(Point(0, 0, 0), Point(length, width, 0.7))
domain_2 = Box(Point(0, 0, 0.7), Point(length, width, 0.8))
domain_3 = Box(Point(0, 0, 0.8), Point(length, width, height))
domain = domain_1 + domain_2 + domain_3

cylinder_hole = Cylinder(Point(1, 0, 0.5), Point(1, 0.5, 0.5), radius, radius)

domain_with_hole = domain - cylinder_hole

mesh = generate_mesh(domain_with_hole, 50)

tol = 1E-14
class Omega_0(SubDomain):
    def inside(self, x, on_boundary):
        return x[2] <= 0.7 + tol
    
class Omega_1(SubDomain):
    def inside(self, x, on_boundary):
        return 0.8 + tol >= x[2] >= 0.7 - tol
    
class Omega_2(SubDomain):
    def inside(self, x, on_boundary):
        return 1 + tol >= x[2] >= 0.8 - tol
    
materials = MeshFunction('size_t', mesh, mesh.topology().dim())
subdomain_0 = Omega_0()
subdomain_1 = Omega_1()
subdomain_2 = Omega_2()
subdomain_0.mark(materials, 0)
subdomain_1.mark(materials, 1)
subdomain_2.mark(materials, 2)

print("Number of vertices: ", mesh.num_vertices())
print("Number of cells: ", mesh.num_cells())

#Save materials to VTK file
vtkfile = File("mesh.pvd")
vtkfile << mesh
# Save materials to VTK file
vtkfile = File("materials.pvd")
vtkfile << materials

1709991988779
1709991988779489×515 15 KB

2

Two suggestions:

  1. try with a larger tolerance, e.g. 1e-6
  2. avoid assigning value equal to zero to a subdomain. That is the default value for MeshFunction, and if you use it also in a call to Subdomain.mark you’ll forever be in doubt if a cell was marked as zero by the call or if it was zero because of the default assignment.
3

I understand what you mean, but without adding mesh, it still seems impossible to distinguish different materials clearly on the plane z=0.7 and 0.8, because the mesh is not evenly distributed

4

I wanted a grid like this, so that the different materials could be clearly separated on a certain plane.

1709996005535
1709996005535912×554 159 KB

5

I would suggest that use use GMSH for this. Mshr is long deprecated.

6

If you prefer graphical user interface, try to use Gmsh in FreeCAD, there is a direct mesh export for FEniCS functionality.

Another method is Salome by using for example NetGen inside, export med and then convert to FEniCS, have a look in this post for more details on this: