GMSH guide #875
Labels
Milestone
Comments
|
I haven't had chance to fully test it yet but I think my cube mesh should work - when visualising in Paraview all the markers look to be correct. I'm just having issues getting more complicated meshes to compile. The GUI workflow has issues when it comes to exporting files, but is fairly okay to use to create the base mesh, I think it's to do with what meshio can and can't read and how GMSH exports things when using the actual program. |
Have you tried importing CAD files into GMSH?
What you could do is try to mesh the exact same geometry in the GUI and with the python API. They should produce the exact same .msh files. If the .msh file previously produced by the python API is readable by meshio, then so should the one produced by the GUI. |
|
Not yet. I'm finding that converting from a .msh to an XDMF that retains all the necessary boundaries and tags isn't as straight forward as I thought. The issue with the GUI vs Python API is the method of exporting is that the GUI has problems with the formatting of the file compared to what meshio can read due to version errors |
|
It looks like this answer from @jorgensd works with the latest meshio version for msh-to-xdmf https://fenicsproject.discourse.group/t/transitioning-from-mesh-xml-to-mesh-xdmf-from-dolfin-convert-to-meshio/412/79 Is this what you are using? |
|
Hi,
I’ve been using workflows from this thread but when I then try to use
MeshFromXDMF it has trouble reading the volume and surface files. If I
visualise my surface and volume XDMF files in Paraview however they all
look ok
…On Tue, 13 Aug 2024 at 09:32, Rémi Delaporte-Mathurin < ***@***.***> wrote:
It looks like this answer from @jorgensd <https://github.com/jorgensd>
works with the latest meshio version for msh-to-xdmf
https://fenicsproject.discourse.group/t/transitioning-from-mesh-xml-to-mesh-xdmf-from-dolfin-convert-to-meshio/412/79
Is this what you are using?
—
Reply to this email directly, view it on GitHub
<#875 (comment)>,
or unsubscribe
<https://github.com/notifications/unsubscribe-auth/BJPSYFVROM3EQ5WWIAYES6LZRHADHAVCNFSM6AAAAABMMJEWNOVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDEOBVGY3TGOJWHE>
.
You are receiving this because you commented.Message ID:
***@***.***>
|
|
@celyneira I managed to get the following to work: Mesh file (taken from the thread) Conversion and FESTIM sim import meshio
import numpy as np
import festim as F
# convert mesh to xdmf
msh = meshio.read("mesh.msh")
for cell in msh.cells:
if cell.type == "triangle":
triangle_cells = cell.data
elif cell.type == "tetra":
tetra_cells = cell.data
for key in msh.cell_data_dict["gmsh:physical"].keys():
if key == "triangle":
triangle_data = msh.cell_data_dict["gmsh:physical"][key]
elif key == "tetra":
tetra_data = msh.cell_data_dict["gmsh:physical"][key]
tetra_mesh = meshio.Mesh(
points=msh.points,
cells={"tetra": tetra_cells},
cell_data={"f": [tetra_data]},
)
triangle_mesh = meshio.Mesh(
points=msh.points,
cells=[("triangle", triangle_cells)],
cell_data={"f": [triangle_data]},
)
print("surface ids ", np.unique(triangle_data))
print("volume ids ", np.unique(tetra_data))
meshio.write("mesh.xdmf", tetra_mesh)
meshio.write("mf.xdmf", triangle_mesh)
# FESTIM simulation
my_model = F.Simulation()
my_model.mesh = F.MeshFromXDMF(volume_file="mesh.xdmf", boundary_file="mf.xdmf")
my_model.boundary_conditions = [
F.DirichletBC(surfaces=[40], value=10, field=0),
F.DirichletBC(surfaces=[10], value=2, field=0),
F.DirichletBC(surfaces=[20], value=5, field=0),
]
my_model.materials = F.Material(id=50, D_0=1, E_D=0)
my_model.T = 500
my_model.settings = F.Settings(
absolute_tolerance=1e-10,
relative_tolerance=1e-10,
transient=False,
)
my_model.exports = [F.XDMFExport(field="solute")]
my_model.initialise()
my_model.run()It produces the following field:
|
|
Please note that you can use |
|
Yep will do that for FESTIM2 but we are still trying to get a workflow for legacy-fenics. But I got it to work!!! Meshing script: import gmsh
def main():
gmsh.initialize()
# Define the first material cube
cube1 = gmsh.model.occ.addBox(x=0, y=0, z=0, dx=1, dy=1, dz=1)
# Define the second material cube
# Adjust the position to make the cubes touch each other
cube2 = gmsh.model.occ.addBox(x=1, y=0, z=0, dx=1, dy=1, dz=1)
gmsh.model.occ.synchronize()
# Perform a Boolean union to merge the cubes and their shared surfaces
gmsh.model.occ.fragment([(3, cube1)], [(3, cube2)])
gmsh.model.occ.synchronize()
# Create physical groups
gmsh.model.addPhysicalGroup(3, [cube1], 1)
gmsh.model.addPhysicalGroup(3, [cube2], 2)
# Tagging the surfaces
# Get the surfaces of cube 1 and cube 2
surfaces_cube1 = gmsh.model.getBoundary(
[(3, cube1)], oriented=False, recursive=False
)
surfaces_cube2 = gmsh.model.getBoundary(
[(3, cube2)], oriented=False, recursive=False
)
# Tag the left surface of cube 1
left_surface_cube1 = None
for surface in surfaces_cube1:
com = gmsh.model.occ.getCenterOfMass(surface[0], surface[1])
if com[0] == 0:
left_surface_cube1 = surface[1]
break
# Tag the right surface of cube 2
right_surface_cube2 = None
for surface in surfaces_cube2:
com = gmsh.model.occ.getCenterOfMass(surface[0], surface[1])
if com[0] == 2:
right_surface_cube2 = surface[1]
break
if left_surface_cube1 is not None:
gmsh.model.addPhysicalGroup(
2, [left_surface_cube1], 3
) # Tag ID 3 for left surface of cube 1
if right_surface_cube2 is not None:
gmsh.model.addPhysicalGroup(
2, [right_surface_cube2], 4
) # Tag ID 4 for right surface of cube 2
# Generate the mesh
gmsh.model.mesh.generate(3)
# Save the mesh
gmsh.write("two_material_cubes_touching.msh")
# If you want to visualize the mesh
gmsh.fltk.run()
gmsh.finalize()
if __name__ == "__main__":
main()Conversion and FESTIM: import meshio
import numpy as np
import festim as F
# Convert mesh to XDMF
filename = "two_material_cubes_touching.msh"
msh = meshio.read(filename)
# Initialize lists to store cells and their corresponding data
triangle_cells_list = []
tetra_cells_list = []
triangle_data_list = []
tetra_data_list = []
# Extract cell data for all types
for cell in msh.cells:
if cell.type == "triangle":
triangle_cells_list.append(cell.data)
elif cell.type == "tetra":
tetra_cells_list.append(cell.data)
# Extract physical tags
for key, data in msh.cell_data_dict["gmsh:physical"].items():
if key == "triangle":
triangle_data_list.append(data)
elif key == "tetra":
tetra_data_list.append(data)
# Concatenate all tetrahedral cells and their data
tetra_cells = np.concatenate(tetra_cells_list)
tetra_data = np.concatenate(tetra_data_list)
# Concatenate all triangular cells and their data
triangle_cells = np.concatenate(triangle_cells_list)
triangle_data = np.concatenate(triangle_data_list)
# Create the tetrahedral mesh
tetra_mesh = meshio.Mesh(
points=msh.points,
cells=[("tetra", tetra_cells)],
cell_data={"f": [tetra_data]},
)
# Create the triangular mesh for the surface
triangle_mesh = meshio.Mesh(
points=msh.points,
cells=[("triangle", triangle_cells)],
cell_data={"f": [triangle_data]},
)
# Print unique surface and volume IDs
print("Surface IDs: ", np.unique(triangle_data))
print("Volume IDs: ", np.unique(tetra_data))
# Write the mesh files
meshio.write("mesh.xdmf", tetra_mesh)
meshio.write("mf.xdmf", triangle_mesh)
# FESTIM simulation
my_model = F.Simulation()
my_model.mesh = F.MeshFromXDMF(volume_file="mesh.xdmf", boundary_file="mf.xdmf")
my_model.boundary_conditions = [
F.DirichletBC(surfaces=[3], value=10, field=0),
F.DirichletBC(surfaces=[4], value=1, field=0),
]
my_model.materials = [
F.Material(id=1, D_0=1, E_D=0),
F.Material(id=2, D_0=4, E_D=0),
]
my_model.T = 500
my_model.settings = F.Settings(
absolute_tolerance=1e-10,
relative_tolerance=1e-10,
transient=False,
)
my_model.exports = [F.XDMFExport(field="solute")]
my_model.initialise()
my_model.run()It produces this concentration field: |
|
ok yes your second method works for me, I think I had an issue with the first method due to there being two volumes? I haven't had a lot of time today to look at this but I'll continue to mess with things |
Merged
10 tasks
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Following #739 we would need to add a guide/demo for a GMSH workflow.
@celyn has a working example for a cube mesh with the python API of GMSH.
I don't know if the GUI workflow is needed.
I'm opening this issue to start the discussion
The text was updated successfully, but these errors were encountered: