OGS
NodeReordering.cpp File Reference

Detailed Description

2013/13/06 KR Initial implementation

Definition in file NodeReordering.cpp.

#include <tclap/CmdLine.h>
#include <algorithm>
#include <array>
#include <memory>
#include <vector>
#include "BaseLib/Algorithm.h"
#include "BaseLib/Logging.h"
#include "BaseLib/MPI.h"
#include "BaseLib/TCLAPArguments.h"
#include "InfoLib/GitInfo.h"
#include "MeshLib/Elements/Element.h"
#include "MeshLib/IO/readMeshFromFile.h"
#include "MeshLib/IO/writeMeshToFile.h"
#include "MeshLib/Mesh.h"
#include "MeshLib/Node.h"
Include dependency graph for NodeReordering.cpp:

Go to the source code of this file.

Functions

void reverseNodeOrder (std::vector< MeshLib::Element * > &elements, bool const forced)
 Reverses order of nodes. In particular, this fixes issues between OGS5 and OGS6 meshes.
 
void fixVtkInconsistencies (std::vector< MeshLib::Element * > &elements)
 
void reorderNonlinearNodes (MeshLib::Mesh &mesh)
 Orders the base nodes of each elements before its non-linear nodes.
 
int main (int argc, char *argv[])
 

Function Documentation

◆ fixVtkInconsistencies()

void fixVtkInconsistencies ( std::vector< MeshLib::Element * > & elements)

Fixes inconsistencies between VTK's and OGS' node order for prism elements. In particular, this fixes issues between OGS6 meshes with and without InSitu-Lib

Definition at line 98 of file NodeReordering.cpp.

99{
100 std::size_t const nElements(elements.size());
101 for (std::size_t i = 0; i < nElements; ++i)
102 {
103 const unsigned nElemNodes(elements[i]->getNumberOfBaseNodes());
104 std::vector<MeshLib::Node*> nodes(elements[i]->getNodes(),
105 elements[i]->getNodes() + nElemNodes);
106
107 for (std::size_t j = 0; j < nElemNodes; ++j)
108 {
109 if (elements[i]->getGeomType() == MeshLib::MeshElemType::PRISM)
110 {
111 for (std::size_t k = 0; k < 3; ++k)
112 {
113 elements[i]->setNode(k, nodes[k + 3]);
114 elements[i]->setNode(k + 3, nodes[k]);
115 }
116 break;
117 }
118 }
119 }
120}
std::vector< std::size_t > getNodes(GeoLib::Point const &pnt, std::vector< MeshLib::Node * > const &nodes, MeshLib::PropertyVector< int > const &mat_ids, std::pair< int, int > const &mat_limits, std::pair< double, double > const &elevation_limits, MeshLib::Mesh const &mesh)

References getNodes(), and MeshLib::PRISM.

Referenced by main().

◆ main()

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

Definition at line 157 of file NodeReordering.cpp.

158{
159 TCLAP::CmdLine cmd(
160 "Reorders mesh nodes in elements to make old or incorrectly ordered "
161 "meshes compatible with OGS6.\n"
162 "Three options are available:\n"
163 "Method 0: Reversing order of nodes for all elements.\n"
164 "Method 1: Reversing order of nodes unless it's perceived correct by "
165 "OGS6 standards. This is the default selection.\n"
166 "Method 2: Fixing node ordering issues between VTK and OGS6 (only "
167 "applies to prism-elements)\n"
168 "Method 3: Re-ordering of mesh node vector such that all base nodes "
169 "are sorted before all nonlinear nodes.\n\n"
170 "OpenGeoSys-6 software, version " +
172 ".\n"
173 "Copyright (c) 2012-2025, OpenGeoSys Community "
174 "(http://www.opengeosys.org)",
176
177 std::vector<int> method_ids{0, 1, 2, 3};
178 TCLAP::ValuesConstraint<int> allowed_values(method_ids);
179 TCLAP::ValueArg<int> method_arg("m", "method",
180 "reordering method selection", false, 1,
181 &allowed_values);
182 cmd.add(method_arg);
183 TCLAP::ValueArg<std::string> output_mesh_arg(
184 "o", "output_mesh", "Output (.vtu). The name of the output mesh file",
185 true, "", "OUTPUT_FILE");
186 cmd.add(output_mesh_arg);
187 TCLAP::ValueArg<std::string> input_mesh_arg(
188 "i", "input_mesh",
189 "Input (.vtu | .vtk | .msh). The name of the input mesh file", true, "",
190 "INPUT_FILE");
191 cmd.add(input_mesh_arg);
192 auto log_level_arg = BaseLib::makeLogLevelArg();
193 cmd.add(log_level_arg);
194 cmd.parse(argc, argv);
195
196 BaseLib::MPI::Setup mpi_setup(argc, argv);
197 BaseLib::initOGSLogger(log_level_arg.getValue());
198
199 std::unique_ptr<MeshLib::Mesh> mesh(
200 MeshLib::IO::readMeshFromFile(input_mesh_arg.getValue()));
201
202 if (!mesh)
203 {
204 return EXIT_FAILURE;
205 }
206
207 INFO("Reordering nodes... ");
208 if (!method_arg.isSet() || method_arg.getValue() < 2)
209 {
210 bool const forced = (method_arg.getValue() == 0);
212 const_cast<std::vector<MeshLib::Element*>&>(mesh->getElements()),
213 forced);
214 }
215 else if (method_arg.getValue() == 2)
216 {
218 const_cast<std::vector<MeshLib::Element*>&>(mesh->getElements()));
219 }
220 else if (method_arg.getValue() == 3)
221 {
223 }
224
225 MeshLib::IO::writeMeshToFile(*mesh, output_mesh_arg.getValue());
226
227 INFO("VTU file written.");
228
229 return EXIT_SUCCESS;
230}
void INFO(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:36
void fixVtkInconsistencies(std::vector< MeshLib::Element * > &elements)
void reverseNodeOrder(std::vector< MeshLib::Element * > &elements, bool const forced)
Reverses order of nodes. In particular, this fixes issues between OGS5 and OGS6 meshes.
void reorderNonlinearNodes(MeshLib::Mesh &mesh)
Orders the base nodes of each elements before its non-linear nodes.
TCLAP::ValueArg< std::string > makeLogLevelArg()
void initOGSLogger(std::string const &log_level)
Definition Logging.cpp:64
GITINFOLIB_EXPORT const std::string ogs_version
MeshLib::Mesh * readMeshFromFile(const std::string &file_name, bool const compute_element_neighbors)
int writeMeshToFile(const MeshLib::Mesh &mesh, std::filesystem::path const &file_path, std::set< std::string > variable_output_names)

References fixVtkInconsistencies(), INFO(), BaseLib::initOGSLogger(), BaseLib::makeLogLevelArg(), GitInfoLib::GitInfo::ogs_version, MeshLib::IO::readMeshFromFile(), reorderNonlinearNodes(), reverseNodeOrder(), and MeshLib::IO::writeMeshToFile().

◆ reorderNonlinearNodes()

void reorderNonlinearNodes ( MeshLib::Mesh & mesh)

Orders the base nodes of each elements before its non-linear nodes.

Definition at line 123 of file NodeReordering.cpp.

124{
125 std::vector<MeshLib::Node*> base_nodes;
126 std::vector<MeshLib::Node*> nonlinear_nodes;
127 for (MeshLib::Element const* e : mesh.getElements())
128 {
129 for (unsigned i = 0; i < e->getNumberOfBaseNodes(); i++)
130 {
131 base_nodes.push_back(const_cast<MeshLib::Node*>(e->getNode(i)));
132 }
133 for (unsigned i = e->getNumberOfBaseNodes(); i < e->getNumberOfNodes();
134 i++)
135 {
136 nonlinear_nodes.push_back(
137 const_cast<MeshLib::Node*>(e->getNode(i)));
138 }
139 }
140
141 BaseLib::makeVectorUnique(base_nodes,
143 BaseLib::makeVectorUnique(nonlinear_nodes,
145
146 std::vector<MeshLib::Node*>& allnodes =
147 const_cast<std::vector<MeshLib::Node*>&>(mesh.getNodes());
148 allnodes.clear();
149
150 allnodes.insert(allnodes.end(), base_nodes.begin(), base_nodes.end());
151 allnodes.insert(allnodes.end(), nonlinear_nodes.begin(),
152 nonlinear_nodes.end());
153
154 mesh.resetNodeIDs();
155}
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition Mesh.h:108
void resetNodeIDs()
Resets the IDs of all mesh-nodes to their position in the node vector.
Definition Mesh.cpp:160
void makeVectorUnique(std::vector< T > &v)
Definition Algorithm.h:180
bool idsComparator(T const a, T const b)
Definition Mesh.h:208

References MeshLib::Mesh::getElements(), MeshLib::Mesh::getNodes(), MeshLib::idsComparator(), BaseLib::makeVectorUnique(), and MeshLib::Mesh::resetNodeIDs().

Referenced by main().

◆ reverseNodeOrder()

void reverseNodeOrder ( std::vector< MeshLib::Element * > & elements,
bool const forced )

Reverses order of nodes. In particular, this fixes issues between OGS5 and OGS6 meshes.

Parameters
elementsMesh elements whose nodes should be reordered
forcedIf true, nodes are reordered for all elements, if false it is first checked if the node order is correct according to OGS6 element definitions.

Definition at line 39 of file NodeReordering.cpp.

41{
42 std::size_t n_corrected_elements = 0;
43 std::size_t const nElements(elements.size());
44 for (std::size_t i = 0; i < nElements; ++i)
45 {
46 if (!forced && elements[i]->testElementNodeOrder())
47 {
48 continue;
49 }
50 n_corrected_elements++;
51
52 const unsigned nElemNodes(elements[i]->getNumberOfBaseNodes());
53 std::vector<MeshLib::Node*> nodes(elements[i]->getNodes(),
54 elements[i]->getNodes() + nElemNodes);
55
56 switch (elements[i]->getGeomType())
57 {
59 for (std::size_t j = 0; j < 4; ++j)
60 {
61 elements[i]->setNode(j, nodes[(j + 1) % 4]);
62 }
63 break;
65 elements[i]->setNode(0, nodes[1]);
66 elements[i]->setNode(1, nodes[0]);
67 elements[i]->setNode(2, nodes[3]);
68 elements[i]->setNode(3, nodes[2]);
69 break;
71 for (std::size_t j = 0; j < 3; ++j)
72 {
73 elements[i]->setNode(j, nodes[j + 3]);
74 elements[i]->setNode(j + 3, nodes[j]);
75 }
76 break;
78 for (std::size_t j = 0; j < 4; ++j)
79 {
80 elements[i]->setNode(j, nodes[j + 4]);
81 elements[i]->setNode(j + 4, nodes[j]);
82 }
83 break;
84 default:
85 for (std::size_t j = 0; j < nElemNodes; ++j)
86 {
87 elements[i]->setNode(j, nodes[nElemNodes - j - 1]);
88 }
89 }
90 }
91
92 INFO("Corrected {:d} elements.", n_corrected_elements);
93}

References getNodes(), MeshLib::HEXAHEDRON, INFO(), MeshLib::PRISM, MeshLib::PYRAMID, and MeshLib::TETRAHEDRON.

Referenced by main().