Create a circular pattern axisymetric mesh, with GMSH

38 views Asked by At

I'm cumming from this issue (https://gitlab.onelab.info/gmsh/gmsh/-/issues/1283) and I use the https://gitlab.onelab.info/gmsh/gmsh/-/blob/master/examples/api/mirror_mesh.py example to perform my circular pattern.

My goal is to create a axi-symetric rotatory machine mesh by copying its blade mesh N times as required. The problem here is that the common part of the mesh (blade base), seems to be a superposed unconnected meshes and not a single mesh. How can I do to re-mesh that part of the mesh to connect all blades bases, or some kind of mesh boolean union operation?, even if I loose the mesh symetry in that mesh part.

I provide a "simple" working example of my code:

The geometry file used: blade.geo

Python code:

import math

import gmsh
import numpy as np


def rotate3D(p, v, angle):
    angle = np.deg2rad(angle)
    v = v / np.linalg.norm(v)
    R = np.array(
        [
            [
                np.cos(angle) + v[0] ** 2 * (1 - np.cos(angle)),
                v[0] * v[1] * (1 - np.cos(angle)) - v[2] * np.sin(angle),
                v[0] * v[2] * (1 - np.cos(angle)) + v[1] * np.sin(angle),
            ],
            [
                v[1] * v[0] * (1 - np.cos(angle)) + v[2] * np.sin(angle),
                np.cos(angle) + v[1] ** 2 * (1 - np.cos(angle)),
                v[1] * v[2] * (1 - np.cos(angle)) - v[0] * np.sin(angle),
            ],
            [
                v[2] * v[0] * (1 - np.cos(angle)) - v[1] * np.sin(angle),
                v[2] * v[1] * (1 - np.cos(angle)) + v[0] * np.sin(angle),
                np.cos(angle) + v[2] ** 2 * (1 - np.cos(angle)),
            ],
        ]
    )
    return (*np.dot(R, p),)


def rotate_mesh(
    ref: int,
    entities,
    offset_entity,
    offset_node,
    offset_element,
    ax=(1, 0, 0),
    angle=180,
):
    offset_entity = offset_entity * ref
    offset_node = offset_node * ref
    offset_element = offset_element * ref
    m = {}
    for e in entities:
        bnd = gmsh.model.getBoundary([e])
        nod = gmsh.model.mesh.getNodes(e[0], e[1])
        ele = gmsh.model.mesh.getElements(e[0], e[1])
        m[e] = (bnd, nod, ele)

    for e in sorted(m):
        gmsh.model.addDiscreteEntity(
            e[0],
            e[1] + offset_entity,
            [(abs(b[1]) + offset_entity) * math.copysign(1, b[1]) for b in m[e][0]],
        )
        coord = []
        for i in range(0, len(m[e][1][1]), 3):
            x = m[e][1][1][i]
            y = m[e][1][1][i + 1]
            z = m[e][1][1][i + 2]
            x, y, z = rotate3D((x, y, z), ax, angle)
            coord.append(x)
            coord.append(y)
            coord.append(z)
        gmsh.model.mesh.addNodes(
            e[0], e[1] + offset_entity, m[e][1][0] + offset_node, coord
        )
        gmsh.model.mesh.addElements(
            e[0],
            e[1] + offset_entity,
            m[e][2][0],
            [t + offset_element for t in m[e][2][1]],
            [n + offset_node for n in m[e][2][2]],
        )

# ROTOR SETTINGS
N_BLADE = 3
GEO_FILE = "blade.geo"

# MESHING
gmsh.initialize()
gmsh.merge(GEO_FILE)

# get the Physical Groups before the circular pattern
physical_groups = gmsh.model.getPhysicalGroups()
physical_groups_names = [gmsh.model.getPhysicalName(*g) for g in physical_groups]
physical_groups_entities = [
    gmsh.model.getEntitiesForPhysicalGroup(*g) for g in physical_groups
]

# perform circular pattern
entities = gmsh.model.getEntities()
offset_entity = len(entities)
offset_node = gmsh.model.mesh.getMaxNodeTag()
offset_element = gmsh.model.mesh.getMaxElementTag()
for i in range(1, N_BLADE):
    ang = i * 360 / N_BLADE
    rotate_mesh(i, entities, offset_entity, offset_node, offset_element, (1, 0, 0), ang)

# add Physical Groups to the new mesh
gmsh.model.remove_physical_groups(physical_groups)
for group in zip(physical_groups, physical_groups_names, physical_groups_entities):
    dim = group[0][0]
    tags = group[2]
    name = group[1]
    gmsh.model.add_physical_group(
        dim,
        tags.tolist()
        + [s + offset_entity * i for s in tags for i in range(1, N_BLADE)],
        name=name,
    )

gmsh.model.mesh.removeDuplicateElements()
gmsh.model.mesh.removeDuplicateNodes()
gmsh.model.mesh.reclassifyNodes()
gmsh.model.mesh.optimize()

gmsh.write("rotor.msh")
gmsh.write("rotor.stl")
gmsh.write("rotor.vtk")

gmsh.fltk.run()

gmsh.finalize()
0

There are 0 answers