The 3D export of the unstructured grid to VTK is done for all VTK (b568a87e) · Commits · Tomáš Jakubec / GTMesh

Unstructured_mesh/MeshWriter.h

+54 −9

Original line number	Diff line number	Diff line
		#ifndef MESHWRITER_H
		#define MESHWRITER_H
		#include "MeshNativeType.h"
		#include "Vertex.h"
		#include "MeshElement.h"

		template<unsigned int MeshDimension, typename IndexType, typename Real>
		class MeshWriter{

		};

		protected:
		/**
		* @brief The MeshHash struct<HR>
		* A container to store cumulative data about a mesh
		* to recognize changes in the mesh.
		* It uses number of elements and sum of vertices coordinations.
		*/
		struct MeshHash {
		IndexType numberOfElements = 0;
		Real totalVert = 0;

		template <typename IndexType, typename Real>
		class MeshWriter<2, IndexType, Real> {
		public:
		using type = MeshNativeType<2>;
		bool operator== (const MeshHash& rhs) const {
		return numberOfElements == rhs.numberOfElements &&
		totalVert == rhs.totalVert;
		}

		bool operator!= (const MeshHash& rhs) const {
		return !((*this)==rhs);
		}
		};

		template <typename IndexType, typename Real>
		class MeshWriter<3, IndexType, Real> {
		private:
		template<unsigned int Dim, typename Dummy = void>
		struct sumOfMeshElements{
		template<unsigned int ...Reserve>
		static IndexType sum(MeshElements<MeshDimension, IndexType, Real, Reserve...>& mesh){
		return mesh.template getElements<Dim>().size() + sumOfMeshElements<Dim - 1>::sum(mesh);
		}
		};
		template<typename Dummy>
		struct sumOfMeshElements<0, Dummy>{
		template<unsigned int ...Reserve>
		static IndexType sum(MeshElements<MeshDimension, IndexType, Real, Reserve...>& mesh){
		return mesh.template getElements<0>().size();
		}
		};
		public:
		using type = MeshNativeType<3>;
		using type = MeshNativeType<MeshDimension>;

		template<unsigned int ...Reserve>
		static MeshHash computeHash(MeshElements<MeshDimension, IndexType, Real, Reserve...>& mesh){
		MeshHash res;
		// A vector of ones that simplifies the sum of coordinates
		Vertex<MeshDimension, Real> ones;
		for(unsigned int dim = 0; dim < MeshDimension; dim++){
		ones[dim] = 1.0;
		}

		for(IndexType i = 0; i < mesh.getVertices().size(); i++) {
		res.totalVert += mesh.getVertices().at(i) * ones;
		}

		res.numberOfElements = sumOfMeshElements<MeshDimension>::sum(mesh);
		return res;
		}
		};


		#endif // MESHWRITER_H

Unstructured_mesh/Unstructured_mesh.pro.user

+1 −1

Original line number	Diff line number	Diff line
		<?xml version="1.0" encoding="UTF-8"?>
		<!DOCTYPE QtCreatorProject>
		<!-- Written by QtCreator 4.10.0, 2019-10-09T22:44:25. -->
		<!-- Written by QtCreator 4.10.0, 2019-10-10T18:51:05. -->
		<qtcreator>
		<data>
		<variable>EnvironmentId</variable>

Unstructured_mesh/VTKMeshWriter.h

+201 −86

Original line number	Diff line number	Diff line
		@@ -31,69 +31,130 @@ public:
		dataName << "\nASCII\nDATASET UNSTRUCTURED_GRID" << std::endl;
		}

		void writeToStream(std::ostream& ost,
		MeshElements<2, IndexType, Real, Reserve...>& mesh,
		MeshDataContainer<typename writer::type::ElementType, 2> cellTypes){
		// first write verices
		ost << "POINTS " << mesh.getVertices().size() <<
		" double" << std::endl;
		/**
		* @brief lastHash<HR>
		* The hash of the last written mesh.
		*/
		typename writer::MeshHash lastHash;
		/**
		* @brief cellVert<HR>
		* Vertices of all cells in correct order for vtk export.
		*/
		MeshDataContainer<std::vector<IndexType>, 3> cellVert;

		/**
		* @brief totalNumberOfWrittenElements<HR>
		* Information required by VTK format.
		*/
		IndexType totalNumberOfWrittenElements = 0;


		/**
		* @brief indexCell<HR>
		* This funcion stores indexes of vertices in correct order
		* in output vector verticesIndexed
		* @param mesh the structure of the mesh
		* @param face the face of the mesh to be indexed
		* @param cell the cell which are the vertices being indexed to
		* @param faceEdgeOri the orientation of the edges to the faces
		* @param verticesIndexed output vector of indexed vertices
		*/
		void indexCell(MeshElements<2, IndexType, Real, Reserve...>& mesh,
		typename MeshElements<2, IndexType, Real, Reserve...>::Cell& cell,
		std::vector<IndexType>& verticesIndexed){


		IndexType tmpEdgeIndex = cell.getBoundaryElementIndex();

		for(auto vert : mesh.getVertices()) {
		ost << vert[0] << " " << vert[1] << " 0.0\n";
		IndexType nextVertex = INVALID_INDEX(IndexType);
		if (mesh.getEdges().at(tmpEdgeIndex).getCellLeftIndex() == cell.getIndex()){
		verticesIndexed.push_back(mesh.getEdges().at(tmpEdgeIndex).getVertexAIndex());
		nextVertex = mesh.getEdges().at(tmpEdgeIndex).getVertexBIndex();
		} else {
		verticesIndexed.push_back(mesh.getEdges().at(tmpEdgeIndex).getVertexBIndex());
		nextVertex = mesh.getEdges().at(tmpEdgeIndex).getVertexAIndex();
		}
		ost << std::endl;

		auto cellVert = ::MeshConnections<2,0>::connections(mesh);
		int cnt = 0;
		for (auto verts : cellVert.template getDataByPos<0>()){
		cnt += verts.size() + 1;
		}
		ost << "CELLS " << mesh.getCells().size() << ' ' << cnt << std::endl;
		size_t lastWrittenEdge = tmpEdgeIndex;

		for (auto cell : mesh.getCells()){
		ost << cellVert.at(cell).size() << " ";
		do {
		if ((lastWrittenEdge != tmpEdgeIndex)){
		if (mesh.getEdges().at(tmpEdgeIndex).getVertexBIndex() == nextVertex) {

		size_t tmpEdgeIndex = cell.getBoundaryElementIndex();
		nextVertex = mesh.getEdges().at(tmpEdgeIndex).getVertexAIndex();

		size_t lastWrittenVertex = INVALID_INDEX(IndexType);
		verticesIndexed.push_back(mesh.getEdges().at(tmpEdgeIndex).getVertexBIndex());

		if (mesh.getEdges().at(tmpEdgeIndex).getCellLeftIndex() == cell.getIndex()){
		ost << mesh.getEdges().at(tmpEdgeIndex).getVertexAIndex() << ' ';
		lastWrittenVertex = mesh.getEdges().at(tmpEdgeIndex).getVertexBIndex();
		} else {
		ost << mesh.getEdges().at(tmpEdgeIndex).getVertexBIndex() << ' ';
		lastWrittenVertex = mesh.getEdges().at(tmpEdgeIndex).getVertexAIndex();
		lastWrittenEdge = tmpEdgeIndex;

		} else if (mesh.getEdges().at(tmpEdgeIndex).getVertexAIndex() == nextVertex){

		nextVertex = mesh.getEdges().at(tmpEdgeIndex).getVertexBIndex();

		verticesIndexed.push_back(mesh.getEdges().at(tmpEdgeIndex).getVertexAIndex());

		lastWrittenEdge = tmpEdgeIndex;
		}

		size_t lastWrittenEdge = tmpEdgeIndex;
		size_t verticesWritten = 1;
		do {
		if ((lastWrittenEdge != tmpEdgeIndex \|\| lastWrittenEdge != INVALID_INDEX(IndexType) ) &&
		(lastWrittenVertex == INVALID_INDEX(IndexType) \|\| (lastWrittenVertex == mesh.getEdges().at(tmpEdgeIndex).getVertexBIndex() \|\| lastWrittenVertex == mesh.getEdges().at(tmpEdgeIndex).getVertexAIndex()))) {
		if (mesh.getEdges().at(tmpEdgeIndex).getCellLeftIndex() == cell.getIndex()) {

		lastWrittenVertex = mesh.getEdges().at(tmpEdgeIndex).getVertexBIndex();
		}
		tmpEdgeIndex = mesh.getEdges().at(tmpEdgeIndex).getNextBElem(cell.getIndex());

		ost << mesh.getEdges().at(tmpEdgeIndex).getVertexBIndex();
		} while (nextVertex != verticesIndexed.at(0));

		} else {
		}

		lastWrittenVertex = mesh.getEdges().at(tmpEdgeIndex).getVertexAIndex();
		void indexMesh(MeshElements<2, IndexType, Real, Reserve...>& mesh){
		typename writer::MeshHash curHash = writer::computeHash(mesh);

		ost << mesh.getEdges().at(tmpEdgeIndex).getVertexAIndex();
		// if the mesh is the same as it was, return
		if (lastHash == curHash){
		return;
		}
		lastHash = curHash;

		std::vector<IndexType> vertIndex;
		for (auto& cell : mesh.getCells()) {
		vertIndex.clear();
		vertIndex.reserve(4);
		indexCell(mesh, cell, vertIndex);
		cellVert.template getDataByPos<0>().push_back(vertIndex);
		}
		cellVert.template getDataByPos<0>().shrink_to_fit();

		lastWrittenEdge = tmpEdgeIndex;
		totalNumberOfWrittenElements = 0;
		for (auto verts : cellVert.template getDataByPos<0>()){
		totalNumberOfWrittenElements += verts.size() + 1;
		}

		verticesWritten++;

		ost << ' ';
		}

		void writeToStream(std::ostream& ost,
		MeshElements<2, IndexType, Real, Reserve...>& mesh,
		MeshDataContainer<typename writer::type::ElementType, 2> cellTypes){

		indexMesh(mesh);
		// first write verices
		ost << "POINTS " << mesh.getVertices().size() <<
		" double" << std::endl;

		for(auto vert : mesh.getVertices()) {
		ost << vert[0] << " " << vert[1] << " 0.0\n";
		}
		tmpEdgeIndex = mesh.getEdges().at(tmpEdgeIndex).getNextBElem(cell.getIndex());
		ost << std::endl;


		ost << "CELLS " << mesh.getCells().size() << ' ' << totalNumberOfWrittenElements << std::endl;

		} while (verticesWritten < cellVert.at(cell).size());
		for (const auto& verts : cellVert.template getDataByPos<0>()){
		ost << verts.size() << ' ';
		for (IndexType i = 0; i < verts.size(); i++) {
		ost << verts[i];
		if (i < verts.size() - 1) {
		ost << ' ';
		}
		}
		ost << std::endl;
		}
		ost << std::endl;
		@@ -126,29 +187,52 @@ public:
		dataName << "\nASCII\nDATASET UNSTRUCTURED_GRID" << std::endl;
		}


		std::vector<IndexType> writeFace(std::ostream& ost,
		MeshElements<3, IndexType, Real, Reserve...>& mesh,
		// indexing of the mesh
		/**
		* @brief lastHash<HR>
		* The hash of the last written mesh.
		*/
		typename writer::MeshHash lastHash;
		/**
		* @brief cellVert<HR>
		* Vertices of all cells in correct order for vtk export.
		*/
		MeshDataContainer<std::vector<IndexType>, 3> cellVert;

		/**
		* @brief totalNumberOfWrittenElements<HR>
		* Information required by VTK format.
		*/
		IndexType totalNumberOfWrittenElements = 0;

		/**
		* @brief indexFace<HR>
		* This funcion stores indexes of vertices in correct order
		* in output vector verticesIndexed
		* @param mesh the structure of the mesh
		* @param face the face of the mesh to be indexed
		* @param cell the cell which are the vertices being indexed to
		* @param faceEdgeOri the orientation of the edges to the faces
		* @param verticesIndexed output vector of indexed vertices
		*/
		void indexFace(MeshElements<3, IndexType, Real, Reserve...>& mesh,
		typename MeshElements<3, IndexType, Real, Reserve...>::Face& face,
		typename MeshElements<3, IndexType, Real, Reserve...>::Cell& cell,
		MeshDataContainer<std::vector<bool>, 2>& faceEdgeOri){
		MeshDataContainer<std::vector<bool>, 2>& faceEdgeOri,
		std::vector<IndexType>& verticesIndexed){

		std::vector<IndexType> verticesWritten;
		verticesWritten.reserve(face.getSubelements().getNumberOfSubElements());

		IndexType startVertex = INVALID_INDEX(IndexType);
		IndexType nextVertex = INVALID_INDEX(IndexType);
		if (faceEdgeOri[face][0] == true && cell.getIndex() == face.getCellLeftIndex()){ // the edge is left to the face
		startVertex = mesh.getEdges().at(face.getSubelements()[0].index).getVertexBIndex();
		ost << startVertex << ' ';
		nextVertex = mesh.getEdges().at(face.getSubelements()[0].index).getVertexAIndex();
		} else {
		startVertex = mesh.getEdges().at(face.getSubelements()[0].index).getVertexAIndex();
		ost << startVertex << ' ';
		nextVertex = mesh.getEdges().at(face.getSubelements()[0].index).getVertexBIndex();
		}

		verticesWritten.push_back(startVertex);
		verticesIndexed.push_back(startVertex);

		IndexType lastWrittenEdge = face.getSubelements()[0].index;
		while (startVertex != nextVertex){
		@@ -158,51 +242,52 @@ public:
		if (edge.getIndex() != lastWrittenEdge) {
		if (edge.getVertexAIndex() == nextVertex) {
		lastWrittenEdge = edge.getIndex();
		ost << edge.getVertexAIndex() << ' ';
		verticesWritten.push_back(edge.getVertexAIndex());
		verticesIndexed.push_back(edge.getVertexAIndex());
		nextVertex = edge.getVertexBIndex();
		} else if (edge.getVertexBIndex() == nextVertex) {
		lastWrittenEdge = edge.getIndex();
		ost << edge.getVertexBIndex() << ' ';
		verticesWritten.push_back(edge.getVertexBIndex());
		verticesIndexed.push_back(edge.getVertexBIndex());
		nextVertex = edge.getVertexAIndex();
		}
		}
		}
		}
		return verticesWritten;

		}

		void writeToStream(std::ostream& ost,
		MeshElements<3, IndexType, Real, Reserve...>& mesh,

		void indexMesh(MeshElements<3, IndexType, Real, Reserve...>& mesh,
		MeshDataContainer<typename writer::type::ElementType, 3> cellTypes){
		// first write verices
		ost << "POINTS " << mesh.getVertices().size() <<
		" double" << std::endl;
		typename writer::MeshHash curHash = writer::computeHash(mesh);

		for(auto vert : mesh.getVertices()) {
		ost << vert[0] << ' ' << vert[1] << ' ' << vert[2] <<"\n";
		// if the mesh is the same as it was, return
		if (lastHash == curHash){
		return;
		}
		ost << std::endl;

		lastHash = curHash;

		DBGMSG("indexing mesh");
		// write cells of the mesh
		// prepare connections
		auto cellVert = MeshConnections<3,0>::connections(mesh);
		auto cellFace = MeshConnections<3,2>::connections(mesh);
		auto faceVert = MeshConnections<2,0>::connections(mesh);
		// prepare orientation for correct export
		// this is very expensive procedure
		auto faceEdgeOri = edgesOrientation(mesh);

		int cnt = 0;
		for (auto verts : cellVert.template getDataByPos<0>()){
		cnt += verts.size() + 1;
		}
		ost << "CELLS " << mesh.getCells().size() << ' ' << cnt << std::endl;

		auto faceEdgeOri = edgesOrientation(mesh);
		totalNumberOfWrittenElements = cnt;


		std::vector<IndexType> vertWrit;
		for (auto cell : mesh.getCells()){
		ost << cellVert.at(cell).size() << " ";
		vertWrit.clear();
		vertWrit.reserve(cellVert[cell].size());

		DBGVAR(cell.getIndex());

		@@ -212,7 +297,7 @@ public:
		// write vertices of one face
		auto& face = mesh.getFaces().at(cell.getBoundaryElementIndex());

		std::vector<IndexType> vertWrit = writeFace(ost,mesh,face, cell, faceEdgeOri);
		indexFace(mesh,face, cell, faceEdgeOri, vertWrit);

		for (IndexType index : cellVert[cell]) {
		bool vertOK = true;
		@@ -222,7 +307,6 @@ public:
		}
		}
		if (vertOK){
		ost << index;
		vertWrit.push_back(index);
		}
		}
		@@ -233,12 +317,12 @@ public:
		typename MeshElements<3, IndexType, Real, Reserve...>::Face* face = nullptr;
		// search for the base face
		for (IndexType faceIndex : cellFace[cell]){
		if (faceVert.template getDataByPos<0>().at(faceIndex).size() == 4){
		if (faceVert.template getDataByPos<0>().at(faceIndex).size() > 3){
		face = &mesh.getFaces().at(faceIndex);
		}
		}

		std::vector<IndexType> vertWrit = writeFace(ost, mesh, *face, cell, faceEdgeOri);
		indexFace(mesh, *face, cell, faceEdgeOri, vertWrit);
		// write the last vertex
		for (IndexType index : cellVert.at(cell)) {
		bool vertOK = true;
		@@ -248,7 +332,6 @@ public:
		}
		}
		if (vertOK){
		ost << index;
		vertWrit.push_back(index);
		}
		}
		@@ -266,12 +349,14 @@ public:
		}
		}
		DBGVAR(face->getIndex());
		std::vector<IndexType> vertWrit = writeFace(ost, mesh, *face, cell, faceEdgeOri);
		indexFace(mesh, *face, cell, faceEdgeOri, vertWrit);
		// write vertices of the oposite triangular side
		DBGVAR(vertWrit);
		for (IndexType index : vertWrit) {

		for (IndexType i = 0; i < vertWrit.size(); i++) {
		IndexType index = vertWrit[i];

		MeshRun<3,3,1,3,false, true>::run(mesh, cell.getIndex(),cell.getIndex(),
		[&ost,&mesh,&index,&vertWrit](IndexType, IndexType edgeIndex){
		[&mesh,&index,&vertWrit](IndexType, IndexType edgeIndex){
		auto& edge = mesh.getEdges().at(edgeIndex);

		if (edge.getVertexAIndex() == index){
		@@ -282,7 +367,6 @@ public:
		}
		}
		if(edgeOK){
		ost << edge.getVertexBIndex() << ' ';
		vertWrit.push_back(edge.getVertexBIndex());
		}
		}
		@@ -295,13 +379,13 @@ public:
		}
		}
		if(edgeOK){
		ost << edge.getVertexAIndex() << ' ';
		vertWrit.push_back(edge.getVertexAIndex());
		}
		}
		}
		);
		if (vertWrit.size() == 6) {
		DBGVAR(vertWrit);
		break;
		}
		}
		@@ -311,11 +395,12 @@ public:
		// write vertices of one face
		auto& face = mesh.getFaces().at(cell.getBoundaryElementIndex());

		std::vector<IndexType> vertWrit = writeFace(ost, mesh, face, cell, faceEdgeOri);
		indexFace(mesh, face, cell, faceEdgeOri, vertWrit);
		// write vertices of the oposite triangular side
		for (IndexType index : vertWrit) {
		for (IndexType i = 0; i < vertWrit.size(); i++) {
		IndexType index = vertWrit[i];
		MeshRun<3,3,1,3,false, true>::run(mesh, cell.getIndex(),cell.getIndex(),
		[&ost,&mesh,&index,&vertWrit](IndexType, IndexType edgeIndex){
		[&mesh,&index,&vertWrit](IndexType, IndexType edgeIndex){
		auto& edge = mesh.getEdges().at(edgeIndex);

		if (edge.getVertexAIndex() == index){
		@@ -326,7 +411,6 @@ public:
		}
		}
		if(edgeOK){
		ost << edge.getVertexBIndex() << ' ';
		vertWrit.push_back(edge.getVertexBIndex());
		}
		}
		@@ -339,7 +423,6 @@ public:
		}
		}
		if(edgeOK){
		ost << edge.getVertexAIndex() << ' ';
		vertWrit.push_back(edge.getVertexAIndex());
		}
		}
		@@ -352,9 +435,41 @@ public:
		}break;
		default: throw std::runtime_error("it is not possible yet to write any object into VTK");
		}
		this->cellVert.template getDataByPos<0>().push_back(vertWrit);
		}
		}



		void writeToStream(std::ostream& ost,
		MeshElements<3, IndexType, Real, Reserve...>& mesh,
		MeshDataContainer<typename writer::type::ElementType, 3> cellTypes){
		// first write verices
		ost << "POINTS " << mesh.getVertices().size() <<
		" double" << std::endl;

		for(auto vert : mesh.getVertices()) {
		ost << vert[0] << ' ' << vert[1] << ' ' << vert[2] <<"\n";
		}
		ost << std::endl;


		// write cells of the mesh
		indexMesh(mesh, cellTypes);
		ost << "CELLS " << mesh.getCells().size() << ' ' << totalNumberOfWrittenElements << std::endl;


		for (const auto& verts : cellVert.template getDataByPos<0>()){
		ost << verts.size() << ' ';
		for (IndexType i = 0; i < verts.size(); i++) {
		ost << verts[i];
		if (i < verts.size() - 1) {
		ost << ' ';
		}
		}
		ost << std::endl;
		}

		ost << std::endl;

		ost << "CELL_TYPES " << mesh.getCells().size() << std::endl;

Unstructured_mesh/main.cpp

+14 −3

Original line number	Diff line number	Diff line
		@@ -552,6 +552,9 @@ void testMesh3D() {
		MeshDataContainer<MeshNativeType<3>::ElementType,3> types(mesh3, MeshNativeType<3>::ElementType::WEDGE);

		VTKMeshWriter<3, size_t, double, 6> writer;
		writer.indexMesh(mesh3, types);

		DBGVAR(writer.cellVert.getDataByPos<0>());
		ofstream out3D("3D_test_mesh_two_prisms.vtk");
		writer.writeHeader(out3D, "test data");
		writer.writeToStream(out3D, mesh3, types);
		@@ -590,6 +593,14 @@ void test3DMeshDeformedPrisms() {
		for(double cellM : measures.getDataByDim<3>()) {
		DBGVAR(cellM);
		}


		MeshDataContainer<MeshNativeType<3>::ElementType,3> types(mesh3, MeshNativeType<3>::ElementType::WEDGE);

		VTKMeshWriter<3, size_t, double, 6> writer;
		ofstream out3D("3D_test_mesh_two_deformed_prisms.vtk");
		writer.writeHeader(out3D, "test data");
		writer.writeToStream(out3D, mesh3, types);
		}


		@@ -694,13 +705,13 @@ DBGVAR(mesh.getVertices().size(),mesh.getEdges().size(), mesh.getFaces().size(),
		int main()
		{
		//testMesh2D();
		//testMesh2DLoadAndWrite();
		//testMesh3D();
		testMesh2DLoadAndWrite();
		testMesh3D();
		//test3DMeshDeformedPrisms();
		//testMeshDataContainer();
		//UnstructuredMesh<5, size_t, double, 6,5,4> m;
		//m.ComputeElementMeasures();
		test3DMeshLoad();
		//test3DMeshLoad();


		}