mesh_gen Namespace Reference

Namespaces
namespace	anonymous_namespace{annulusMesh2D.cpp}
namespace	anonymous_namespace{annulusMesh3D.cpp}
namespace	anonymous_namespace{builtinGmshGeometry.cpp}
namespace	anonymous_namespace{meshGenerator.cpp}
namespace	anonymous_namespace{openBoundaryWalls3D.cpp}

Functions
void	buildAnnulus2DInCurrentModel (const geom::AnnulusGeomObject &a, double h)
void	buildAnnulus3DInCurrentModel (const geom::AnnulusGeomObject &a, double h)
void	buildGmshGeometryInCurrentModel (const geom::GeomObject &g, double h)
int	gmshMeshGenerateDim (const inp::ModelDeck *modelDeck)
void	fillMeshFromActiveGmshModel (mesh::Mesh *mesh_p, const inp::MeshDeck *meshDeck, const inp::ModelDeck *modelDeck)
int	addDiskOcc (double xc, double yc, double zc, double rx, double ry, double theta=0.)
	OCC disk compatible with older Gmsh (no zAxis/xAxis) and newer APIs.
void	generateBuiltinParticleMeshGmsh (const std::shared_ptr< geom::GeomObject > &geomObj, double h, const std::string &filenameStem, bool vtk_out, bool write_mesh_file, mesh::Mesh *out_mesh, const inp::MeshDeck *meshDeck, const inp::ModelDeck *modelDeck)
	In-process Gmsh mesh from geom::GeomObject (after createGeomObject on deck / Particle data).
void	physicalGroupsWallOpenFromY2D (int surfaceTag, double yOpen, double tol, const std::string &physWall, const std::string &physOpen)
void	physicalGroupsWallOpenFromFace3D (int volumeTag, int openFace, const util::Point &lo, const util::Point &hi, double t, double tol, const std::string &physWall, const std::string &physOpen)
std::shared_ptr< mesh::Mesh >	createParticleMesh (const inp::MeshDeck &zmeshDeck, const geom::GeomData &zgeomDeck, const inp::ModelDeck *modelDeck, const std::string &modelName)
	Build a reference-particle mesh: file, uniform rectangle, or in-process Gmsh.
void	buildCylinderOcc (const geom::Cylinder &c)
void	buildEllipseOcc (const geom::Ellipse &e, double h)
void	buildEllipsoidOcc (const geom::Ellipsoid &e)

Function Documentation

addDiskOcc()

int mesh_gen::addDiskOcc ( double  xc, double  yc, double  zc, double  rx, double  ry, double  theta = 0.  )

inline

OCC disk compatible with older Gmsh (no zAxis/xAxis) and newer APIs.

Gmsh before ~4.11 only accepts addDisk(xc,yc,zc,rx,ry,tag). Later versions add optional axis vectors for in-plane rotation. On older Gmsh, apply a z-rotation about the disk center when theta is nonzero.

Definition at line 26 of file gmshOccCompat.h.

                                         {
#if (GMSH_API_VERSION_MAJOR > 4) ||                                        \
    (GMSH_API_VERSION_MAJOR == 4 && GMSH_API_VERSION_MINOR >= 11)
  const std::vector<double> zAxis = {0., 0., 1.};
  const std::vector<double> xAxis = {std::cos(theta), std::sin(theta), 0.};
  return gmsh::model::occ::addDisk(xc, yc, zc, rx, ry, -1, zAxis, xAxis);
#else
  const int tag = gmsh::model::occ::addDisk(xc, yc, zc, rx, ry, -1);
  if (std::abs(theta) > 1.0e-15)
    gmsh::model::occ::rotate({{2, tag}}, xc, yc, zc, 0., 0., 1., theta);
  return tag;
#endif
}

Referenced by mesh_gen::anonymous_namespace{annulusMesh2D.cpp}::addEllipseDiskOcc(), mesh_gen::anonymous_namespace{annulusMesh2D.cpp}::buildCircleAnnulus2D(), and buildEllipseOcc().

Here is the caller graph for this function:

buildAnnulus2DInCurrentModel()

void mesh_gen::buildAnnulus2DInCurrentModel	(	const geom::AnnulusGeomObject &	a,
		double	h
	)

OCC: outer rectangle minus inner rectangle (coplanar XY), then mesh-size point embed.

Definition at line 201 of file annulusMesh2D.cpp.

                                                                            {
 
  if (a.d_dim != 2)
    throw std::runtime_error("buildAnnulus2DInCurrentModel: expected d_dim == 2.");
  if (!a.d_inObj_p || !a.d_outObj_p)
    throw std::runtime_error("buildAnnulus2DInCurrentModel: null inner or outer geometry.");
 
  if (a.d_inObj_p->d_name == "rectangle" && a.d_outObj_p->d_name == "rectangle") {
    buildRectangleAnnulus2D(a, h);
    return;
  }
  if (a.d_inObj_p->d_name == "circle" && a.d_outObj_p->d_name == "circle") {
    buildCircleAnnulus2D(a, h);
    return;
  }
  if (a.d_inObj_p->d_name == "ellipse" && a.d_outObj_p->d_name == "ellipse") {
    buildEllipseAnnulus2D(a, h);
    return;
  }
 
  throw std::runtime_error(
      "buildAnnulus2DInCurrentModel: unsupported 2D pair (supported: rectangle−rectangle, "
      "circle−circle, ellipse−ellipse).");
}

References geom::AnnulusGeomObject::d_dim, geom::AnnulusGeomObject::d_inObj_p, geom::GeomObject::d_name, and geom::AnnulusGeomObject::d_outObj_p.

Referenced by buildGmshGeometryInCurrentModel().

Here is the caller graph for this function:

buildAnnulus3DInCurrentModel()

void mesh_gen::buildAnnulus3DInCurrentModel	(	const geom::AnnulusGeomObject &	a,
		double	h
	)

OCC: outer cuboid minus inner cuboid (axis-aligned shell), then mesh-size point embed.

Definition at line 143 of file annulusMesh3D.cpp.

                                                                            {
 
  if (a.d_dim != 3)
    throw std::runtime_error("buildAnnulus3DInCurrentModel: expected d_dim == 3.");
  if (!a.d_inObj_p || !a.d_outObj_p)
    throw std::runtime_error("buildAnnulus3DInCurrentModel: null inner or outer geometry.");
 
  if (a.d_inObj_p->d_name == "cuboid" && a.d_outObj_p->d_name == "cuboid") {
    buildCuboidAnnulus3D(a, h);
    return;
  }
  if (a.d_inObj_p->d_name == "sphere" && a.d_outObj_p->d_name == "sphere") {
    buildSphereAnnulus3D(a, h);
    return;
  }
 
  throw std::runtime_error(
      "buildAnnulus3DInCurrentModel: unsupported 3D pair (supported: cuboid−cuboid, sphere−sphere).");
}

References geom::AnnulusGeomObject::d_dim, geom::AnnulusGeomObject::d_inObj_p, geom::GeomObject::d_name, and geom::AnnulusGeomObject::d_outObj_p.

Referenced by buildGmshGeometryInCurrentModel().

Here is the caller graph for this function:

buildCylinderOcc()

void mesh_gen::buildCylinderOcc

(

const geom::Cylinder &

c

)

OCC cylinder matching geom::Cylinder (base center, axis × length, radius).

Definition at line 23 of file primitiveOccMesh.cpp.

                                             {
 
  if (c.d_r <= 0.)
    throw std::runtime_error("buildCylinderOcc: radius must be positive.");
  if (c.d_l <= 0.)
    throw std::runtime_error("buildCylinderOcc: length must be positive.");
 
  const double dx = c.d_l * c.d_xa.d_x;
  const double dy = c.d_l * c.d_xa.d_y;
  const double dz = c.d_l * c.d_xa.d_z;
 
  gmsh::model::occ::addCylinder(c.d_xBegin.d_x, c.d_xBegin.d_y, c.d_xBegin.d_z, dx, dy, dz, c.d_r);
  gmsh::model::occ::synchronize();
}

References geom::Cylinder::d_l, geom::Cylinder::d_r, util::Point::d_x, geom::Cylinder::d_xa, geom::Cylinder::d_xBegin, util::Point::d_y, and util::Point::d_z.

Referenced by buildGmshGeometryInCurrentModel().

Here is the caller graph for this function:

buildEllipseOcc()

void mesh_gen::buildEllipseOcc	(	const geom::Ellipse &	e,
		double	h
	)

OCC disk / ellipse in the plane z = center.d_z (see geom::Ellipse).

Definition at line 38 of file primitiveOccMesh.cpp.

                                                     {
 
  if (e.d_a <= 0. || e.d_b <= 0.)
    throw std::runtime_error("buildEllipseOcc: semi-axes must be positive.");
 
  const double rx = std::max(e.d_a, e.d_b);
  const double ry = std::min(e.d_a, e.d_b);
  double theta = e.d_theta;
  if (e.d_a < e.d_b)
    theta += M_PI / 2.;
 
  const int s = addDiskOcc(e.d_x.d_x, e.d_x.d_y, e.d_x.d_z, rx, ry, theta);
  gmsh::model::occ::synchronize();
  const int p = gmsh::model::occ::addPoint(e.d_x.d_x, e.d_x.d_y, e.d_x.d_z, h);
  gmsh::model::occ::synchronize();
  gmsh::model::mesh::embed(0, {p}, 2, s);
}

References addDiskOcc(), geom::Ellipse::d_a, geom::Ellipse::d_b, geom::Ellipse::d_theta, geom::Ellipse::d_x, util::Point::d_x, util::Point::d_y, and util::Point::d_z.

Referenced by buildGmshGeometryInCurrentModel().

Here is the call graph for this function:

Here is the caller graph for this function:

buildEllipsoidOcc()

void mesh_gen::buildEllipsoidOcc

(

const geom::Ellipsoid &

e

)

OCC ellipsoid: unit sphere + affine map R diag(r1,r2,r3) and translation.

Definition at line 56 of file primitiveOccMesh.cpp.

                                               {
 
  if (e.d_a <= 0. || e.d_b <= 0. || e.d_c <= 0.)
    throw std::runtime_error("buildEllipsoidOcc: semi-axes must be positive.");
 
  // Unit sphere at origin, then x' = R diag(r1,r2,r3) x + c (row-major 4×4 for gmsh).
  double R[9];
  geom::ellipsoidRotationMatrix(e, R);
 
  std::vector<double> mat(16);
  for (int i = 0; i < 3; ++i) {
    mat[static_cast<size_t>(i * 4 + 0)] = R[i * 3 + 0] * e.d_a;
    mat[static_cast<size_t>(i * 4 + 1)] = R[i * 3 + 1] * e.d_b;
    mat[static_cast<size_t>(i * 4 + 2)] = R[i * 3 + 2] * e.d_c;
    mat[static_cast<size_t>(i * 4 + 3)] = (i == 0) ? e.d_x.d_x : (i == 1) ? e.d_x.d_y : e.d_x.d_z;
  }
  mat[12] = mat[13] = mat[14] = 0.;
  mat[15] = 1.;
 
  const int v = gmsh::model::occ::addSphere(0., 0., 0., 1.);
  gmsh::model::occ::synchronize();
  gmsh::model::occ::affineTransform({{3, v}}, mat);
  gmsh::model::occ::synchronize();
  gmsh::model::occ::removeAllDuplicates();
  gmsh::model::occ::synchronize();
}

References geom::Ellipsoid::d_a, geom::Ellipsoid::d_b, geom::Ellipsoid::d_c, geom::Ellipsoid::d_x, util::Point::d_x, util::Point::d_y, util::Point::d_z, and geom::ellipsoidRotationMatrix().

Referenced by buildGmshGeometryInCurrentModel().

Here is the call graph for this function:

Here is the caller graph for this function:

buildGmshGeometryInCurrentModel()

void mesh_gen::buildGmshGeometryInCurrentModel	(	const geom::GeomObject &	geom,
		double	h
	)

Build OCC/geo entities in the current Gmsh model from a concrete geom::GeomObject.

New shapes: add a branch here or in helpers (e.g. primitiveOccMesh.cpp), and register the name in geom::acceptable_geometries / createGeomObject. AnnulusGeomObject uses d_name annulus_object; meshGenerator also allows that name (deck types are circle_minus_circle / rectangle_minus_rectangle / …). open_rect_channel_2d: U-channel polygon (geom::OpenRectChannel2D) + optional physical groups via openBoundaryWalls2D. open_cuboid_channel_3d: cuboidal shell with one outer face slab removed (geom::OpenCuboidChannel3D) + physical groups via openBoundaryWalls3D.

Definition at line 293 of file builtinGmshGeometry.cpp.

                                                                        {
  const std::string &n = g.d_name;
 
  if (n == "annulus_object") {
    const auto &ag = static_cast<const geom::AnnulusGeomObject &>(g);
    if (ag.d_dim == 2)
      return buildAnnulus2DInCurrentModel(ag, h);
    if (ag.d_dim == 3)
      return buildAnnulus3DInCurrentModel(ag, h);
    throw std::runtime_error("buildGmshGeometryInCurrentModel: annulus_object has invalid d_dim.");
  }
 
  if (n == "circle")
    return buildCircleOcc(static_cast<const geom::Circle &>(g), h);
  if (n == "ellipse")
    return buildEllipseOcc(static_cast<const geom::Ellipse &>(g), h);
  if (n == "sphere")
    return buildSphereOcc(static_cast<const geom::Sphere &>(g));
  if (n == "ellipsoid")
    return buildEllipsoidOcc(static_cast<const geom::Ellipsoid &>(g));
  if (n == "cube")
    return buildOccBoxFromAabb(g);
  if (n == "cuboid")
    return buildOccBoxFromAabb(g);
  if (n == "cylinder")
    return buildCylinderOcc(static_cast<const geom::Cylinder &>(g));
 
  if (n == "square")
    return meshPolygon2DGeoFromVertices(static_cast<const geom::Square &>(g).d_vertices, h);
  if (n == "rectangle")
    return meshPolygon2DGeoFromVertices(static_cast<const geom::Rectangle &>(g).d_vertices, h);
  if (n == "triangle")
    return meshPolygon2DGeoFromVertices(static_cast<const geom::Triangle &>(g).d_vertices, h);
  if (n == "hexagon")
    return meshPolygon2DGeoFromVertices(static_cast<const geom::Hexagon &>(g).d_vertices, h);
  if (n == "drum2d")
    return meshPolygon2DGeoFromVertices(drum2dPolygonBoundary(static_cast<const geom::Drum2D &>(g)),
                                        h);
  if (n == "open_rect_channel_2d")
    return buildOpenRectChannel2DGeo(static_cast<const geom::OpenRectChannel2D &>(g), h);
  if (n == "open_cuboid_channel_3d")
    return buildOpenCuboidChannel3DGeo(static_cast<const geom::OpenCuboidChannel3D &>(g), h);
  if (n == "complex")
    return buildComplexRectangles2DOcc(static_cast<const geom::ComplexGeomObject &>(g), h);
 
  throw std::runtime_error("buildGmshGeometryInCurrentModel: no Gmsh recipe for geometry \"" + n +
                           "\".");
}

References buildAnnulus2DInCurrentModel(), buildAnnulus3DInCurrentModel(), buildCylinderOcc(), buildEllipseOcc(), buildEllipsoidOcc(), and geom::GeomObject::d_name.

Referenced by generateBuiltinParticleMeshGmsh().

Here is the call graph for this function:

Here is the caller graph for this function:

createParticleMesh()

std::shared_ptr< mesh::Mesh > mesh_gen::createParticleMesh	(	const inp::MeshDeck &	zmeshDeck,
		const geom::GeomData &	zgeomDeck,
		const inp::ModelDeck *	modelDeck,
		const std::string &	modelName
	)

Build a reference-particle mesh: file, uniform rectangle, or in-process Gmsh.

Definition at line 40 of file particleMesh.cpp.

                                                                  {
  if (!zmeshDeck.d_createMesh) {
    return std::make_shared<mesh::Mesh>(&zmeshDeck, modelDeck);
  }
 
  if (zmeshDeck.d_createMeshInfo == "uniform" &&
      (zgeomDeck.d_geomName == "rectangle" || zgeomDeck.d_geomName == "cuboid" ||
       zgeomDeck.d_geomName == "square" || zgeomDeck.d_geomName == "cube")) {
 
    if (!zgeomDeck.d_geom_p)
      throw std::runtime_error(
          modelName + ": uniform mesh on a box requires particle geometry object "
                      "(geom::createGeomObject on Particle.Set_i).");
 
    const size_t dim = modelDeck->d_dim;
    const auto bb = zgeomDeck.d_geom_p->box();
    std::pair<std::vector<double>, std::vector<double>> box;
    box.first.reserve(dim);
    box.second.reserve(dim);
    for (size_t i = 0; i < dim; ++i) {
      const double lo_i =
          (i == 0) ? bb.first.d_x : (i == 1) ? bb.first.d_y : bb.first.d_z;
      const double hi_i =
          (i == 0) ? bb.second.d_x : (i == 1) ? bb.second.d_y : bb.second.d_z;
      box.first.push_back(lo_i);
      box.second.push_back(hi_i);
    }
 
    std::vector<size_t> nGrid(dim);
    for (size_t i = 0; i < dim; ++i) {
      const double span = box.second[i] - box.first[i];
      if (span <= 0.)
        throw std::runtime_error(
            modelName + ": uniform mesh: non-positive axis extent from geom box "
                        "(axis " +
            std::to_string(i) + ").");
      nGrid[i] = static_cast<size_t>(span / zmeshDeck.d_hMeshing);
    }
 
    mesh::Mesh temp_mesh;
    mesh::createUniformMesh(&temp_mesh, dim, box, nGrid);
    if (!zmeshDeck.d_voidRegions.empty())
      mesh::removeNodesInBoxes(&temp_mesh, zmeshDeck.d_voidRegions);
    maybeWriteMeshFile(zmeshDeck, temp_mesh);
    return std::make_shared<mesh::Mesh>(temp_mesh);
  }
 
  if (zmeshDeck.d_createMeshInfo == "gmsh_builtin_mesh") {
    if (!zgeomDeck.d_geom_p)
      throw std::runtime_error(
          modelName + ": gmsh_builtin_mesh requires particle geometry object "
                      "(geom::createGeomObject on Particle.Set_i).");
    mesh::Mesh temp_mesh;
    const std::string mesh_stem =
        zmeshDeck.d_filename.empty()
            ? std::string()
            : util::io::removeExtensionFromFile(zmeshDeck.d_filename);
    mesh_gen::generateBuiltinParticleMeshGmsh(
        zgeomDeck.d_geom_p, zmeshDeck.d_hMeshing, mesh_stem, false,
        zmeshDeck.d_writeMeshFile, &temp_mesh, &zmeshDeck, modelDeck);
    return std::make_shared<mesh::Mesh>(temp_mesh);
  }
 
  throw std::runtime_error(
      "Error: Unsupported in-built mesh: CreateMesh.Info = " +
      zmeshDeck.d_createMeshInfo + " with geometry = " + zgeomDeck.d_geomName);
}

References mesh::createUniformMesh(), inp::MeshDeck::d_createMesh, inp::MeshDeck::d_createMeshInfo, inp::ModelDeck::d_dim, inp::MeshDeck::d_filename, geom::GeomData::d_geom_p, geom::GeomData::d_geomName, inp::MeshDeck::d_hMeshing, inp::MeshDeck::d_voidRegions, inp::MeshDeck::d_writeMeshFile, generateBuiltinParticleMeshGmsh(), anonymous_namespace{particleMesh.cpp}::maybeWriteMeshFile(), util::io::removeExtensionFromFile(), and mesh::removeNodesInBoxes().

Referenced by particle::createReferenceParticles().

Here is the call graph for this function:

Here is the caller graph for this function:

fillMeshFromActiveGmshModel()

void mesh_gen::fillMeshFromActiveGmshModel	(	mesh::Mesh *	mesh_p,
		const inp::MeshDeck *	meshDeck,
		const inp::ModelDeck *	modelDeck
	)

Fill PeriDEM Mesh from the active Gmsh model (after mesh::generate).

Definition at line 31 of file gmshMeshPipeline.cpp.

                                                                {
 
  std::vector<std::size_t> nodeTags;
  std::vector<double> coord;
  std::vector<double> paramCoord;
  gmsh::model::mesh::getNodes(nodeTags, coord, paramCoord);
 
  std::unordered_map<std::size_t, std::size_t> tagToIdx;
  tagToIdx.reserve(nodeTags.size());
  for (size_t i = 0; i < nodeTags.size(); ++i)
    tagToIdx[nodeTags[i]] = i;
 
  std::vector<util::Point> nodes(nodeTags.size());
  for (size_t i = 0; i < nodeTags.size(); ++i)
    nodes[i] = util::Point(coord[3 * i], coord[3 * i + 1], coord[3 * i + 2]);
 
  std::vector<int> elementTypes;
  std::vector<std::vector<std::size_t>> elementTags, elementNodeTags;
  gmsh::model::mesh::getElements(elementTypes, elementTags, elementNodeTags, -1, -1);
 
  bool hasTetra = false;
  for (size_t t = 0; t < elementTypes.size(); ++t) {
    if (elementTypes[t] == util::msh_type_tetrahedron) {
      hasTetra = true;
      break;
    }
  }
 
  if (hasTetra) {
    std::vector<size_t> enc;
    for (size_t t = 0; t < elementTypes.size(); ++t) {
      if (elementTypes[t] != util::msh_type_tetrahedron)
        continue;
      const auto &nt = elementNodeTags[t];
      for (size_t j = 0; j < nt.size(); j += 4) {
        enc.push_back(tagToIdx.at(nt[j]));
        enc.push_back(tagToIdx.at(nt[j + 1]));
        enc.push_back(tagToIdx.at(nt[j + 2]));
        enc.push_back(tagToIdx.at(nt[j + 3]));
      }
    }
    if (enc.empty())
      throw std::runtime_error("fillMeshFromActiveGmshModel: no Gmsh tetrahedron elements found.");
    mesh_p->loadFromTetraElements3D(std::move(nodes), std::move(enc), meshDeck, modelDeck);
    return;
  }
 
  std::vector<size_t> enc;
  for (size_t t = 0; t < elementTypes.size(); ++t) {
    if (elementTypes[t] != util::msh_type_triangle)
      continue;
    const auto &nt = elementNodeTags[t];
    for (size_t j = 0; j < nt.size(); j += 3) {
      enc.push_back(tagToIdx.at(nt[j]));
      enc.push_back(tagToIdx.at(nt[j + 1]));
      enc.push_back(tagToIdx.at(nt[j + 2]));
    }
  }
 
  if (enc.empty())
    throw std::runtime_error("fillMeshFromActiveGmshModel: no Gmsh triangle elements found.");
 
  mesh_p->loadFromTriangleElements2D(std::move(nodes), std::move(enc), meshDeck, modelDeck);
}

References mesh::Mesh::loadFromTetraElements3D(), mesh::Mesh::loadFromTriangleElements2D(), util::msh_type_tetrahedron, and util::msh_type_triangle.

Referenced by generateBuiltinParticleMeshGmsh().

Here is the call graph for this function:

Here is the caller graph for this function:

generateBuiltinParticleMeshGmsh()

void mesh_gen::generateBuiltinParticleMeshGmsh	(	const std::shared_ptr< geom::GeomObject > &	geomObj,
		double	h,
		const std::string &	filenameStem,
		bool	vtk_out,
		bool	write_mesh_file,
		mesh::Mesh *	out_mesh,
		const inp::MeshDeck *	meshDeck,
		const inp::ModelDeck *	modelDeck
	)

In-process Gmsh mesh from geom::GeomObject (after createGeomObject on deck / Particle data).

Single entry point: geometry is only the shared object; mesh size, VTK, file stem, and write flags are explicit. Simulation code unpacks MeshDeck + GeomData.d_geom_p at the call site.

Definition at line 66 of file meshGenerator.cpp.

                                                                                                 {
  // Geometry is only `geomObj`; decks supply mesh/model metadata for `out_mesh` fill — no parallel param arrays.
 
  if (!geomObj)
    throw std::runtime_error("generateBuiltinParticleMeshGmsh: GeomObject pointer is null.");
 
  const std::string &geomName = geomObj->d_name;
 
  // AnnulusGeomObject always uses d_name == "annulus_object" (deck may say
  // circle_minus_circle / rectangle_minus_rectangle / etc.).
  bool known = (geomName == "annulus_object");
  if (!known) {
    for (const auto &g : geom::getAcceptableGeometries()) {
      if (g == geomName) {
        known = true;
        break;
      }
    }
  }
  if (!known)
    throw std::runtime_error(
        "generateBuiltinParticleMeshGmsh: geometry type is not in geom::acceptable_geometries.");
 
  gmsh::initialize();
  gmsh::option::setNumber("Mesh.MshFileVersion", 2.2);
  gmsh::clear();
 
  try {
    buildGmshGeometryInCurrentModel(*geomObj, h);
    // OCC-built 3D shapes (e.g. sphere, ellipsoid) do not attach `h` to geometry; without this,
    // Gmsh uses a default size and can produce far too few tets so mesh spacing > PD horizon.
    if (h > 0.) {
      gmsh::option::setNumber("Mesh.MeshSizeMin", h);
      gmsh::option::setNumber("Mesh.MeshSizeMax", h);
    }
    const int genDim = gmshMeshGenerateDim(modelDeck);
    if (genDim == 2) {
      // Frontal-Delaunay (6) is Gmsh's 2D default and, with MeshSizeMin=Max=h,
      // reproduces the v0.1.0 circle disk (hmin ≈ 0.142 mm for h = R/5).
      // Delaunay (5) on the same circle is ~15% finer in min nodal spacing and
      // was only needed for GEO-kernel polygons (frontal can fail on those).
      const std::string &n = geomObj->d_name;
      const bool geo_polygon =
          n == "square" || n == "rectangle" || n == "triangle" || n == "hexagon" ||
          n == "drum2d" || n == "open_rect_channel_2d";
      // complex (OCC cuts) and annuli use Frontal-Delaunay (6)
      gmsh::option::setNumber("Mesh.Algorithm", geo_polygon ? 5 : 6);
    }
    gmsh::model::mesh::generate(genDim);
  } catch (...) {
    gmsh::finalize();
    throw;
  }
 
  if (!write_mesh_file && out_mesh == nullptr)
    throw std::runtime_error(
        "generateBuiltinParticleMeshGmsh: when write_mesh_file is false, out_mesh is required.");
 
  if (out_mesh != nullptr) {
    if (meshDeck == nullptr || modelDeck == nullptr)
      throw std::runtime_error(
          "generateBuiltinParticleMeshGmsh: meshDeck and modelDeck are required when out_mesh is set.");
    fillMeshFromActiveGmshModel(out_mesh, meshDeck, modelDeck);
  }
 
  if (write_mesh_file) {
    if (filenameStem.empty())
      throw std::runtime_error(
          "generateBuiltinParticleMeshGmsh: non-empty filename stem is required when writing a .msh file.");
    const std::string mshPath = filenameStem + ".msh";
    if (modelDeck != nullptr && modelDeck->d_dim == 2 && out_mesh != nullptr)
      writeGmshMsh22From2DTriangleMesh(*out_mesh, mshPath);
    else
      gmsh::write(mshPath);
  }
 
  if (vtk_out) {
    if (filenameStem.empty())
      throw std::runtime_error(
          "generateBuiltinParticleMeshGmsh: non-empty filename stem is required for VTK output.");
    gmsh::write(filenameStem + ".vtk");
  }
 
  gmsh::finalize();
}

References buildGmshGeometryInCurrentModel(), inp::ModelDeck::d_dim, fillMeshFromActiveGmshModel(), geom::getAcceptableGeometries(), and gmshMeshGenerateDim().

Referenced by anonymous_namespace{testMeshGen.cpp}::builtinMeshCheckForGeometry(), createParticleMesh(), getInputJson(), and anonymous_namespace{testMeshGen.cpp}::testCircleMeshWritesFiles().

Here is the call graph for this function:

Here is the caller graph for this function:

gmshMeshGenerateDim()

int mesh_gen::gmshMeshGenerateDim

(

const inp::ModelDeck *

modelDeck

)

Meshing dimension for Gmsh::generate (1–3), aligned with ModelDeck::d_dim.

Definition at line 22 of file gmshMeshPipeline.cpp.

                                                       {
  if (modelDeck == nullptr)
    return 3;
  const int d = static_cast<int>(modelDeck->d_dim);
  if (d < 1 || d > 3)
    return 3;
  return d;
}

References inp::ModelDeck::d_dim.

Referenced by generateBuiltinParticleMeshGmsh().

Here is the caller graph for this function:

physicalGroupsWallOpenFromFace3D()

void mesh_gen::physicalGroupsWallOpenFromFace3D	(	int	volumeTag,
		int	openFace,
		const util::Point &	lo,
		const util::Point &	hi,
		double	t,
		double	tol,
		const std::string &	physWall = "wall",
		const std::string &	physOpen = "open"
	)

After a 3D volume is built in the current Gmsh model, partition its boundary surfaces into physical groups "wall" vs "open" for an open cuboidal shell (one outer face slab removed).

Parameters

volumeTag	Gmsh volume tag (dimension 3).
openFace	Same convention as geom::OpenCuboidChannel3D::d_openFace (0..5).
lo	Outer AABB low corner (world axes).
hi	Outer AABB high corner.
t	Wall thickness (used to locate the opening plane).
tol	Length tolerance for classifying nearly planar faces.
physWall	Physical name for wall surfaces.
physOpen	Physical name for opening / cavity-facing rim surfaces.

Definition at line 32 of file openBoundaryWalls3D.cpp.

                                                                                            {
 
  const double x0 = lo.d_x, y0 = lo.d_y, z0 = lo.d_z;
  const double x1 = hi.d_x, y1 = hi.d_y, z1 = hi.d_z;
  const double Lx = x1 - x0;
  const double Ly = y1 - y0;
  const double Lz = z1 - z0;
  const double horizTol = std::max(tol, 1.0e-9 * std::max({Lx, Ly, Lz, 1.0}));
 
  std::vector<std::pair<int, int>> bnd;
  // Request unoriented boundary entities so tags are always valid entity IDs.
  gmsh::model::getBoundary({{3, volumeTag}}, bnd, false, false);
 
  std::vector<int> wallTags;
  std::vector<int> openTags;
  wallTags.reserve(bnd.size());
  openTags.reserve(bnd.size());
 
  for (const auto &pr : bnd) {
    if (pr.first != 2)
      continue;
    const int surfTag = std::abs(pr.second);
    double xmin = 0., ymin = 0., zmin = 0., xmax = 0., ymax = 0., zmax = 0.;
    gmsh::model::getBoundingBox(pr.first, surfTag, xmin, ymin, zmin, xmax, ymax, zmax);
    const double dx = xmax - xmin;
    const double dy = ymax - ymin;
    const double dz = zmax - zmin;
    const double cx = 0.5 * (xmin + xmax);
    const double cy = 0.5 * (ymin + ymax);
    const double cz = 0.5 * (zmin + zmax);
    const double emin = minEdge2D(dx, dy, dz);
    const double emax = maxEdge2D(dx, dy, dz);
    const bool thin = emin < horizTol * std::max(1.0, 0.1 * emax);
 
    bool isOpen = false;
    if (thin) {
      switch (openFace) {
      case 0: // +x
        isOpen = cx > x1 - t - 2. * horizTol && std::abs(dx - t) < 0.25 * t + horizTol;
        break;
      case 1: // -x
        isOpen = cx < x0 + t + 2. * horizTol && std::abs(dx - t) < 0.25 * t + horizTol;
        break;
      case 2: // +y
        isOpen = cy > y1 - t - 2. * horizTol && std::abs(dy - t) < 0.25 * t + horizTol;
        break;
      case 3: // -y
        isOpen = cy < y0 + t + 2. * horizTol && std::abs(dy - t) < 0.25 * t + horizTol;
        break;
      case 4: // +z — opening at top: rim / horizontal faces near z1 or cavity top
        isOpen = cz > z1 - t - 2. * horizTol && dz < horizTol * 10.;
        break;
      case 5: // -z
        isOpen = cz < z0 + t + 2. * horizTol && dz < horizTol * 10.;
        break;
      default:
        break;
      }
    }
 
    if (isOpen)
      openTags.push_back(surfTag);
    else
      wallTags.push_back(surfTag);
  }
 
  if (!wallTags.empty()) {
    const int g = gmsh::model::addPhysicalGroup(2, wallTags, -1);
    gmsh::model::setPhysicalName(2, g, physWall);
  }
  if (!openTags.empty()) {
    const int g = gmsh::model::addPhysicalGroup(2, openTags, -1);
    gmsh::model::setPhysicalName(2, g, physOpen);
  }
}

References util::Point::d_x, util::Point::d_y, and util::Point::d_z.

Referenced by mesh_gen::anonymous_namespace{builtinGmshGeometry.cpp}::buildOpenCuboidChannel3DGeo().

Here is the caller graph for this function:

physicalGroupsWallOpenFromY2D()

void mesh_gen::physicalGroupsWallOpenFromY2D	(	int	surfaceTag,
		double	yOpen,
		double	tol,
		const std::string &	physWall = "wall",
		const std::string &	physOpen = "open"
	)

After a 2D plane surface is built in the current Gmsh model, partition its boundary curves into physical groups: "wall" vs "open" for a top opening at +y.

Curves classified as "open" are nearly horizontal segments whose elevation is near yOpen (outer top of a closed rectangular annulus). U-channels have no such edge (the opening is a gap), so all curves are typically "wall".

Parameters

surfaceTag	Gmsh surface tag (dimension 2).
yOpen	Reference y-coordinate of the intended open (top) side.
tol	Length tolerance for "horizontal" vs slanted edges.

Definition at line 18 of file openBoundaryWalls2D.cpp.

                                                              {
 
  std::vector<std::pair<int, int>> bnd;
  gmsh::model::getBoundary({{2, surfaceTag}}, bnd, false);
 
  std::vector<int> wallTags;
  std::vector<int> openTags;
  wallTags.reserve(bnd.size());
  openTags.reserve(bnd.size());
 
  const double horizEps = std::max(tol * 1.0e-3, 1.0e-12 * (1.0 + std::abs(yOpen)));
 
  for (const auto &pr : bnd) {
    if (pr.first != 1)
      continue;
    double xmin = 0., ymin = 0., zmin = 0., xmax = 0., ymax = 0., zmax = 0.;
    gmsh::model::getBoundingBox(pr.first, pr.second, xmin, ymin, zmin, xmax, ymax, zmax);
    const double dy = ymax - ymin;
    const double yc = 0.5 * (ymin + ymax);
    const bool nearlyHorizontal = dy < horizEps;
    const bool nearTopOpen = yc >= yOpen - tol && nearlyHorizontal;
    if (nearTopOpen)
      openTags.push_back(pr.second);
    else
      wallTags.push_back(pr.second);
  }
 
  if (!wallTags.empty()) {
    const int g = gmsh::model::addPhysicalGroup(1, wallTags, -1);
    gmsh::model::setPhysicalName(1, g, physWall);
  }
  if (!openTags.empty()) {
    const int g = gmsh::model::addPhysicalGroup(1, openTags, -1);
    gmsh::model::setPhysicalName(1, g, physOpen);
  }
}

Referenced by mesh_gen::anonymous_namespace{builtinGmshGeometry.cpp}::buildOpenRectChannel2DGeo().

Here is the caller graph for this function: