From 05a7f742b70c576b6ea31cf7a278e15b25510c6e Mon Sep 17 00:00:00 2001
From: Norihiro Watanabe <norihiro.watanabe@ufz.de>
Date: Tue, 25 Oct 2016 16:54:15 +0200
Subject: [PATCH] [App/Utils] add the new tool createQuadraticMesh
---
Applications/Utils/MeshEdit/CMakeLists.txt | 5 +
.../Utils/MeshEdit/createQuadraticeMesh.cpp | 163 ++++++++++++++++++
2 files changed, 168 insertions(+)
create mode 100644 Applications/Utils/MeshEdit/createQuadraticeMesh.cpp
diff --git a/Applications/Utils/MeshEdit/CMakeLists.txt b/Applications/Utils/MeshEdit/CMakeLists.txt
index a44be8c211e..35ef24556d8 100644
--- a/Applications/Utils/MeshEdit/CMakeLists.txt
+++ b/Applications/Utils/MeshEdit/CMakeLists.txt
@@ -81,6 +81,10 @@ add_executable(swapNodeCoordinateAxes swapNodeCoordinateAxes.cpp)
target_link_libraries(swapNodeCoordinateAxes MeshLib)
set_target_properties(swapNodeCoordinateAxes PROPERTIES FOLDER Utilities)
+add_executable(createQuadraticeMesh createQuadraticeMesh.cpp )
+target_link_libraries(createQuadraticeMesh MeshLib)
+set_target_properties(createQuadraticeMesh PROPERTIES FOLDER Utilities)
+
####################
### Installation ###
####################
@@ -90,6 +94,7 @@ install(TARGETS
checkMesh
CreateBoundaryConditionsAlongPolylines
createLayeredMeshFromRasters
+ createQuadraticeMesh
editMaterialID
ExtractSurface
MapGeometryToMeshSurface
diff --git a/Applications/Utils/MeshEdit/createQuadraticeMesh.cpp b/Applications/Utils/MeshEdit/createQuadraticeMesh.cpp
new file mode 100644
index 00000000000..47486c281f5
--- /dev/null
+++ b/Applications/Utils/MeshEdit/createQuadraticeMesh.cpp
@@ -0,0 +1,163 @@
+/**
+ * @copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/LICENSE.txt
+ *
+ */
+
+#include <memory>
+#include <string>
+
+#include <tclap/CmdLine.h>
+
+#include "Applications/ApplicationsLib/LogogSetup.h"
+
+#include "BaseLib/BuildInfo.h"
+#include "BaseLib/makeVectorUnique.h"
+
+#include "MeshLib/Mesh.h"
+#include "MeshLib/Node.h"
+#include "MeshLib/Elements/Element.h"
+#include "MeshLib/Elements/Line.h"
+#include "MeshLib/Elements/Quad.h"
+#include "MeshLib/MeshEditing/DuplicateMeshComponents.h"
+#include "MeshLib/Properties.h"
+#include "MeshLib/PropertyVector.h"
+
+#include "MeshLib/IO/readMeshFromFile.h"
+#include "MeshLib/IO/writeMeshToFile.h"
+
+namespace
+{
+
+struct Edge
+{
+ Edge(std::size_t id1, std::size_t id2)
+ : _id1(std::min(id1, id2)), _id2(std::max(id1, id2))
+ {}
+
+ bool operator==(Edge const& r) const
+ {
+ return (_id1 == r._id1 && _id2 == r._id2);
+ }
+
+ std::size_t _id1;
+ std::size_t _id2;
+
+ std::size_t _edge_id = 0;
+};
+
+bool operator< (Edge const& l, Edge const& r)
+{
+ return (l._id1 != r._id1) ? (l._id1 < r._id1) : l._id2 < r._id2;
+}
+
+template <typename T_ELEMENT>
+T_ELEMENT* createQuadraticElement(MeshLib::Element const* e, std::vector<Edge> const& vec_edges,
+ std::vector<MeshLib::Node*> const& vec_new_nodes, const std::size_t n_mesh_base_nodes)
+{
+ auto const n_all_nodes = T_ELEMENT::n_all_nodes;
+ auto const n_base_nodes = T_ELEMENT::n_base_nodes;
+ MeshLib::Node** nodes = new MeshLib::Node*[n_all_nodes];
+ for (unsigned i=0; i<e->getNumberOfBaseNodes(); i++)
+ nodes[i] = const_cast<MeshLib::Node*>(e->getNode(i));
+ for (unsigned i=0; i<e->getNumberOfEdges(); i++)
+ {
+ auto itr = std::find(vec_edges.begin(), vec_edges.end(),
+ Edge(e->getEdgeNode(i, 0)->getID(),
+ e->getEdgeNode(i, 1)->getID()));
+ assert(itr != vec_edges.end());
+ nodes[n_base_nodes + i] =
+ vec_new_nodes[n_mesh_base_nodes + itr->_edge_id];
+ }
+ return new T_ELEMENT(nodes);
+}
+
+
+static std::unique_ptr<MeshLib::Mesh> createQuadraticOrderMesh(MeshLib::Mesh const& org_mesh)
+{
+ std::vector<MeshLib::Node*> vec_new_nodes = MeshLib::copyNodeVector(org_mesh.getNodes());
+
+ // collect edges
+ std::vector<Edge> vec_edges;
+ for (MeshLib::Element const* e : org_mesh.getElements())
+ {
+ for (unsigned i=0; i<e->getNumberOfEdges(); i++)
+ {
+ auto node0 = e->getEdgeNode(i, 0);
+ auto node1 = e->getEdgeNode(i, 1);
+ vec_edges.emplace_back(node0->getID(), node1->getID());
+ }
+ }
+ BaseLib::makeVectorUnique(vec_edges);
+ for (std::size_t i=0; i<vec_edges.size(); i++)
+ vec_edges[i]._edge_id = i;
+ INFO("Found %d edges in the mesh", vec_edges.size());
+
+ // create mid-point nodes
+ double coords[3];
+ for (Edge const& edge : vec_edges)
+ {
+ auto const& node0 = *vec_new_nodes[edge._id1];
+ auto const& node1 = *vec_new_nodes[edge._id2];
+ for (unsigned i=0; i<3; i++)
+ coords[i] = (node0[i] + node1[i]) * 0.5;
+ vec_new_nodes.push_back(new MeshLib::Node(coords));
+ }
+
+ // create new elements with the quadratic nodes
+ std::vector<MeshLib::Element*> vec_new_eles;
+ for (MeshLib::Element const* e : org_mesh.getElements())
+ {
+ if (e->getCellType() == MeshLib::CellType::LINE2)
+ {
+ vec_new_eles.push_back(createQuadraticElement<MeshLib::Line3>(
+ e, vec_edges, vec_new_nodes, org_mesh.getNumberOfNodes()));
+ }
+ else if (e->getCellType() == MeshLib::CellType::QUAD4)
+ {
+ vec_new_eles.push_back(createQuadraticElement<MeshLib::Quad8>(
+ e, vec_edges, vec_new_nodes, org_mesh.getNumberOfNodes()));
+ }
+ else
+ {
+ OGS_FATAL("Mesh element type %s is not supported", MeshLib::CellType2String(e->getCellType()).c_str());
+ }
+ }
+
+ std::unique_ptr<MeshLib::Mesh> new_mesh(
+ new MeshLib::Mesh(org_mesh.getName(), vec_new_nodes, vec_new_eles,
+ org_mesh.getProperties().excludeCopyProperties(
+ std::vector<MeshLib::MeshItemType>(
+ 1, MeshLib::MeshItemType::Node))));
+ return new_mesh;
+}
+
+} // no named namespace
+
+int main(int argc, char *argv[])
+{
+ ApplicationsLib::LogogSetup logog_setup;
+
+ TCLAP::CmdLine cmd("Create quadratic order mesh", ' ', BaseLib::BuildInfo::git_describe);
+ TCLAP::ValueArg<std::string> input_arg("i", "input-mesh-file","input mesh file",true,"","string");
+ cmd.add( input_arg );
+ TCLAP::ValueArg<std::string> output_arg("o", "output-mesh-file","output mesh file",true,"","string");
+ cmd.add( output_arg );
+ cmd.parse( argc, argv );
+
+ std::unique_ptr<MeshLib::Mesh> mesh(
+ MeshLib::IO::readMeshFromFile(input_arg.getValue()));
+ if (!mesh)
+ return EXIT_FAILURE;
+
+ INFO("Create a quadratic order mesh");
+ auto new_mesh(createQuadraticOrderMesh(*mesh));
+
+ INFO("Save the new mesh into a file");
+ MeshLib::IO::writeMeshToFile(*new_mesh, output_arg.getValue());
+
+ return EXIT_SUCCESS;
+}
--
GitLab