Also, FYI, this is my effort towards generating a mesh for the problem stated in #1588
Some points, which are used to generate the curve are very slightly off the plane (4e-2 error)
For the plane shown in the image above the points from which the CurveLoops are generated are as follows:
Side.ZX_PLUS
Closing Surface Curve for ZX_PLUS
Curve Points:
[[ 7.85122204 98.99953461 15.79696178]
[10.50467396 98.96764374 18.242239 ]
[11.13861942 99. 21. ]
[11.19631767 99. 22. ]
[10. 99. 24.06012344]
[ 8.81718636 99. 25.70767975]
[ 7.2372365 98.9969635 28.11213493]
[ 5. 99. 30.43367195]
[ 3.0911572 98.98747253 30.77775955]
[ 0. 99. 29.88615227]
[ 0. 99. 99. ]
[15.5063715 99. 99. ]
[15.45397758 99. 98. ]
[15.9881649 99. 96.98873901]
[18. 99. 96.31934357]
[19. 99. 96.7681427 ]
[19.44369698 99. 99. ]
[99. 99. 99. ]
[99. 99. 82.83720398]
[98.31617737 99. 85.5 ]
[97.46531677 98.99739838 86.51204681]
[95.24176788 98.99988556 87.8412323 ]
[91. 99. 89.32992554]
[89.69467926 98.98996735 88.77893829]
[85.61408234 98.99964905 85.97666931]
[83.31599426 99. 82.88822937]
[80.98503113 99. 80.02915192]
[80.56729889 99. 77. ]
[84.14086914 99. 73.1754837 ]
[87.97346497 98.99175262 69.81996918]
[92.1038208 98.98532867 68.32574463]
[95.54121399 99. 70.15272522]
[98.10284424 98.99835968 74.36509705]
[99. 99. 75.42941284]
[99. 99. 0. ]
[75.53888702 99. 0. ]
[75.51417542 99. 1.5 ]
[75.08796692 99. 4.91558027]
[74. 99. 6.96033144]
[73. 99. 7.47651482]
[71.95001221 99. 7. ]
[71.49875641 99. 3. ]
[70.37236023 99. 2. ]
[70.15048981 99. 1. ]
[70.42929077 99. 0. ]
[ 0. 99. 0. ]
[ 0. 99. 15.47007656]
[ 1. 99. 15.21899414]
[ 5.74206448 98.97620392 14.74265862]]
Boundary Curve 9
Curve Points:
[[15. 99. 22.24702454]
[16. 99. 22.19179726]
[16.99629974 99. 23.00666809]
[17.33810616 99. 24. ]
[16.54361153 99. 26. ]
[16.17150879 99. 31. ]
[14. 99. 32.59423065]
[13. 99. 32.23978424]
[11.55433369 99. 30.5 ]
[11.93755245 99. 25.62646484]
[14.04186153 99. 23.0480175 ]]
Boundary Curve 14
Curve Points:
[[38. 99. 28.3058033 ]
[39. 99. 28.35775566]
[41.41404724 99. 32.14810562]
[41.08234024 99. 35. ]
[39.8121376 99. 36.16167068]
[38. 99. 36.6192894 ]
[36.88687897 99. 36.11312103]
[35.34084702 99. 34. ]
[37.02362061 99. 29.02130318]]
Boundary Curve 16
Curve Points:
[[47.03936005 99. 40.03447723]
[51. 99. 38.49594879]
[52.04640579 99. 39. ]
[53.83101273 99. 42. ]
[53.39136124 99. 45.27084351]
[51.16548157 99. 49. ]
[50. 99. 49.36804199]
[47.57462311 98.99409485 47.95490646]
[44. 99. 45.6236496 ]
[43.78265381 99. 44. ]]
Additionally for Side YZ_Minus, the same issue occurs(not shown in images above)
Side.YZ_MINUS
Closing Surface Curve for YZ_MINUS
Curve Points:
[[0.00000000e+00 9.90000000e+01 0.00000000e+00]
[0.00000000e+00 9.90000000e+01 1.54700766e+01]
[2.03750879e-02 9.70255966e+01 1.63947048e+01]
[0.00000000e+00 9.37203140e+01 1.90000000e+01]
[0.00000000e+00 9.34260101e+01 2.05000000e+01]
[0.00000000e+00 9.37338409e+01 2.50000000e+01]
[0.00000000e+00 9.47352982e+01 2.69131451e+01]
[0.00000000e+00 9.72763748e+01 2.93700123e+01]
[0.00000000e+00 9.90000000e+01 2.98861523e+01]
[0.00000000e+00 9.90000000e+01 9.90000000e+01]
[0.00000000e+00 0.00000000e+00 9.90000000e+01]
[0.00000000e+00 0.00000000e+00 0.00000000e+00]]
Boundary Curve 2
Curve Points:
[[0.00000000e+00 3.00000000e+01 7.21812897e+01]
[0.00000000e+00 2.83452549e+01 7.12331314e+01]
[0.00000000e+00 2.55975952e+01 6.88677063e+01]
[0.00000000e+00 2.50059986e+01 6.60000000e+01]
[0.00000000e+00 2.80000000e+01 6.34713058e+01]
[3.84895341e-03 3.22473373e+01 6.43671646e+01]
[1.99757033e-05 3.72287598e+01 6.51782684e+01]
[0.00000000e+00 4.17634621e+01 6.66981277e+01]
[0.00000000e+00 4.38651924e+01 6.86455841e+01]
[0.00000000e+00 4.43637009e+01 7.05000000e+01]
[0.00000000e+00 4.29237442e+01 7.39319611e+01]
[1.17071001e-02 3.81151123e+01 7.46424713e+01]
[0.00000000e+00 3.41208916e+01 7.35708313e+01]]
Boundary Curve 4
Curve Points:
[[ 0. 81. 39.7156601 ]
[ 0. 82. 40.23688507]
[ 0. 83.29307556 42.5 ]
[ 0. 84.3600769 47. ]
[ 0. 84.08524323 48. ]
[ 0. 83. 48.82769012]
[ 0. 82. 48.82404327]
[ 0. 80.91921997 48.08078766]
[ 0. 80. 47.60902405]
[ 0. 78. 47.86508179]
[ 0. 76.63981628 47. ]
[ 0. 76.21346283 42. ]
[ 0. 79.22987366 39.88877106]]
Dear friends of GMSH,
I am currently using GMSH to generate a mesh for finite element analysis. Unfortunately, when meshing a plane surface (that is created from a curveloop, which is created from set of ordered points) the mesh is not on the plane, or has an offset angle (see image below). This does not however always occur (see second image).
Please ignore the cut-outs at the edges of the plane surface. They are part of a bigger 3d model that I am trying to setup like shown below
Is there anyway that this can be fixed, or am I doing something wrong here?
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:
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)
Thank you for your response Prof. Geuzaine. I used the following modified code to solve my problem
import gmsh
import numpy as np
factory = gmsh.model.geo
gmsh.initialize()
gmsh.option.setNumber("General.Verbosity", 2)
gmsh.model.add("RVE")
gmsh.merge("graphite_boundary_mesh_optimized_copy.stl")
# the STL is just a bunch of triangles - create (discrete) surfaces, curves and
# points to have full topology
gmsh.model.mesh.createTopology()
# create curve loops with the curves
holes = factory.addCurveLoops([c[1] for c in gmsh.model.getEntities(1)])
graphite_boundary = [c[1] for c in gmsh.model.getEntities(2)]
factory.synchronize()
print("Holes: ", holes)
hole_surfaces = []
for i in range(len(holes)):
hole_surface = factory.addPlaneSurface([holes[i]])
hole_surfaces.append(hole_surface)
factory.synchronize()
sl2 = factory.addSurfaceLoop(np.append(graphite_boundary, [h for h in hole_surfaces]))
gmsh.model.mesh.setCompound(2, np.append(graphite_boundary, [h for h in hole_surfaces]))
# Create a volume from all the surfaces
v2 = factory.addVolume([sl2])
factory.synchronize()
gmsh.fltk.run()
# create other cad entities
p0 = factory.addPoint(0, 0, 0)
p1 = factory.addPoint(0, 0, 99)
p2 = factory.addPoint(0, 99, 0)
p3 = factory.addPoint(0, 99, 99)
p4 = factory.addPoint(99, 0, 0)
p5 = factory.addPoint(99, 0, 99)
p6 = factory.addPoint(99, 99, 0)
p7 = factory.addPoint(99, 99, 99)
l0 = factory.addLine(p0, p4)
l1 = factory.addLine(p4, p6)
l2 = factory.addLine(p6, p2)
l3 = factory.addLine(p2, p0)
l4 = factory.addLine(p1, p5)
l5 = factory.addLine(p5, p7)
l6 = factory.addLine(p7, p3)
l7 = factory.addLine(p3, p1)
l8 = factory.addLine(p0, p1)
l9 = factory.addLine(p4, p5)
l10 = factory.addLine(p6, p7)
l11 = factory.addLine(p2, p3)
c0 = factory.addCurveLoop([l0, l1, l2, l3])
c1 = factory.addCurveLoop([l0, l9, -l4, -l8])
c2 = factory.addCurveLoop([l1, l10, -l5, -l9])
c3 = factory.addCurveLoop([l2, l11, -l6, -l10])
c4 = factory.addCurveLoop([-l3, l11, l7, -l8])
c5 = factory.addCurveLoop([l4, l5, l6, l7])
s0 = factory.addPlaneSurface([c0])
s1 = factory.addPlaneSurface([c1])
s2 = factory.addPlaneSurface([c2])
s3 = factory.addPlaneSurface(np.append(c3, [h for h in holes])) # hybrid!
s4 = factory.addPlaneSurface([c4])
s5 = factory.addPlaneSurface([c5])
factory.synchronize()
gmsh.fltk.run()
# create the (hybrid) volume
sl1 = factory.addSurfaceLoop([s[1] for s in gmsh.model.getEntities(2)])
v1 = factory.addVolume([sl1])
factory.synchronize()
gmsh.model.mesh.generate(2)
gmsh.fltk.run()
gmsh.model.mesh.generate(3)
gmsh.fltk.run()
However this works for the specific case that I had posted before where there are arbitrary number of holes only on the top surface. Is there a way to generalize the procedure for when there are holes on all six surfaces and further even when topological holes are across two surfaces of the cube? For example look at the image below:
My idea is to identify the hole curves that belong to a particular closing surface of the cube and create hybrid planes for all the six sides of the cube with the corresponding hole curves. However I cannot think of any good solution for when the hole curves are across two sides of the cube..
Any help would be appreciated.
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:
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)