VerticalSliceFromLayers.cpp File Reference

Detailed Description

Copyright

Copyright (c) 2012-2024, OpenGeoSys Community (http://www.opengeosys.org) Distributed under a Modified BSD License. See accompanying file LICENSE.txt or http://www.opengeosys.org/project/license

Definition in file VerticalSliceFromLayers.cpp.

#include <algorithm>
#include <cmath>
#include <memory>
#include <string>
#include <vector>
#include <tclap/CmdLine.h>
#include <QCoreApplication>
#include "Applications/FileIO/Gmsh/GMSHInterface.h"
#include "Applications/FileIO/Gmsh/GmshReader.h"
#include "BaseLib/FileTools.h"
#include "BaseLib/IO/readStringListFromFile.h"
#include "BaseLib/MPI.h"
#include "GeoLib/AABB.h"
#include "GeoLib/AnalyticalGeometry.h"
#include "GeoLib/GEOObjects.h"
#include "GeoLib/IO/XmlIO/Qt/XmlGmlInterface.h"
#include "GeoLib/Point.h"
#include "GeoLib/Polygon.h"
#include "GeoLib/Polyline.h"
#include "InfoLib/GitInfo.h"
#include "MathLib/MathTools.h"
#include "MathLib/Point3d.h"
#include "MeshGeoToolsLib/GeoMapper.h"
#include "MeshLib/Elements/Element.h"
#include "MeshLib/IO/VtkIO/VtuInterface.h"
#include "MeshLib/IO/readMeshFromFile.h"
#include "MeshLib/Mesh.h"
#include "MeshLib/MeshSearch/ElementSearch.h"
#include "MeshLib/Node.h"
#include "MeshToolsLib/MeshEditing/MeshRevision.h"
#include "MeshToolsLib/MeshEditing/RemoveMeshComponents.h"
#include "MeshToolsLib/MeshSurfaceExtraction.h"

Include dependency graph for VerticalSliceFromLayers.cpp:

Go to the source code of this file.

Functions
std::vector< GeoLib::Point * >	createPoints (MathLib::Point3d const &start, MathLib::Point3d const &end, std::size_t const n_intervals)
	creates a vector of sampling points based on the specified resolution
GeoLib::Polyline *	createPolyline (std::vector< GeoLib::Point * > const &points)
	creates a polyline to be mapped on a mesh layer
std::vector< std::string >	createGeometries (GeoLib::GEOObjects &geo, std::vector< std::string > const &layer_names, MathLib::Point3d const &pnt_start, MathLib::Point3d const &pnt_end, double const resolution)
	creates a mapped line for each of the mesh layers
static GeoLib::Point *	getMinElevationPoint (GeoLib::Point *a, GeoLib::Point *b)
void	mergeGeometries (GeoLib::GEOObjects &geo, std::vector< std::string > const &geo_names, std::string &merged_geo_name)
std::pair< Eigen::Matrix3d, double >	rotateGeometryToXY (std::vector< GeoLib::Point * > &points)
	rotates the merged geometry into the XY-plane
void	consolidateGeometry (GeoLib::GEOObjects &geo, std::string const &output_name, std::string &merged_geo_name, bool const keep_gml_file)
MeshLib::Mesh *	generateMesh (GeoLib::GEOObjects &geo, std::string const &geo_name, std::string const &output_name, double res)
	converts geometry into GMSH format and creates mesh
void	rotateMesh (MeshLib::Mesh &mesh, Eigen::Matrix3d const &rot_mat, double const z_shift)
	inverse rotation of the mesh, back into original position
MeshLib::Mesh *	removeLineElements (MeshLib::Mesh const &mesh)
	removes line elements from mesh such that only triangles remain
void	writeBoundary (MeshLib::Mesh const &mesh, std::vector< std::size_t > const &idx_array, std::string const &file_name)
void	extractBoundaries (MeshLib::Mesh const &mesh, std::string const &output_name, MathLib::Point3d const &pnt_start)
int	main (int argc, char *argv[])

Function Documentation

consolidateGeometry()

void consolidateGeometry	(	GeoLib::GEOObjects &	geo,
		std::string const &	output_name,
		std::string &	merged_geo_name,
		bool const	keep_gml_file )

This encapsulates a workaround: For unknown reasons, the GML->GEO converter will not work correctly when inputting the merged geometry directly. However, if the geometry is saved to a file and immediately loaded again, everything works fine.

Definition at line 208 of file VerticalSliceFromLayers.cpp.

{
    std::string const filename(output_name + ".gml");
    GeoLib::IO::XmlGmlInterface xml(geo);
    xml.export_name = merged_geo_name;
    BaseLib::IO::writeStringToFile(xml.writeToString(), filename);
 
    geo.removePolylineVec(merged_geo_name);
    geo.removePointVec(merged_geo_name);
 
    xml.readFile(filename);
 
    if (!keep_gml_file)
    {
        BaseLib::removeFile(filename);
        BaseLib::removeFile(filename + ".md5");
    }
}

References BaseLib::IO::XMLInterface::export_name, GeoLib::IO::XmlGmlInterface::readFile(), BaseLib::removeFile(), GeoLib::GEOObjects::removePointVec(), GeoLib::GEOObjects::removePolylineVec(), BaseLib::IO::writeStringToFile(), and BaseLib::IO::Writer::writeToString().

Referenced by main().

createGeometries()

std::vector< std::string > createGeometries	(	GeoLib::GEOObjects &	geo,
		std::vector< std::string > const &	layer_names,
		MathLib::Point3d const &	pnt_start,
		MathLib::Point3d const &	pnt_end,
		double const	resolution )

creates a mapped line for each of the mesh layers

Definition at line 82 of file VerticalSliceFromLayers.cpp.

{
    std::vector<std::string> geo_name_list;
    std::size_t const n_layers = layer_names.size();
    for (std::size_t i = 0; i < n_layers; ++i)
    {
        std::unique_ptr<MeshLib::Mesh> const layer(
            MeshLib::IO::readMeshFromFile(layer_names[i]));
        if (layer == nullptr)
        {
            ERR("Could not read file {:s}. Skipping layer...", layer_names[i]);
            continue;
        }
        if (layer->getDimension() != 2)
        {
            ERR("Layer {:d} is not a 2D mesh. Skipping layer...", i);
            continue;
        }
 
        std::string geo_name(std::to_string(i));
        auto points(createPoints(pnt_start, pnt_end, resolution));
        geo.addPointVec(std::move(points), geo_name,
                        GeoLib::PointVec::NameIdMap{});
 
        std::vector<GeoLib::Polyline*> lines;
        lines.push_back(createPolyline(*geo.getPointVec(geo_name)));
        geo.addPolylineVec(std::move(lines), geo_name,
                           GeoLib::PolylineVec::NameIdMap{});
 
        MeshGeoToolsLib::GeoMapper mapper(geo, geo_name);
        mapper.mapOnMesh(layer.get());
        geo_name_list.push_back(geo_name);
    }
    return geo_name_list;
}

References GeoLib::GEOObjects::addPointVec(), GeoLib::GEOObjects::addPolylineVec(), createPoints(), createPolyline(), ERR(), GeoLib::GEOObjects::getPointVec(), MeshGeoToolsLib::GeoMapper::mapOnMesh(), and MeshLib::IO::readMeshFromFile().

Referenced by main().

createPoints()

std::vector< GeoLib::Point * > createPoints	(	MathLib::Point3d const &	start,
		MathLib::Point3d const &	end,
		std::size_t const	n_intervals )

creates a vector of sampling points based on the specified resolution

Definition at line 48 of file VerticalSliceFromLayers.cpp.

{
    std::vector<GeoLib::Point*> points;
    points.push_back(new GeoLib::Point(start, 0));
    double const length = std::sqrt(MathLib::sqrDist(start, end));
    double const interval = length / n_intervals;
 
    GeoLib::Point const vec((end[0] - start[0]), (end[1] - start[1]),
                            (end[2] - start[2]));
    for (std::size_t i = 1; i < n_intervals; ++i)
    {
        points.push_back(new GeoLib::Point(
            start[0] + ((i * interval) / length * vec[0]),
            start[1] + ((i * interval) / length * vec[1]), 0, i));
    }
    points.push_back(new GeoLib::Point(end, n_intervals));
    return points;
}

References MathLib::sqrDist().

Referenced by createGeometries().

createPolyline()

GeoLib::Polyline * createPolyline

(

std::vector< GeoLib::Point * > const &

points

)

creates a polyline to be mapped on a mesh layer

Definition at line 70 of file VerticalSliceFromLayers.cpp.

{
    GeoLib::Polyline* line = new GeoLib::Polyline(points);
    std::size_t const length = points.size();
    for (std::size_t i = 0; i < length; ++i)
    {
        line->addPoint(i);
    }
    return line;
}

References GeoLib::Polyline::addPoint().

Referenced by createGeometries().

extractBoundaries()

void extractBoundaries	(	MeshLib::Mesh const &	mesh,
		std::string const &	output_name,
		MathLib::Point3d const &	pnt_start )

Definition at line 305 of file VerticalSliceFromLayers.cpp.

{
    auto const edge_length = minMaxEdgeLength(mesh.getElements());
    double const eps = edge_length.first / 100.0;
    std::unique_ptr<MeshLib::Mesh> boundary_mesh(
        MeshToolsLib::BoundaryExtraction::getBoundaryElementsAsMesh(
            mesh,
            MeshLib::getBulkIDString(MeshLib::MeshItemType::Node),
            MeshLib::getBulkIDString(MeshLib::MeshItemType::Cell),
            MeshLib::getBulkIDString(MeshLib::MeshItemType::Face)));
 
    auto const& elems = boundary_mesh->getElements();
    std::vector<std::size_t> left_bound_idx, right_bound_idx, top_bound_idx,
        bottom_bound_idx;
    Eigen::Vector2d const anchor(pnt_start[0], pnt_start[1]);
    for (auto e : elems)
    {
        Eigen::Vector2d const n1((*e->getNode(0))[0], (*e->getNode(0))[1]);
        Eigen::Vector2d const n2((*e->getNode(1))[0], (*e->getNode(1))[1]);
        Eigen::Vector2d const dist1(n1 - anchor);
        Eigen::Vector2d const dist2(n2 - anchor);
        std::size_t const id = e->getID();
        // elements located at left or right side
        if ((dist1 - dist2).squaredNorm() < eps)
        {
            top_bound_idx.push_back(id);
            bottom_bound_idx.push_back(id);
 
            if (dist1.squaredNorm() < eps)
            {
                right_bound_idx.push_back(id);
            }
            else
            {
                left_bound_idx.push_back(id);
            }
            continue;
        }
        // elements located at top or bottom
        if (dist2.squaredNorm() < dist1.squaredNorm())
        {
            top_bound_idx.push_back(id);
            left_bound_idx.push_back(id);
            right_bound_idx.push_back(id);
        }
        else
        {
            bottom_bound_idx.push_back(id);
            left_bound_idx.push_back(id);
            right_bound_idx.push_back(id);
        }
    }
 
    writeBoundary(*boundary_mesh, left_bound_idx, output_name + "_left");
    writeBoundary(*boundary_mesh, right_bound_idx, output_name + "_right");
    writeBoundary(*boundary_mesh, top_bound_idx, output_name + "_top");
    writeBoundary(*boundary_mesh, bottom_bound_idx, output_name + "_bottom");
}

References MeshLib::Cell, MeshLib::Face, MeshToolsLib::BoundaryExtraction::getBoundaryElementsAsMesh(), MeshLib::getBulkIDString(), MeshLib::Mesh::getElements(), MeshLib::Node, and writeBoundary().

Referenced by main().

generateMesh()

MeshLib::Mesh * generateMesh	(	GeoLib::GEOObjects &	geo,
		std::string const &	geo_name,
		std::string const &	output_name,
		double	res )

converts geometry into GMSH format and creates mesh

Definition at line 231 of file VerticalSliceFromLayers.cpp.

{
    std::string const gmsh_geo_name(output_name + ".geo");
    std::vector<std::string> gmsh_geo;
    gmsh_geo.push_back(geo_name);
    FileIO::GMSH::GMSHInterface gmsh_io(
        geo, true, FileIO::GMSH::MeshDensityAlgorithm::FixedMeshDensity, res, 0,
        0, gmsh_geo, false, false);
    gmsh_io.writePhysicalGroups(true);
    if (!BaseLib::IO::writeStringToFile(gmsh_io.writeToString(), gmsh_geo_name))
    {
        ERR("Writing gmsh geo file '{:s}' failed.", gmsh_geo_name);
    }
 
    std::string const gmsh_mesh_name = output_name + ".msh";
    std::string gmsh_command = "gmsh -2 -algo meshadapt " + gmsh_geo_name;
    gmsh_command += " -o " + gmsh_mesh_name + " -format msh22";
    int const return_value = std::system(gmsh_command.c_str());
    if (return_value != 0)
    {
        ERR("Execution of gmsh command returned non-zero status, %d",
            return_value);
    }
    return FileIO::GMSH::readGMSHMesh(gmsh_mesh_name,
                                      false /*is_created_with_gmsh2*/);
}

References ERR(), FileIO::GMSH::FixedMeshDensity, FileIO::GMSH::readGMSHMesh(), FileIO::GMSH::GMSHInterface::writePhysicalGroups(), BaseLib::IO::writeStringToFile(), and BaseLib::IO::Writer::writeToString().

Referenced by main().

getMinElevationPoint()

static GeoLib::Point * getMinElevationPoint ( GeoLib::Point * a, GeoLib::Point * b )

static

Definition at line 123 of file VerticalSliceFromLayers.cpp.

{
    return ((*a)[2] < (*b)[2]) ? a : b;
}

Referenced by mergeGeometries().

main()

int main	(	int	argc,
		char *	argv[] )

Definition at line 366 of file VerticalSliceFromLayers.cpp.

{
    QCoreApplication a(argc, argv);
    TCLAP::CmdLine cmd(
        "Creates a triangle-mesh of a vertical slice out of a list of input "
        "layer meshes. The slice is defined by a start- and end-point. In "
        "addition, the resolution for meshing the extracted slice (i.e. the "
        "maximum edge length of the domain discretisation) needs to be "
        "specified. The utility requires access to the meshing utility GMSH to "
        "work correctly.\n\n"
        "OpenGeoSys-6 software, version " +
            GitInfoLib::GitInfo::ogs_version +
            ".\n"
            "Copyright (c) 2012-2024, OpenGeoSys Community "
            "(http://www.opengeosys.org)",
        ' ', GitInfoLib::GitInfo::ogs_version);
    TCLAP::SwitchArg test_arg("t", "testdata", "keep test data", false);
    cmd.add(test_arg);
    TCLAP::SwitchArg bound_arg("b", "bounds", "save mesh boundaries", false);
    cmd.add(bound_arg);
    TCLAP::ValueArg<double> res_arg(
        "r", "resolution",
        "desired edge length of triangles in the resulting slice.", true, 0,
        "floating point number");
    cmd.add(res_arg);
    TCLAP::ValueArg<double> end_y_arg(
        "", "end-y", "y-coordinates of the end point defining the slice", true,
        0, "floating point number");
    cmd.add(end_y_arg);
    TCLAP::ValueArg<double> end_x_arg(
        "", "end-x", "x-coordinates of the end point defining the slice", true,
        0, "floating point number");
    cmd.add(end_x_arg);
    TCLAP::ValueArg<double> start_y_arg(
        "", "start-y", "y-coordinates of the start point defining the slice",
        true, 0, "floating point number");
    cmd.add(start_y_arg);
    TCLAP::ValueArg<double> start_x_arg(
        "", "start-x", "x-coordinates of the start point defining the slice",
        true, 0, "floating point number");
    cmd.add(start_x_arg);
    TCLAP::ValueArg<std::string> output_arg(
        "o", "output", "name of output mesh (*.vtu)", true, "", "string");
    cmd.add(output_arg);
    TCLAP::ValueArg<std::string> input_arg(
        "i", "input",
        "name of the input file list containing the paths the all input layers "
        "in correct order from top to bottom",
        true, "", "string");
    cmd.add(input_arg);
    cmd.parse(argc, argv);
 
    BaseLib::MPI::Setup mpi_setup(argc, argv);
 
    std::string const input_name = input_arg.getValue();
    std::string const output_name = output_arg.getValue();
 
    MathLib::Point3d const pnt_start{
        {start_x_arg.getValue(), start_y_arg.getValue(), 0.0}};
    MathLib::Point3d const pnt_end{
        {end_x_arg.getValue(), end_y_arg.getValue(), 0.0}};
    double const length = std::sqrt(MathLib::sqrDist(pnt_start, pnt_end));
 
    std::size_t const res = std::ceil(length / res_arg.getValue());
    double const interval_length = length / res;
 
    std::vector<std::string> const layer_names =
        BaseLib::IO::readStringListFromFile(input_name);
    if (layer_names.size() < 2)
    {
        ERR("At least two layers are required to extract a slice.");
        return EXIT_FAILURE;
    }
 
    GeoLib::GEOObjects geo;
    std::vector<std::string> const geo_name_list =
        createGeometries(geo, layer_names, pnt_start, pnt_end, res);
 
    if (geo_name_list.size() < 2)
    {
        ERR("Less than two geometries could be created from layers. Aborting "
            "extraction...");
        return EXIT_FAILURE;
    }
 
    std::string merged_geo_name = "merged_geometries";
    mergeGeometries(geo, geo_name_list, merged_geo_name);
    std::vector<GeoLib::Point*> points = *geo.getPointVec(merged_geo_name);
    auto const [rot_mat, z_shift] = rotateGeometryToXY(points);
    consolidateGeometry(geo, output_name, merged_geo_name, test_arg.getValue());
 
    std::unique_ptr<MeshLib::Mesh> mesh(
        generateMesh(geo, merged_geo_name, output_name, interval_length));
    if (!test_arg.getValue())
    {
        BaseLib::removeFile(output_name + ".geo");
        BaseLib::removeFile(output_name + ".msh");
    }
    if (mesh == nullptr)
    {
        ERR("Error generating mesh... (GMSH was unable to output mesh)");
        return EXIT_FAILURE;
    }
    rotateMesh(*mesh, rot_mat, z_shift);
    std::unique_ptr<MeshLib::Mesh> new_mesh(removeLineElements(*mesh));
    if (new_mesh == nullptr)
    {
        ERR("Error generating mesh... (GMSH created line mesh)");
        return EXIT_FAILURE;
    }
 
    // collapse all nodes that might have been created due to gmsh physically
    // separating layers
    MeshToolsLib::MeshRevision rev(*new_mesh);
    auto const edge_length = minMaxEdgeLength(new_mesh->getElements());
    double const minimum_node_distance = edge_length.first / 100.0;
 
    std::unique_ptr<MeshLib::Mesh> revised_mesh(
        rev.simplifyMesh("RevisedMesh", minimum_node_distance));
 
    if (bound_arg.getValue())
    {
        extractBoundaries(*revised_mesh, output_name, pnt_start);
    }
    MeshLib::IO::VtuInterface vtu(revised_mesh.get());
    vtu.writeToFile(output_name + ".vtu");
    return EXIT_SUCCESS;
}

References consolidateGeometry(), createGeometries(), ERR(), extractBoundaries(), generateMesh(), GeoLib::GEOObjects::getPointVec(), mergeGeometries(), GitInfoLib::GitInfo::ogs_version, BaseLib::IO::readStringListFromFile(), BaseLib::removeFile(), removeLineElements(), rotateGeometryToXY(), rotateMesh(), MeshToolsLib::MeshRevision::simplifyMesh(), MathLib::sqrDist(), and MeshLib::IO::VtuInterface::writeToFile().

mergeGeometries()

void mergeGeometries	(	GeoLib::GEOObjects &	geo,
		std::vector< std::string > const &	geo_names,
		std::string &	merged_geo_name )

Merges all layer geometries into one. Each layer is specified by one polygon. This will ensure that GMSH deals with the subsequent assignment of material groups automatically.

Definition at line 131 of file VerticalSliceFromLayers.cpp.

{
    std::vector<GeoLib::Point*> points;
    std::vector<GeoLib::Polyline*> lines;
 
    auto DEM_pnts = *geo.getPointVec(geo_names[0]);
    std::size_t const pnts_per_line = DEM_pnts.size();
    std::size_t const n_layers = geo_names.size();
    std::vector<std::size_t> last_line_idx(pnts_per_line, 0);
 
    auto layer_pnts = *geo.getPointVec(geo_names.back());
    for (std::size_t i = 0; i < pnts_per_line; ++i)
    {
        points.push_back(new GeoLib::Point(
            *getMinElevationPoint(layer_pnts[i], DEM_pnts[i]), i));
        last_line_idx[i] = i;
    }
    for (int j = n_layers - 2; j >= 0; --j)
    {
        GeoLib::Polyline* line = new GeoLib::Polyline(points);
        for (std::size_t i = 0; i < pnts_per_line; ++i)
        {
            line->addPoint(last_line_idx[pnts_per_line - i - 1]);
        }
        layer_pnts = *geo.getPointVec(geo_names[j]);
        for (std::size_t i = 0; i < pnts_per_line; ++i)
        {
            // check if for current point the upper layer boundary is actually
            // located above the upper boundary
            std::size_t idx = last_line_idx[i];
            if ((*points[idx])[2] < (*layer_pnts[i])[2])
            {
                idx = points.size();
                // check current point against DEM
                points.push_back(new GeoLib::Point(
                    *getMinElevationPoint(layer_pnts[i], DEM_pnts[i]), idx));
                last_line_idx[i] = idx;
            }
            line->addPoint(idx);
        }
        // close polygon
        line->addPoint(line->getPointID(0));
        lines.push_back(line);
    }
 
    geo.addPointVec(std::move(points), merged_geo_name,
                    GeoLib::PointVec::NameIdMap{});
    geo.addPolylineVec(std::move(lines), merged_geo_name,
                       GeoLib::PolylineVec::NameIdMap{});
}

References GeoLib::Polyline::addPoint(), GeoLib::GEOObjects::addPointVec(), GeoLib::GEOObjects::addPolylineVec(), getMinElevationPoint(), GeoLib::Polyline::getPointID(), and GeoLib::GEOObjects::getPointVec().

Referenced by main().

removeLineElements()

MeshLib::Mesh * removeLineElements

(

MeshLib::Mesh const &

mesh

)

removes line elements from mesh such that only triangles remain

Definition at line 271 of file VerticalSliceFromLayers.cpp.

{
    std::vector<std::size_t> line_idx;
    std::vector<MeshLib::Element*> const& elems = mesh.getElements();
    std::for_each(elems.begin(), elems.end(),
                  [&](auto e)
                  {
                      if (e->getGeomType() == MeshLib::MeshElemType::LINE)
                      {
                          line_idx.push_back(e->getID());
                      }
                  });
    if (line_idx.size() == mesh.getNumberOfElements())
    {
        return nullptr;
    }
    return MeshToolsLib::removeElements(mesh, line_idx, "mesh");
}

References MeshLib::Mesh::getElements().

Referenced by main().

rotateGeometryToXY()

std::pair< Eigen::Matrix3d, double > rotateGeometryToXY

(

std::vector< GeoLib::Point * > &

points

)

rotates the merged geometry into the XY-plane

Definition at line 185 of file VerticalSliceFromLayers.cpp.

{
    // compute the plane normal
    auto const [plane_normal, d] =
        GeoLib::getNewellPlane(points.begin(), points.end());
    // rotate points into x-y-plane
    Eigen::Matrix3d const rotation_matrix =
        GeoLib::computeRotationMatrixToXY(plane_normal);
    GeoLib::rotatePoints(rotation_matrix, points.begin(), points.end());
 
    GeoLib::AABB aabb(points.begin(), points.end());
    double const z_shift =
        (aabb.getMinPoint()[2] + aabb.getMaxPoint()[2]) / 2.0;
    std::for_each(points.begin(), points.end(),
                  [z_shift](GeoLib::Point* p) { (*p)[2] -= z_shift; });
    return {rotation_matrix, z_shift};
}

References GeoLib::computeRotationMatrixToXY(), GeoLib::AABB::getMaxPoint(), GeoLib::AABB::getMinPoint(), GeoLib::getNewellPlane(), and GeoLib::rotatePoints().

Referenced by main().

rotateMesh()

void rotateMesh	(	MeshLib::Mesh &	mesh,
		Eigen::Matrix3d const &	rot_mat,
		double const	z_shift )

inverse rotation of the mesh, back into original position

Definition at line 261 of file VerticalSliceFromLayers.cpp.

{
    std::vector<MeshLib::Node*> const& nodes = mesh.getNodes();
    std::for_each(nodes.begin(), nodes.end(),
                  [z_shift](MeshLib::Node* n) { (*n)[2] += z_shift; });
    GeoLib::rotatePoints(rot_mat.transpose(), nodes.begin(), nodes.end());
}

References MeshLib::Mesh::getNodes(), and GeoLib::rotatePoints().

Referenced by main().

writeBoundary()

void writeBoundary	(	MeshLib::Mesh const &	mesh,
		std::vector< std::size_t > const &	idx_array,
		std::string const &	file_name )

Definition at line 290 of file VerticalSliceFromLayers.cpp.

{
    std::unique_ptr<MeshLib::Mesh> const boundary(
        MeshToolsLib::removeElements(mesh, idx_array, "mesh"));
    if (boundary == nullptr)
    {
        ERR("Error extracting boundary '{:s}'", file_name);
        return;
    }
    MeshLib::IO::VtuInterface vtu(boundary.get());
    vtu.writeToFile(file_name + ".vtu");
}

References ERR(), MeshToolsLib::removeElements(), and MeshLib::IO::VtuInterface::writeToFile().

Referenced by extractBoundaries().