OGS
NodeReordering.cpp
Go to the documentation of this file.
1// SPDX-FileCopyrightText: Copyright (c) OpenGeoSys Community (opengeosys.org)
2// SPDX-License-Identifier: BSD-3-Clause
3
4#include <tclap/CmdLine.h>
5
6#include <Eigen/Dense>
7#include <algorithm>
8#include <array>
9#include <concepts>
10#include <functional>
11#include <vector>
12
13#include "BaseLib/Algorithm.h"
14#include "BaseLib/Logging.h"
15#include "BaseLib/MPI.h"
16#include "BaseLib/TCLAPArguments.h"
17#include "InfoLib/GitInfo.h"
18#include "MeshLib/ElementCoordinatesMappingLocal.h"
19#include "MeshLib/Elements/Element.h"
20#include "MeshLib/IO/readMeshFromFile.h"
21#include "MeshLib/IO/writeMeshToFile.h"
22#include "MeshLib/Mesh.h"
23#include "MeshLib/MeshEnums.h"
24#include "MeshLib/Node.h"
25#include "MeshLib/Utils/GetSpaceDimension.h"
26#include "NumLib/Fem/ShapeFunction/ShapeHex20.h"
27#include "NumLib/Fem/ShapeFunction/ShapeHex8.h"
28#include "NumLib/Fem/ShapeFunction/ShapeLine2.h"
29#include "NumLib/Fem/ShapeFunction/ShapeLine3.h"
30#include "NumLib/Fem/ShapeFunction/ShapePoint1.h"
31#include "NumLib/Fem/ShapeFunction/ShapePrism15.h"
32#include "NumLib/Fem/ShapeFunction/ShapePrism6.h"
33#include "NumLib/Fem/ShapeFunction/ShapePyra13.h"
34#include "NumLib/Fem/ShapeFunction/ShapePyra5.h"
35#include "NumLib/Fem/ShapeFunction/ShapeQuad4.h"
36#include "NumLib/Fem/ShapeFunction/ShapeQuad8.h"
37#include "NumLib/Fem/ShapeFunction/ShapeQuad9.h"
38#include "NumLib/Fem/ShapeFunction/ShapeTet10.h"
39#include "NumLib/Fem/ShapeFunction/ShapeTet4.h"
40#include "NumLib/Fem/ShapeFunction/ShapeTri3.h"
41#include "NumLib/Fem/ShapeFunction/ShapeTri6.h"
42
43template <typename GradShapeFunction, int Dim>
44bool checkJacobianDeterminant(MeshLib::Element const& e,
45 int const mesh_space_dimension,
46 std::array<double, Dim> const& xi,
47 bool const check_reordered = false)
48{
49 Eigen::Matrix<double, GradShapeFunction::DIM, GradShapeFunction::NPOINTS,
50 Eigen::RowMajor>
51 dNdxi;
52 Eigen::Map<Eigen::VectorXd> dN_vec(dNdxi.data(), dNdxi.size());
53 GradShapeFunction::computeGradShapeFunction(xi.data(), dN_vec);
54 MeshLib::ElementCoordinatesMappingLocal ele_local_coord(
55 e, mesh_space_dimension);
56
57 Eigen::MatrixXd J = Eigen::MatrixXd::Zero(Dim, Dim);
58 assert(e.getNumberOfNodes() == GradShapeFunction::NPOINTS);
59 for (unsigned k = 0; k < GradShapeFunction::NPOINTS; k++)
60 {
61 const MathLib::Point3d& mapped_pt =
62 ele_local_coord.getMappedCoordinates(k);
63 J += dNdxi.col(k) * mapped_pt.asEigenVector3d().head(Dim).transpose();
64 }
65
66 if (!check_reordered)
67 {
68 return !(J.determinant() < 0);
69 }
70
71 if (J.determinant() < 0)
72 {
73 OGS_FATAL(
74 "Element {:d} has negative Jacobian determinant {:g}. "
75 "NodeReordering fails.",
76 e.getID(), J.determinant());
77 }
78
79 return true;
80}
81
82template <typename T>
83concept ShapeFunction = requires(double* xi, Eigen::Map<Eigen::VectorXd> dN) {
84 { T::DIM } -> std::convertible_to<int>;
85 { T::NPOINTS } -> std::convertible_to<int>;
86 T::computeGradShapeFunction(xi, dN);
87};
88
89struct ElementReorderConfigBase
90{
91 std::function<bool(MeshLib::Element&, int, bool)> is_node_ordering_correct;
92 std::function<void(MeshLib::Element&, const std::vector<MeshLib::Node*>&)>
93 reorder_element_nodes;
94};
95
96template <ShapeFunction ShapeFunc, int Dim>
97ElementReorderConfigBase makeElementConfig(
98 std::array<double, Dim> xi,
99 std::function<void(MeshLib::Element&, const std::vector<MeshLib::Node*>&)>
100 reorder)
101{
102 return ElementReorderConfigBase{
103 [xi](MeshLib::Element& e, int mesh_space_dimension,
104 bool check_reordered)
105 {
106 return checkJacobianDeterminant<ShapeFunc, Dim>(
107 e, mesh_space_dimension, xi, check_reordered);
108 },
109 reorder};
110}
111
112static const std::array<ElementReorderConfigBase,
113 static_cast<int>(MeshLib::CellType::enum_length)>
114 element_configs_array = []
115{
116 auto swap_nodes_i_j = [](MeshLib::Element& element,
117 const std::vector<MeshLib::Node*>& nodes,
118 unsigned const i, unsigned const j)
119 {
120 element.setNode(i, nodes[j]);
121 element.setNode(j, nodes[i]);
122 };
123
124 auto order_nodes_quadratic_quad =
125 [swap_nodes_i_j](MeshLib::Element& element,
126 const std::vector<MeshLib::Node*>& nodes)
127 {
128 swap_nodes_i_j(element, nodes, 1, 3);
129 swap_nodes_i_j(element, nodes, 5, 6);
130 swap_nodes_i_j(element, nodes, 4, 7);
131 };
132
133 std::array<ElementReorderConfigBase,
134 static_cast<int>(MeshLib::CellType::enum_length)>
135 arr{};
136
137 arr[static_cast<int>(MeshLib::CellType::LINE2)] =
138 makeElementConfig<NumLib::ShapeLine2, 1>(
139 {0.5},
140 [&swap_nodes_i_j](MeshLib::Element& element,
141 const std::vector<MeshLib::Node*>& nodes)
142 { swap_nodes_i_j(element, nodes, 0, 1); });
143
144 arr[static_cast<int>(MeshLib::CellType::LINE3)] =
145 makeElementConfig<NumLib::ShapeLine3, 1>(
146 {0.5},
147 [&swap_nodes_i_j](MeshLib::Element& element,
148 const std::vector<MeshLib::Node*>& nodes)
149 { swap_nodes_i_j(element, nodes, 0, 1); });
150
151 arr[static_cast<int>(MeshLib::CellType::TRI3)] =
152 makeElementConfig<NumLib::ShapeTri3, 2>(
153 {1.0 / 3.0, 1.0 / 3.0},
154 [&swap_nodes_i_j](MeshLib::Element& element,
155 const std::vector<MeshLib::Node*>& nodes)
156 { swap_nodes_i_j(element, nodes, 1, 2); });
157
158 arr[static_cast<int>(MeshLib::CellType::TRI6)] =
159 makeElementConfig<NumLib::ShapeTri6, 2>(
160 {1.0 / 3.0, 1.0 / 3.0},
161 [&swap_nodes_i_j](MeshLib::Element& element,
162 const std::vector<MeshLib::Node*>& nodes)
163 {
164 swap_nodes_i_j(element, nodes, 1, 2);
165 swap_nodes_i_j(element, nodes, 3, 5);
166 });
167
168 arr[static_cast<int>(MeshLib::CellType::QUAD4)] =
169 makeElementConfig<NumLib::ShapeQuad4, 2>(
170 {0.0, 0.0},
171 [&swap_nodes_i_j](MeshLib::Element& element,
172 const std::vector<MeshLib::Node*>& nodes)
173 { swap_nodes_i_j(element, nodes, 0, 2); });
174
175 arr[static_cast<int>(MeshLib::CellType::QUAD8)] =
176 makeElementConfig<NumLib::ShapeQuad8, 2>(
177 {0.0, 0.0},
178 [&order_nodes_quadratic_quad](
179 MeshLib::Element& element,
180 const std::vector<MeshLib::Node*>& nodes)
181 { order_nodes_quadratic_quad(element, nodes); });
182 arr[static_cast<int>(MeshLib::CellType::QUAD9)] =
183 makeElementConfig<NumLib::ShapeQuad9, 2>(
184 {0.0, 0.0},
185 [&order_nodes_quadratic_quad](
186 MeshLib::Element& element,
187 const std::vector<MeshLib::Node*>& nodes)
188 { order_nodes_quadratic_quad(element, nodes); });
189 arr[static_cast<int>(MeshLib::CellType::TET4)] =
190 makeElementConfig<NumLib::ShapeTet4, 3>(
191 {0.25, 0.25, 0.25},
192 [&swap_nodes_i_j](MeshLib::Element& element,
193 const std::vector<MeshLib::Node*>& nodes)
194 { swap_nodes_i_j(element, nodes, 1, 2); });
195
196 arr[static_cast<int>(MeshLib::CellType::TET10)] =
197 makeElementConfig<NumLib::ShapeTet10, 3>(
198 {0.25, 0.25, 0.25},
199 [&swap_nodes_i_j](MeshLib::Element& element,
200 const std::vector<MeshLib::Node*>& nodes)
201 {
202 swap_nodes_i_j(element, nodes, 1, 2);
203 swap_nodes_i_j(element, nodes, 4, 6);
204 swap_nodes_i_j(element, nodes, 8, 9);
205 });
206 arr[static_cast<int>(MeshLib::CellType::PRISM6)] =
207 makeElementConfig<NumLib::ShapePrism6, 3>(
208 {1.0 / 3.0, 1.0 / 3.0, 0.5},
209 [&swap_nodes_i_j](MeshLib::Element& element,
210 const std::vector<MeshLib::Node*>& nodes)
211 {
212 swap_nodes_i_j(element, nodes, 1, 2);
213 swap_nodes_i_j(element, nodes, 4, 5);
214 });
215
216 arr[static_cast<int>(MeshLib::CellType::PRISM15)] =
217 makeElementConfig<NumLib::ShapePrism15, 3>(
218 {1.0 / 3.0, 1.0 / 3., 0.5},
219 [&swap_nodes_i_j](MeshLib::Element& element,
220 const std::vector<MeshLib::Node*>& nodes)
221 {
222 swap_nodes_i_j(element, nodes, 0, 3);
223 swap_nodes_i_j(element, nodes, 1, 4);
224 swap_nodes_i_j(element, nodes, 2, 5);
225 swap_nodes_i_j(element, nodes, 6, 9);
226 swap_nodes_i_j(element, nodes, 7, 10);
227 swap_nodes_i_j(element, nodes, 8, 11);
228 });
229
230 arr[static_cast<int>(MeshLib::CellType::PYRAMID5)] =
231 makeElementConfig<NumLib::ShapePyra5, 3>(
232 {0.25, 0.25, 0.5},
233 [&swap_nodes_i_j](MeshLib::Element& element,
234 const std::vector<MeshLib::Node*>& nodes)
235 { swap_nodes_i_j(element, nodes, 0, 2); });
236
237 arr[static_cast<int>(MeshLib::CellType::PYRAMID13)] =
238 makeElementConfig<NumLib::ShapePyra13, 3>(
239 {0.25, 0.25, 0.5},
240 [&swap_nodes_i_j](MeshLib::Element& element,
241 const std::vector<MeshLib::Node*>& nodes)
242 {
243 swap_nodes_i_j(element, nodes, 0, 2);
244 swap_nodes_i_j(element, nodes, 9, 11);
245 swap_nodes_i_j(element, nodes, 5, 6);
246 swap_nodes_i_j(element, nodes, 7, 8);
247 });
248 arr[static_cast<int>(MeshLib::CellType::HEX8)] =
249 makeElementConfig<NumLib::ShapeHex8, 3>(
250 {0.5, 0.5, 0.5},
251 [&swap_nodes_i_j](MeshLib::Element& element,
252 const std::vector<MeshLib::Node*>& nodes)
253 {
254 swap_nodes_i_j(element, nodes, 0, 2);
255 swap_nodes_i_j(element, nodes, 4, 6);
256 });
257
258 arr[static_cast<int>(MeshLib::CellType::HEX20)] =
259 makeElementConfig<NumLib::ShapeHex20, 3>(
260 {0.5, 0.5, 0.5},
261 [&swap_nodes_i_j](MeshLib::Element& element,
262 const std::vector<MeshLib::Node*>& nodes)
263 {
264 swap_nodes_i_j(element, nodes, 0, 2);
265 swap_nodes_i_j(element, nodes, 4, 6);
266 swap_nodes_i_j(element, nodes, 16, 18);
267 swap_nodes_i_j(element, nodes, 8, 9);
268 swap_nodes_i_j(element, nodes, 10, 11);
269 swap_nodes_i_j(element, nodes, 12, 13);
270 swap_nodes_i_j(element, nodes, 14, 15);
271 });
272 return arr;
273}();
274
285void reverseNodeOrder(std::vector<MeshLib::Element*>& elements,
286 int const mesh_space_dimension, bool const forced,
287 bool const volume_check)
288{
289 std::size_t n_corrected_elements = 0;
290
291 for (auto* element : elements)
292 {
293 auto const cell_type = element->getCellType();
294
295 if (cell_type == MeshLib::CellType::INVALID)
296 {
297 OGS_FATAL("Element type `{:s}' does not exist.",
298 MeshLib::CellType2String(cell_type));
299 }
300 // This check is already done in MeshLib::readMeshFromFile(). Therefore,
301 // this is just a safeguard.
302 if (cell_type == MeshLib::CellType::HEX27 ||
303 cell_type == MeshLib::CellType::PRISM18)
304 {
305 OGS_FATAL("Element type `{:s}' is not supported in OGS.",
306 MeshLib::CellType2String(cell_type));
307 }
308
309 if (cell_type == MeshLib::CellType::POINT1)
310 {
311 continue;
312 }
313
314 const auto& element_config =
315 element_configs_array[static_cast<int>(cell_type)];
316
317 if (element_config.is_node_ordering_correct(
318 *element, mesh_space_dimension, false /*check_reordered*/))
319 {
320 continue;
321 }
322
323 // Save nodes before reordering
324 const unsigned nElemNodes = element->getNumberOfNodes();
325 std::vector<MeshLib::Node*> nodes(element->getNodes(),
326 element->getNodes() + nElemNodes);
327
328 if (volume_check)
329 {
330 double const element_volume_origin = element->computeVolume();
331
332 element_config.reorder_element_nodes(*element, nodes);
333
334 // Ensure that the element volume did not change.
335 double const element_volume = element->computeVolume();
336 // We use a fixed tolerance here because for very small elements
337 // the machine epsilon might be too small.
338 if (std::fabs(element_volume - element_volume_origin) /
339 element_volume_origin >
340 100 * std::numeric_limits<double>::epsilon())
341 {
342 OGS_FATAL(
343 "Reordering the nodes of element {:d} failed as its volume "
344 "changed from {:.20g} to {:.20g}, the relative difference "
345 "is {:.20g} and the threshold is {:.20g}.",
346 element->getID(), element_volume_origin, element_volume,
347 std::fabs(element_volume_origin - element_volume) /
348 element_volume_origin,
349 100 * std::numeric_limits<double>::epsilon());
350 }
351 }
352 else
353 {
354 element_config.reorder_element_nodes(*element, nodes);
355 }
356
357 // Element::computeVolume() uses the element vertices to compute
358 // the element volume, so the change of edge nodes are not
359 // considered. Therefore, we need to additionally check the Jacobian
360 // determinant here if forced is true.
361 if (forced)
362 {
363 element_config.is_node_ordering_correct(
364 *element, mesh_space_dimension, true /*check_reordered*/);
365 }
366
367 ++n_corrected_elements;
368 }
369
370 INFO("Corrected {:d} elements.", n_corrected_elements);
371}
372
377void fixVtkInconsistencies(std::vector<MeshLib::Element*>& elements)
378{
379 for (auto* const element : elements)
380 {
381 const unsigned nElemNodes(element->getNumberOfBaseNodes());
382 std::vector<MeshLib::Node*> nodes(element->getNodes(),
383 element->getNodes() + nElemNodes);
384
385 for (std::size_t j = 0; j < nElemNodes; ++j)
386 {
387 if (element->getGeomType() == MeshLib::MeshElemType::PRISM)
388 {
389 for (std::size_t k = 0; k < 3; ++k)
390 {
391 element->setNode(k, nodes[k + 3]);
392 element->setNode(k + 3, nodes[k]);
393 }
394 break;
395 }
396 }
397 }
398}
399
401void reorderNonlinearNodes(MeshLib::Mesh& mesh)
402{
403 std::vector<MeshLib::Node*> base_nodes;
404 std::vector<MeshLib::Node*> nonlinear_nodes;
405 for (MeshLib::Element const* e : mesh.getElements())
406 {
407 for (unsigned i = 0; i < e->getNumberOfBaseNodes(); i++)
408 {
409 base_nodes.push_back(const_cast<MeshLib::Node*>(e->getNode(i)));
410 }
411 for (unsigned i = e->getNumberOfBaseNodes(); i < e->getNumberOfNodes();
412 i++)
413 {
414 nonlinear_nodes.push_back(
415 const_cast<MeshLib::Node*>(e->getNode(i)));
416 }
417 }
418
419 BaseLib::makeVectorUnique(base_nodes,
420 MeshLib::idsComparator<MeshLib::Node*>);
421 BaseLib::makeVectorUnique(nonlinear_nodes,
422 MeshLib::idsComparator<MeshLib::Node*>);
423
424 std::vector<MeshLib::Node*>& allnodes =
425 const_cast<std::vector<MeshLib::Node*>&>(mesh.getNodes());
426 allnodes.clear();
427
428 allnodes.insert(allnodes.end(), base_nodes.begin(), base_nodes.end());
429 allnodes.insert(allnodes.end(), nonlinear_nodes.begin(),
430 nonlinear_nodes.end());
431
432 mesh.resetNodeIDs();
433}
434
435int main(int argc, char* argv[])
436{
437 enum class ExpectedCellType
438 {
439 INVALID = 0,
440 POINT1 = 1,
441 LINE2 = 2,
442 LINE3 = 3,
443 TRI3 = 4,
444 TRI6 = 5,
445 QUAD4 = 6,
446 QUAD8 = 7,
447 QUAD9 = 8,
448 TET4 = 9,
449 TET10 = 10,
450 HEX8 = 11,
451 HEX20 = 12,
452 HEX27 = 13,
453 PRISM6 = 14,
454 PRISM15 = 15,
455 PRISM18 = 16,
456 PYRAMID5 = 17,
457 PYRAMID13 = 18,
458 enum_length
459 };
460
461 constexpr bool are_expected_cell_types =
462 static_cast<int>(ExpectedCellType::enum_length) ==
463 static_cast<int>(MeshLib::CellType::enum_length) &&
464 static_cast<int>(ExpectedCellType::INVALID) ==
465 static_cast<int>(MeshLib::CellType::INVALID) &&
466 static_cast<int>(ExpectedCellType::POINT1) ==
467 static_cast<int>(MeshLib::CellType::POINT1) &&
468 static_cast<int>(ExpectedCellType::LINE2) ==
469 static_cast<int>(MeshLib::CellType::LINE2) &&
470 static_cast<int>(ExpectedCellType::LINE3) ==
471 static_cast<int>(MeshLib::CellType::LINE3) &&
472 static_cast<int>(ExpectedCellType::TRI3) ==
473 static_cast<int>(MeshLib::CellType::TRI3) &&
474 static_cast<int>(ExpectedCellType::TRI6) ==
475 static_cast<int>(MeshLib::CellType::TRI6) &&
476 static_cast<int>(ExpectedCellType::QUAD4) ==
477 static_cast<int>(MeshLib::CellType::QUAD4) &&
478 static_cast<int>(ExpectedCellType::QUAD8) ==
479 static_cast<int>(MeshLib::CellType::QUAD8) &&
480 static_cast<int>(ExpectedCellType::QUAD9) ==
481 static_cast<int>(MeshLib::CellType::QUAD9) &&
482 static_cast<int>(ExpectedCellType::TET4) ==
483 static_cast<int>(MeshLib::CellType::TET4) &&
484 static_cast<int>(ExpectedCellType::TET10) ==
485 static_cast<int>(MeshLib::CellType::TET10) &&
486 static_cast<int>(ExpectedCellType::HEX8) ==
487 static_cast<int>(MeshLib::CellType::HEX8) &&
488 static_cast<int>(ExpectedCellType::HEX20) ==
489 static_cast<int>(MeshLib::CellType::HEX20) &&
490 static_cast<int>(ExpectedCellType::HEX27) ==
491 static_cast<int>(MeshLib::CellType::HEX27) &&
492 static_cast<int>(ExpectedCellType::PRISM6) ==
493 static_cast<int>(MeshLib::CellType::PRISM6) &&
494 static_cast<int>(ExpectedCellType::PRISM15) ==
495 static_cast<int>(MeshLib::CellType::PRISM15) &&
496 static_cast<int>(ExpectedCellType::PRISM18) ==
497 static_cast<int>(MeshLib::CellType::PRISM18) &&
498 static_cast<int>(ExpectedCellType::PYRAMID5) ==
499 static_cast<int>(MeshLib::CellType::PYRAMID5) &&
500 static_cast<int>(ExpectedCellType::PYRAMID13) ==
501 static_cast<int>(MeshLib::CellType::PYRAMID13);
502 // This error should never occur unless someone has changed the
503 // MeshLib::CellType enum class. These assertions ensure that the
504 // array 'element_configs_array' is up to date.
505 if (!are_expected_cell_types)
506 {
507 OGS_FATAL(
508 "The enum class MeshLib::CellType has changed. Please adapt the "
509 "array 'element_configs_array' in NodeReordering.cpp accordingly.");
510 }
511
512 TCLAP::CmdLine cmd(
513 "Reorders mesh nodes in elements to make old or incorrectly ordered "
514 "meshes compatible with OGS6.\n"
515 "Three options are available:\n"
516 "Method 0: Reversing order of nodes and checking again for all "
517 "elements.\n"
518 "Method 1: Reversing order of nodes unless it's perceived correct by "
519 "OGS6 standards. This is the default selection.\n"
520 "Method 2: Fixing node ordering issues between VTK and OGS6 (only "
521 "applies to prism-elements)\n"
522 "Method 3: Re-ordering of mesh node vector such that all base nodes "
523 "are sorted before all nonlinear nodes.\n\n"
524 "OpenGeoSys-6 software, version " +
525 GitInfoLib::GitInfo::ogs_version +
526 ".\n"
527 "Copyright (c) 2012-2026, OpenGeoSys Community "
528 "(http://www.opengeosys.org)",
529 ' ', GitInfoLib::GitInfo::ogs_version);
530
531 TCLAP::SwitchArg no_volume_check_arg(
532 "", "no_volume_check",
533 "By default the volumes of original and reordered elements are "
534 "compared if they are numerically equal, i.e., relative volume "
535 "difference is smaller than a threshold. This switch disables the "
536 "volume comparison.");
537 cmd.add(no_volume_check_arg);
538
539 std::vector<int> method_ids{0, 1, 2, 3};
540 TCLAP::ValuesConstraint<int> allowed_values(method_ids);
541 TCLAP::ValueArg<int> method_arg("m", "method",
542 "reordering method selection", false, 1,
543 &allowed_values);
544 cmd.add(method_arg);
545 TCLAP::ValueArg<std::string> output_mesh_arg(
546 "o", "output_mesh", "Output (.vtu). The name of the output mesh file",
547 true, "", "OUTPUT_FILE");
548 cmd.add(output_mesh_arg);
549 TCLAP::ValueArg<std::string> input_mesh_arg(
550 "i", "input_mesh",
551 "Input (.vtu | .vtk | .msh). The name of the input mesh file", true, "",
552 "INPUT_FILE");
553 cmd.add(input_mesh_arg);
554 auto log_level_arg = BaseLib::makeLogLevelArg();
555 cmd.add(log_level_arg);
556 cmd.parse(argc, argv);
557
558 BaseLib::MPI::Setup mpi_setup(argc, argv);
559 BaseLib::initOGSLogger(log_level_arg.getValue());
560
561 std::unique_ptr<MeshLib::Mesh> mesh(
562 MeshLib::IO::readMeshFromFile(input_mesh_arg.getValue()));
563
564 if (!mesh)
565 {
566 return EXIT_FAILURE;
567 }
568
569 INFO("Reordering nodes... ");
570 if (!method_arg.isSet() || method_arg.getValue() < 2)
571 {
572 bool const forced = (method_arg.getValue() == 0);
573
574 if (forced)
575 {
576 INFO("Method 0: Reversing order of nodes will be checked again.");
577 }
578 INFO(
579 "Method: Reversing order of nodes unless it is considered correct "
580 "by the OGS6 standard, i.e. such that det(J) > 0, where J is the "
581 "Jacobian of the global-to-local coordinate transformation.");
582 int const mesh_space_dimension =
583 MeshLib::getSpaceDimension(mesh->getNodes());
584 reverseNodeOrder(
585 const_cast<std::vector<MeshLib::Element*>&>(mesh->getElements()),
586 mesh_space_dimension, forced, !no_volume_check_arg.isSet());
587 }
588 else if (method_arg.getValue() == 2)
589 {
590 fixVtkInconsistencies(
591 const_cast<std::vector<MeshLib::Element*>&>(mesh->getElements()));
592 }
593 else if (method_arg.getValue() == 3)
594 {
595 reorderNonlinearNodes(*mesh);
596 }
597
598 MeshLib::IO::writeMeshToFile(*mesh, output_mesh_arg.getValue());
599
600 INFO("VTU file written.");
601
602 return EXIT_SUCCESS;
603}
OGS_FATAL
#define OGS_FATAL(...)
Definition Error.h:19
INFO
void INFO(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:28
main
int main(int argc, char *argv[])
Definition NodeReordering.cpp:435
checkJacobianDeterminant
bool checkJacobianDeterminant(MeshLib::Element const &e, int const mesh_space_dimension, std::array< double, Dim > const &xi, bool const check_reordered=false)
Definition NodeReordering.cpp:44
element_configs_array
static const std::array< ElementReorderConfigBase, static_cast< int >(MeshLib::CellType::enum_length)> element_configs_array
Definition NodeReordering.cpp:114
fixVtkInconsistencies
void fixVtkInconsistencies(std::vector< MeshLib::Element * > &elements)
Definition NodeReordering.cpp:377
reorderNonlinearNodes
void reorderNonlinearNodes(MeshLib::Mesh &mesh)
Orders the base nodes of each elements before its non-linear nodes.
Definition NodeReordering.cpp:401
makeElementConfig
ElementReorderConfigBase makeElementConfig(std::array< double, Dim > xi, std::function< void(MeshLib::Element &, const std::vector< MeshLib::Node * > &)> reorder)
Definition NodeReordering.cpp:97
reverseNodeOrder
void reverseNodeOrder(std::vector< MeshLib::Element * > &elements, int const mesh_space_dimension, bool const forced, bool const volume_check)
Reverses order of nodes. In particular, this fixes issues between (Gmsh or OGS5) and OGS6 meshes.
Definition NodeReordering.cpp:285
MathLib::Point3d::asEigenVector3d
Eigen::Vector3d const & asEigenVector3d() const
Definition Point3d.h:55
MeshLib::ElementCoordinatesMappingLocal::getMappedCoordinates
MathLib::Point3d const & getMappedCoordinates(std::size_t node_id) const
return mapped coordinates of the node
Definition ElementCoordinatesMappingLocal.h:37
MeshLib::Element::getNumberOfNodes
virtual unsigned getNumberOfNodes() const =0
MeshLib::Element::setNode
virtual void setNode(unsigned idx, Node *node)=0
MeshLib::Element::getID
std::size_t getID() const
Returns the ID of the element.
Definition Element.h:80
MeshLib::Mesh::getNodes
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition Mesh.h:97
MeshLib::Mesh::getElements
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition Mesh.h:100
MeshLib::Mesh::resetNodeIDs
void resetNodeIDs()
Resets the IDs of all mesh-nodes to their position in the node vector.
Definition Mesh.cpp:149
BaseLib::makeLogLevelArg
TCLAP::ValueArg< std::string > makeLogLevelArg()
Definition TCLAPArguments.cpp:11
BaseLib::initOGSLogger
void initOGSLogger(std::string const &log_level)
Definition Logging.cpp:56
BaseLib::makeVectorUnique
void makeVectorUnique(std::vector< T > &v)
Definition Algorithm.h:173
GitInfoLib::GitInfo::ogs_version
GITINFOLIB_EXPORT const std::string ogs_version
MeshLib::IO::readMeshFromFile
MeshLib::Mesh * readMeshFromFile(const std::string &file_name, bool const compute_element_neighbors)
Definition readMeshFromFile.cpp:58
MeshLib::IO::writeMeshToFile
int writeMeshToFile(const MeshLib::Mesh &mesh, std::filesystem::path const &file_path, std::set< std::string > variable_output_names)
Definition writeMeshToFile.cpp:18
MeshLib::idsComparator
bool idsComparator(T const a, T const b)
Definition Mesh.h:197
MeshLib::CellType2String
std::string CellType2String(const CellType t)
Given a MeshElemType this returns the appropriate string.
Definition MeshEnums.cpp:146
MeshLib::getSpaceDimension
int getSpaceDimension(std::vector< Node * > const &nodes)
Computes dimension of the embedding space containing the set of given points.
Definition GetSpaceDimension.cpp:15
ElementReorderConfigBase::is_node_ordering_correct
std::function< bool(MeshLib::Element &, int, bool)> is_node_ordering_correct
Definition NodeReordering.cpp:91
ElementReorderConfigBase::reorder_element_nodes
std::function< void(MeshLib::Element &, const std::vector< MeshLib::Node * > &)> reorder_element_nodes
Definition NodeReordering.cpp:93