fusion-energy / brep_to_h5m

Converts Brep CAD geometry files to h5m geometry files compatible with DAGMC simulations
MIT License
1 stars 2 forks source link

avoid using intermediate STL files #25

Closed shimwell closed 1 year ago

shimwell commented 2 years ago

Currently each volume is surface meshed and saved as an STL file

The mb.load_file("vol1.stl") pymoab command is then used to add the meshed volume to the moab core.

If we could add the mesh directly without writing and reading an STL file that would be a nice upgrade.

I can see from the pymoab source code that there are other methods that look useful such as ``create_elementandadd_entity``` that might allow a moab object to be built up from groups of coordinates

https://bitbucket.org/fathomteam/moab/src/master/pymoab/pymoab/core.pyx

Perhaps there is an option to describe the mesh in terms of groups of coordinates and then load in this information directly

Screenshot from 2022-04-19 11-35-48

For the above diagram with 4 faces and 4 coordinates

I assume we can do something like this

mb = pymoab.core.Core()  # creates an empty moab object

coords = np.array([
                                [
                                     [0,0,0], [1,0,0], [0,1,0],  # bottom face 
                                ],
                                [
                                     [0,0,1], [1,0,0], [0,1,0],  # front face
                                ],
                                [
                                     [0,0,0], [0,1,0], [0,0,1],  # back face
                                ],
                                [
                                     [0,0,0], [1,0,0], [0,0,1],  # left face
                                ],
                              ], dtype='float64')

surface = self.mb.create_meshset()  # makes an empty surface

for three_coord_face in coord:
    verts = mb.create_vertices(three_coord_face)

    tri_verts = np.array(((verts[0],verts[1],verts[2]),),dtype='uint64')

    triangle = mb.create_element(types.MBTRI,tri_verts)

    mb.add_entity(surface, triangle)

volume = self.mb.create_meshset()
# should be parent child
mb.add_entity(volume, surface)

Perhaps this could be reworked into a function that accepts a list of vertices and produces a volume meshset. The volume meshset can then be written to h5m in the same way stl-to-h5m does already

shimwell commented 2 years ago

the gmsh getNodesForPhysicalGroup function can return the coordinates

import paramak

rotated_straights = paramak.RotateStraightShape(
    points=[
        (400, 100),
        (400, 200),
        (600, 200),
        (600, 100)
           ],
    rotation_angle = 180
)
brep_filename='shape.brep'
rotated_straights.export_brep(brep_filename)

import gmsh

gmsh.initialize()
gmsh.option.setNumber("General.Terminal", 1)
gmsh.model.add("made_with_brep_to_h5m_package")
volumes = gmsh.model.occ.importShapes(brep_filename)
gmsh.model.occ.synchronize()

gmsh.option.setNumber("Mesh.Algorithm", 1)
gmsh.option.setNumber("Mesh.MeshSizeMin", 10)
gmsh.option.setNumber("Mesh.MeshSizeMax", 20)
gmsh.model.mesh.generate(2)

vol_id = 1
entities_in_volume = gmsh.model.getAdjacencies(3, vol_id)
surfaces_in_volume = entities_in_volume[1]
ps = gmsh.model.addPhysicalGroup(2, surfaces_in_volume)
gmsh.model.setPhysicalName(2, ps, f"surfaces_on_volume_{vol_id}")
nodeTags, coords=gmsh.model.mesh.getNodesForPhysicalGroup(dim=2,tag=1)

verts_to_h5m()
shimwell commented 2 years ago

suggestions for implementing

tet1 = [1,2,3,4]
tet2 = [4,5,6,7]

surface set for all 7 surfaces
volume set for surface_set 1,2,3,4
volume set for surface_set 4,5,6,7

in this case surface 4 would have dagmc sense set to the reverse for the one of the volumes

shimwell commented 1 year ago

solved by the vertices to h5m package github.com/fusion-energy/vertices_to_h5m