diff --git a/Applications/Utils/SimpleMeshCreation/generateStructuredMesh.cpp b/Applications/Utils/SimpleMeshCreation/generateStructuredMesh.cpp index 1bfad11a73415ce9e622ebffc21e14e50d24fac7..3e70e9f24e5e0b693e665fba69139351b29a5b1f 100644 --- a/Applications/Utils/SimpleMeshCreation/generateStructuredMesh.cpp +++ b/Applications/Utils/SimpleMeshCreation/generateStructuredMesh.cpp @@ -79,9 +79,10 @@ int main (int argc, char* argv[]) allowed_ele_types.push_back("tri"); allowed_ele_types.push_back("quad"); allowed_ele_types.push_back("hex"); + allowed_ele_types.push_back("tet"); TCLAP::ValuesConstraint<std::string> allowedVals(allowed_ele_types); TCLAP::ValueArg<std::string> eleTypeArg("e", "element-type", - "element type to be created: line | tri | quad | hex", true, "line", &allowedVals); + "element type to be created: line | tri | quad | hex | tet", true, "line", &allowedVals); cmd.add(eleTypeArg); TCLAP::ValueArg<std::string> mesh_out("o", "mesh-output-file", "the name of the file the mesh will be written to", true, @@ -207,6 +208,9 @@ int main (int argc, char* argv[]) case MeshLib::MeshElemType::HEXAHEDRON: mesh.reset(MeshLib::MeshGenerator::generateRegularHexMesh(*vec_div[0], *vec_div[1], *vec_div[2])); break; + case MeshLib::MeshElemType::TETRAHEDRON: + mesh.reset(MeshLib::MeshGenerator::generateRegularTetMesh(*vec_div[0], *vec_div[1], *vec_div[2])); + break; default: ERR("Given element type is not supported."); break; diff --git a/MeshLib/MeshGenerators/MeshGenerator.cpp b/MeshLib/MeshGenerators/MeshGenerator.cpp index 388775d83ef3faf79faf7e622e27bdcaf357ce25..87fddf6c106b6bf21d459099bba29eb9b6dc9395 100644 --- a/MeshLib/MeshGenerators/MeshGenerator.cpp +++ b/MeshLib/MeshGenerators/MeshGenerator.cpp @@ -16,6 +16,7 @@ #include "MeshLib/Elements/Line.h" #include "MeshLib/Elements/Quad.h" #include "MeshLib/Elements/Hex.h" +#include "MeshLib/Elements/Tet.h" #include "MeshLib/Elements/Tri.h" #include "MeshLib/MeshEditing/AddLayerToMesh.h" #include "MeshLib/MeshEditing/RemoveMeshComponents.h" @@ -460,6 +461,143 @@ Mesh* MeshGenerator::generateRegularPrismMesh( return MeshLib::removeElements(*mesh, elem_ids, mesh_name); } +Mesh* MeshGenerator::generateRegularTetMesh( + const double x_length, + const double y_length, + const double z_length, + const std::size_t x_subdivision, + const std::size_t y_subdivision, + const std::size_t z_subdivision, + const MathLib::Point3d& origin, + std::string const& mesh_name) +{ + return generateRegularTetMesh(x_subdivision, y_subdivision, z_subdivision, + x_length/x_subdivision, y_length/y_subdivision, z_length/z_subdivision, + origin, mesh_name); +} + +Mesh* MeshGenerator::generateRegularTetMesh( + const unsigned n_x_cells, + const unsigned n_y_cells, + const unsigned n_z_cells, + const double cell_size_x, + const double cell_size_y, + const double cell_size_z, + MathLib::Point3d const& origin, + std::string const& mesh_name) +{ + return generateRegularTetMesh( + BaseLib::UniformSubdivision(n_x_cells*cell_size_x, n_x_cells), + BaseLib::UniformSubdivision(n_y_cells*cell_size_y, n_y_cells), + BaseLib::UniformSubdivision(n_z_cells*cell_size_z, n_z_cells), + origin, mesh_name); +} + +Mesh* MeshGenerator::generateRegularTetMesh( + const BaseLib::ISubdivision &div_x, + const BaseLib::ISubdivision &div_y, + const BaseLib::ISubdivision &div_z, + MathLib::Point3d const& origin, + std::string const& mesh_name) +{ + std::vector<double> vec_x(div_x()); + std::vector<double> vec_y(div_y()); + std::vector<double> vec_z(div_z()); + std::vector<Node*> nodes(generateRegularNodes(vec_x, vec_y, vec_z, origin)); + + const unsigned n_x_nodes (vec_x.size()); + const unsigned n_y_nodes (vec_y.size()); + const unsigned n_x_cells (vec_x.size()-1); + const unsigned n_y_cells (vec_y.size()-1); + const unsigned n_z_cells (vec_z.size()-1); + + //elements + std::vector<Element*> elements; + elements.reserve(n_x_cells * n_y_cells * n_z_cells * 6); + + for (std::size_t i = 0; i < n_z_cells; i++) + { + const std::size_t offset_z1 = i * n_x_nodes * n_y_nodes; // bottom + const std::size_t offset_z2 = (i + 1) * n_x_nodes * n_y_nodes; // top + for (std::size_t j = 0; j < n_y_cells; j++) + { + const std::size_t offset_y1 = j * n_x_nodes; + const std::size_t offset_y2 = (j + 1) * n_x_nodes; + for (std::size_t k = 0; k < n_x_cells; k++) + { + // tet 1 + { + std::array<Node*, 4> element_nodes; + // bottom + element_nodes[0] = nodes[offset_z1 + offset_y1 + k]; + element_nodes[1] = nodes[offset_z1 + offset_y2 + k + 1]; + element_nodes[2] = nodes[offset_z1 + offset_y2 + k]; + // top + element_nodes[3] = nodes[offset_z2 + offset_y1 + k]; + elements.push_back (new Tet(element_nodes)); + } + // tet 2 + { + std::array<Node*, 4> element_nodes; + // bottom + element_nodes[0] = nodes[offset_z1 + offset_y2 + k + 1]; + element_nodes[1] = nodes[offset_z1 + offset_y2 + k]; + // top + element_nodes[2] = nodes[offset_z2 + offset_y1 + k]; + element_nodes[3] = nodes[offset_z2 + offset_y2 + k + 1]; + elements.push_back (new Tet(element_nodes)); + } + // tet 3 + { + std::array<Node*, 4> element_nodes; + // bottom + element_nodes[0] = nodes[offset_z1 + offset_y2 + k]; + // top + element_nodes[1] = nodes[offset_z2 + offset_y1 + k]; + element_nodes[2] = nodes[offset_z2 + offset_y2 + k + 1]; + element_nodes[3] = nodes[offset_z2 + offset_y2 + k]; + elements.push_back (new Tet(element_nodes)); + } + // tet 4 + { + std::array<Node*, 4> element_nodes; + // bottom + element_nodes[0] = nodes[offset_z1 + offset_y1 + k]; + element_nodes[1] = nodes[offset_z1 + offset_y1 + k + 1]; + element_nodes[2] = nodes[offset_z1 + offset_y2 + k + 1]; + // top + element_nodes[3] = nodes[offset_z2 + offset_y1 + k + 1]; + elements.push_back (new Tet(element_nodes)); + } + // tet 5 + { + std::array<Node*, 4> element_nodes; + // bottom + element_nodes[0] = nodes[offset_z1 + offset_y1 + k]; + element_nodes[1] = nodes[offset_z1 + offset_y2 + k + 1]; + // top + element_nodes[2] = nodes[offset_z2 + offset_y1 + k]; + element_nodes[3] = nodes[offset_z2 + offset_y1 + k + 1]; + elements.push_back (new Tet(element_nodes)); + } + // tet 6 + { + std::array<Node*, 4> element_nodes; + // bottom + element_nodes[0] = nodes[offset_z1 + offset_y2 + k + 1]; + // top + element_nodes[1] = nodes[offset_z2 + offset_y1 + k]; + element_nodes[2] = nodes[offset_z2 + offset_y1 + k + 1]; + element_nodes[3] = nodes[offset_z2 + offset_y2 + k + 1]; + elements.push_back (new Tet(element_nodes)); + } + } + } + } + + return new Mesh(mesh_name, nodes, elements); +} + MeshLib::Mesh* MeshGenerator::createSurfaceMesh(std::string const& mesh_name, MathLib::Point3d const& ll, MathLib::Point3d const& ur, diff --git a/MeshLib/MeshGenerators/MeshGenerator.h b/MeshLib/MeshGenerators/MeshGenerator.h index a08a747e427a0f47f0ff8c5dac728ede5fc00ba8..b10b114a61a0eb81f25e0e981b48dffae57acf8d 100644 --- a/MeshLib/MeshGenerators/MeshGenerator.h +++ b/MeshLib/MeshGenerators/MeshGenerator.h @@ -421,6 +421,69 @@ Mesh* generateRegularPrismMesh(const unsigned n_x_cells, MathLib::Point3d const& origin = MathLib::ORIGIN, std::string const& mesh_name = "mesh"); +/** + * Generate a regular 3D Tet-Element mesh. + * + * This algorithm generates regular grid points and split each grid cell into six tetrahedrals. + * + * \param div_x Subdivision operator in x direction + * \param div_y Subdivision operator in y direction + * \param div_z Subdivision operator in z direction + * \param origin Optional mesh's origin (the bottom lower left corner) with MathLib::ORIGIN default. + * \param mesh_name Name of the new mesh. + */ +Mesh* generateRegularTetMesh(const BaseLib::ISubdivision &div_x, + const BaseLib::ISubdivision &div_y, + const BaseLib::ISubdivision &div_z, + MathLib::Point3d const& origin = MathLib::ORIGIN, + std::string const& mesh_name = "mesh"); + +/** + * Generate a regular 3D Tet-Element mesh. + * + * This algorithm generates regular grid points and split each grid cell into six tetrahedrals. + * + * \param x_length Mesh dimension in x-direction. + * \param y_length Mesh dimension in y-direction. + * \param z_length Mesh dimension in z-direction. + * \param x_subdivision Number of subdivisions in x-direction. + * \param y_subdivision Number of subdivisions in y-direction. + * \param z_subdivision Number of subdivisions in z-direction. + * \param origin Optional mesh's origin (the bottom lower left corner) with MathLib::ORIGIN default. + * \param mesh_name Name of the new mesh. + */ +Mesh* generateRegularTetMesh(const double x_length, + const double y_length, + const double z_length, + const std::size_t x_subdivision, + const std::size_t y_subdivision, + const std::size_t z_subdivision, + MathLib::Point3d const& origin = MathLib::ORIGIN, + std::string const& mesh_name = "mesh"); + +/** + * Generate a regular 3D Tet-Element mesh. + * + * This algorithm generates regular grid points and split each grid cell into six tetrahedrals. + * + * \param n_x_cells Number of cells in x-direction. + * \param n_y_cells Number of cells in y-direction. + * \param n_z_cells Number of cells in z-direction. + * \param cell_size_x Edge length of Tet elements in x-direction. + * \param cell_size_y Edge length of Tet elements in y s-direction. + * \param cell_size_z Edge length of Tet elements in z-direction + * \param origin Optional mesh's origin (the lower left corner) with MathLib::ORIGIN default. + * \param mesh_name Name of the new mesh. + */ +Mesh* generateRegularTetMesh(const unsigned n_x_cells, + const unsigned n_y_cells, + const unsigned n_z_cells, + const double cell_size_x, + const double cell_size_y, + const double cell_size_z, + MathLib::Point3d const& origin = MathLib::ORIGIN, + std::string const& mesh_name = "mesh"); + /// Constructs a surface mesh approximating a surface in the 3d space given by /// a function. /// The surface within the xy-domain \f$[ll[0], ur[0]] \times [ll[1], ur[1]\f$