Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

creating a 3d mesh with meshio #1471

Open
Mirialex opened this issue May 30, 2024 · 0 comments
Open

creating a 3d mesh with meshio #1471

Mirialex opened this issue May 30, 2024 · 0 comments

Comments

@Mirialex
Copy link

hi all,

need help in creating 3d mesh with meshio, basically my goal is to create a box with a cylinder inside of it, starting with creating a box
this is what i tried to do
`
import numpy
import meshio
import gmsh
import pygmsh

resolution = 0.01

Channel parameters

L = 2.2
H = 0.41
c = [0.2, 0.2, 0]
r = 0.05

Initialize empty geometry using the built-in kernel in GMSH

geometry = pygmsh.geo.Geometry()
model = geometry.enter()

Add circle

circle = model.add_circle(c, r, mesh_size=resolution)

Add points with finer resolution on the left side

points = [
model.add_point((0, 0, 0), mesh_size=resolution),
model.add_point((L, 0, 0), mesh_size=5 * resolution),
model.add_point((L, H, 0), mesh_size=5 * resolution),
model.add_point((0, H, 0), mesh_size=resolution)
]

Add lines between all points creating the rectangle

channel_lines = [model.add_line(points[i], points[i + 1]) for i in range(-1, len(points) - 1)]

Create a line loop and plane surface for meshing

channel_loop = model.add_curve_loop(channel_lines)
plane_surface = model.add_plane_surface(channel_loop, holes=[circle.curve_loop])

Synchronize the model

model.synchronize()

Mark different boundaries and the volume mesh

model.add_physical([plane_surface], "Volume")
model.add_physical([channel_lines[0]], "Inflow")
model.add_physical([channel_lines[2]], "Outflow")
model.add_physical([channel_lines[1], channel_lines[3]], "Walls")
model.add_physical(circle.curve_loop.curves, "Obstacle")

Generate the mesh and write to file

geometry.generate_mesh(dim=2)
gmsh.write("mesh.msh")
gmsh.clear()
geometry.exit()

Convert the mesh to XDMF format

mesh_from_file = meshio.read("mesh.msh")

def create_mesh(mesh, cell_type, prune_z=False):
cells = mesh.get_cells_type(cell_type)
cell_data = mesh.get_cell_data("gmsh:physical", cell_type)
points = mesh.points[:, :2] if prune_z else mesh.points
out_mesh = meshio.Mesh(points=points, cells={cell_type: cells}, cell_data={"name_to_read": [cell_data]})
return out_mesh

line_mesh = create_mesh(mesh_from_file, "line", prune_z=True)
meshio.write("facet_mesh.xdmf", line_mesh)

triangle_mesh = create_mesh(mesh_from_file, "triangle", prune_z=True)
meshio.write("mesh.xdmf", triangle_mesh)

3D Mesh generation using OpenCASCADE geometry kernel

L = 1.0
H = 1.0
W = 1.0
resolution = 0.1

geom = pygmsh.occ.Geometry()
model3D = geom.enter()

points = [
model3D.add_point([0.0, 0.0, 0.0], mesh_size=resolution),
model3D.add_point([L, 0.0, 0.0], mesh_size=5 * resolution),
model3D.add_point([L, H, 0.0], mesh_size=5 * resolution),
model3D.add_point([0.0, H, 0.0], mesh_size=resolution),
model3D.add_point([0.0, 0.0, W], mesh_size=resolution),
model3D.add_point([L, 0.0, W], mesh_size=5 * resolution),
model3D.add_point([L, H, W], mesh_size=5 * resolution),
model3D.add_point([0.0, H, W], mesh_size=resolution)
]

lines = [model3D.add_line(points[i], points[i + 1])
for i in range(-1, len(points) - 1)
]

Define surfaces (make sure lines are in the correct order to form closed loops)

surfaces = []
surfaces.append(model3D.add_plane_surface(model3D.add_curve_loop([lines[0], lines[1], lines[2], lines[3]]))) # Bottom
surfaces.append(model3D.add_plane_surface(model3D.add_curve_loop([lines[4], lines[5], lines[6], lines[7]]))) # Top
surfaces.append(model3D.add_plane_surface(model3D.add_curve_loop([lines[0], lines[9], lines[4], lines[8]]))) # Front
surfaces.append(model3D.add_plane_surface(model3D.add_curve_loop([lines[1], lines[10], lines[5], lines[9]]))) # Right
surfaces.append(model3D.add_plane_surface(model3D.add_curve_loop([lines[2], lines[11], lines[6], lines[10]]))) # Back
surfaces.append(model3D.add_plane_surface(model3D.add_curve_loop([lines[3], lines[8], lines[7], lines[11]]))) # Left

Define the volume

surface_loop = model3D.add_surface_loop(surfaces)
volume = model3D.add_volume(surface_loop)

Add Physical Groups

model3D.add_physical(volume, label="MyVolume")
for i, surface in enumerate(surfaces):
model3D.add_physical(surface, label=f"MySurface_{i}")

geom.generate_mesh(dim=3)
gmsh.write("mesh3D.msh")
model3D.exit()

mesh3D_from_file = meshio.read("mesh3D.msh")

triangle_mesh = create_mesh(mesh3D_from_file, "triangle", prune_z=True)
meshio.write("facet_mesh3D.xdmf", triangle_mesh)

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

`

but im not sure how to use this add_curve_loop and how to change it in order for it to work
if you have any documention about it, would appreciate

error i got
/usr/bin/python3.10 /home/mirialex/PycharmProjects/fibers/from_community.py
Traceback (most recent call last):
File "/home/mirialex/PycharmProjects/fibers/from_community.py", line 94, in
surfaces.append(model3D.add_plane_surface(model3D.add_curve_loop([lines[0], lines[1], lines[2], lines[3]]))) # Bottom
File "/usr/lib/python3/dist-packages/pygmsh/common/geometry.py", line 80, in add_curve_loop
return CurveLoop(self.env, *args, **kwargs)
File "/usr/lib/python3/dist-packages/pygmsh/common/curve_loop.py", line 26, in init
assert curves[-1].points[-1] == curves[0].points[0]
AssertionError

Process finished with exit code 1

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

No branches or pull requests

1 participant