Dividing RVE Cube based on a Surface Mesh STL : PLC Error
Hello friends of gmsh,
I am currently currently working on a task where I need to create a RVE for FEM simulations. The main issue I am facing is that the boundaries with which I would like to divide the volumetric mesh are already generated using a marching cube algorithm from a vtk file. Now using gmsh, I've been able to create a 2d mesh for the RVE Cube along with the imported stl surface mesh. See image below
However, after this step when I perform the 3d meshing, I am receiving the PLC Error: A segment and a facet intersect at point.
Alternatively I also tried the method given here. However, in this case, I believe STL with a lot of triangles could not be handled by the OCC engine.
Any Ideas on how to proceed? The goal here is to create a 3d Cube FE mesh, with element groups divided by the surface as can be seen in the image above. Hoping to receive a response because I am trying to fix this issue for the past two weeks.
Files:
- Msh file upto 2d meshing : 2d_surf_mesh.msh
- Stl file of graphite surface mesh : mesh_open_4.stl
The code I'm using is:
factory = gmsh.model.geo
gmsh.initialize()
gmsh.model.add('RVE')
gmsh.merge('mesh_open_4.stl')
gmsh.model.mesh.createEdges()
n = gmsh.model.getDimension()
s = gmsh.model.getEntities()
gmsh.fltk.run()
graphite_surface_loop = gmsh.model.geo.addSurfaceLoop([s[i][1] for i in range(len(s))])
graphite_volume_id = gmsh.model.geo.addVolume([graphite_surface_loop], tag=-1)
gmsh.fltk.run()
p0 = factory.addPoint(0,0,0, meshSize=0, tag=-1)
p1 = factory.addPoint(0,0,99, meshSize=0, tag=-1)
p2 = factory.addPoint(0,99,0, meshSize=0, tag=-1)
p3 = factory.addPoint(0,99,99, meshSize=0, tag=-1)
p4 = factory.addPoint(99,0,0, meshSize=0, tag=-1)
p5 = factory.addPoint(99,0,99, meshSize=0, tag=-1)
p6 = factory.addPoint(99,99,0, meshSize=0, tag=-1)
p7 = factory.addPoint(99,99,99, meshSize=0, tag=-1)
factory.synchronize()
l0 = factory.addLine(p0, p4, tag =-1)
l1 = factory.addLine(p4, p6, tag =-1)
l2 = factory.addLine(p6, p2, tag =-1)
l3 = factory.addLine(p2, p0, tag =-1)
l4 = factory.addLine(p1, p5, tag =-1)
l5 = factory.addLine(p5, p7, tag =-1)
l6 = factory.addLine(p7, p3, tag =-1)
l7 = factory.addLine(p3, p1, tag =-1)
l8 = factory.addLine(p0, p1, tag =-1)
l9 = factory.addLine(p4, p5, tag =-1)
l10 = factory.addLine(p6, p7, tag =-1)
l11 = factory.addLine(p2, p3, tag =-1)
factory.synchronize()
c0 = factory.addCurveLoop([l0, l1, l2, l3], tag = -1)
c1 = factory.addCurveLoop([l0, l9, -l4, -l8], tag = -1)
c2 = factory.addCurveLoop([l1, l10, -l5, -l9], tag = -1)
c3 = factory.addCurveLoop([l2, l11, -l6, -l10], tag = -1)
c4 = factory.addCurveLoop([-l3, l11, l7, -l8], tag = -1)
c5 = factory.addCurveLoop([l4, l5, l6, l7], tag = -1)
factory.synchronize()
s0 = factory.addPlaneSurface([c0], tag=-1)
s1 = factory.addPlaneSurface([c1], tag=-1)
s2 = factory.addPlaneSurface([c2], tag=-1)
s3 = factory.addPlaneSurface([c3], tag=-1)
s4 = factory.addPlaneSurface([c4], tag=-1)
s5 = factory.addPlaneSurface([c5], tag=-1)
factory.synchronize()
rve_surface_loop = factory.addSurfaceLoop([s0, s1, s2, s3, s4, s5])
rve_volume_id = factory.addVolume([rve_surface_loop, graphite_surface_loop], tag = -1)
factory.synchronize()
gmsh.model.mesh.embed(2, [1], 3, rve_volume_id)
factory.synchronize()
# graphite_physical_group = gmsh.model.addPhysicalGroup(3, [graphite_volume_id], tag=-1)
# gmsh.model.setPhysicalName(3, graphite_physical_group, 'GRAPHITE')
# ferrite_physical_group = gmsh.model.addPhysicalGroup(3, [ferrite_volume_id], tag=-1)
# gmsh.model.setPhysicalName(3, ferrite_physical_group, 'FERRITE')
gmsh.model.mesh.generate(2)
gmsh.model.mesh.refine()
gmsh.fltk.run()
gmsh.model.mesh.generate(3)