From 98a0d3a871e713bd5bddb02e0074994a0d83a966 Mon Sep 17 00:00:00 2001
From: Koen Hillewaert <koen.hillewaert@cenaero.be>
Date: Wed, 20 Apr 2016 20:42:36 +0000
Subject: [PATCH] added correction/averaging functions for periodic connections
---
Mesh/CMakeLists.txt | 1 +
Mesh/meshGEdge.cpp | 34 +++++++++++++++++++++++++
Mesh/meshGEntity.cpp | 60 ++++++++++++++++++++++++++++++++++++++++++++
Mesh/meshGFace.cpp | 55 ++++++++++++++++++++++++++++++++++++++++
4 files changed, 150 insertions(+)
create mode 100644 Mesh/meshGEntity.cpp
diff --git a/Mesh/CMakeLists.txt b/Mesh/CMakeLists.txt
index 7cf93c40fc..7088580950 100644
--- a/Mesh/CMakeLists.txt
+++ b/Mesh/CMakeLists.txt
@@ -9,6 +9,7 @@ set(SRC
delaunay_refinement.cpp
Generator.cpp
meshGEdge.cpp
+ meshGEntity.cpp
meshGEdgeExtruded.cpp
meshGFace.cpp
meshGFaceElliptic.cpp
diff --git a/Mesh/meshGEdge.cpp b/Mesh/meshGEdge.cpp
index 59f2a8894b..f290f990c7 100644
--- a/Mesh/meshGEdge.cpp
+++ b/Mesh/meshGEdge.cpp
@@ -303,6 +303,40 @@ void copyMesh(GEdge *from, GEdge *to, int direction)
}
}
+bool correctPeriodicity(MVertex*,MVertex*,const std::vector<double>&);
+bool correctVertexPeriodicity(GEntity*);
+
+
+bool correctPeriodicity(GEdge* tgt) {
+
+ if (correctVertexPeriodicity(tgt)) {
+
+ bool check = true;
+
+ GEdge* src = dynamic_cast<GEdge*> (tgt->meshMaster());
+
+ if (!src) throw;
+
+ const std::vector<double>& tfo = src->affineTransform;
+
+ std::vector<MLine*>::iterator sIter = src->lines.begin();
+ std::vector<MLine*>::iterator tIter = tgt->lines.begin();
+
+ for (;sIter!=src->lines.end();++sIter,++tIter) {
+
+ MLine* src = *sIter;
+ MLine* tgt = *sIter;
+
+ for (int iVtx=2;iVtx<tgt->getNumVertices();iVtx++) {
+ check = check &&
+ correctPeriodicity(tgt->getVertex(iVtx),src->getVertex(iVtx),tfo);
+ }
+ }
+ }
+ return false;
+}
+
+
void deMeshGEdge::operator() (GEdge *ge)
{
if(ge->geomType() == GEntity::DiscreteCurve && !CTX::instance()->meshDiscrete)
diff --git a/Mesh/meshGEntity.cpp b/Mesh/meshGEntity.cpp
new file mode 100644
index 0000000000..2a8880829e
--- /dev/null
+++ b/Mesh/meshGEntity.cpp
@@ -0,0 +1,60 @@
+// Gmsh - Copyright (C) 1997-2016 C. Geuzaine, J.-F. Remacle
+//
+// See the LICENSE.txt file for license information. Please report all
+// bugs and problems to the public mailing list <gmsh@onelab.info>.
+//
+// Contributor(s):
+// Koen Hillewaert
+
+#include "GEntity.h"
+#include "MVertex.h"
+
+// -----------------------------------------------------------------------------
+
+bool correctPeriodicity(MVertex* tgt,
+ MVertex* src,
+ const std::vector<double>& tfo) {
+
+ double ps[4] = {src->x(),src->y(),src->z(),1};
+ double pt[4] = {tgt->x(),tgt->y(),tgt->z(),1};
+ int idx = 0;
+ for(int i = 0; i < 4; i++)
+ for(int j = 0; j < 4; j++)
+ ps[j] += tfo[idx++] * pt[i];
+
+ for (int i=0;i<4;i++) ps[i] *= 0.5;
+ for (int i=0;i<4;i++) pt[i] = 0;
+
+ src->x() = ps[0];
+ src->y() = ps[1];
+ src->z() = ps[2];
+
+ idx = 0;
+ for(int i = 0; i < 4; i++)
+ for(int j = 0; j < 4; j++)
+ pt[i] += tfo[idx++] * ps[j];
+
+ tgt->x() = pt[0];
+ tgt->y() = pt[1];
+ tgt->z() = pt[2];
+
+ return true;
+}
+
+// -----------------------------------------------------------------------------
+
+bool correctVertexPeriodicity(GEntity* tgt) {
+
+ if (tgt->meshMaster() == NULL) return false;
+ bool check = true;
+ std::map<MVertex*,MVertex*>::iterator vIter = tgt->correspondingVertices.begin();
+ for (;vIter!=tgt->correspondingVertices.end();++vIter) {
+ check = check && correctPeriodicity(vIter->first,
+ vIter->second,
+ tgt->affineTransform);
+
+ }
+ return check;
+}
+
+// -----------------------------------------------------------------------------
diff --git a/Mesh/meshGFace.cpp b/Mesh/meshGFace.cpp
index 8b2f3cbd60..5cd6ec89ec 100644
--- a/Mesh/meshGFace.cpp
+++ b/Mesh/meshGFace.cpp
@@ -25,6 +25,7 @@
#include "MLine.h"
#include "MTriangle.h"
#include "MQuadrangle.h"
+#include "MElementCut.h"
#include "CenterlineField.h"
#include "meshGFaceElliptic.h"
#include "Context.h"
@@ -379,6 +380,60 @@ static void copyMesh(GFace *source, GFace *target)
}
}
+
+bool correctPeriodicity(MVertex*,MVertex*,const std::vector<double>&);
+bool correctVertexPeriodicity(GEntity*);
+
+template<class MElementType>
+bool correctPeriodicity(std::vector<MElementType*>& tgtList,
+ std::vector<MElementType*>& srcList,
+ const std::vector<double>& tfo,
+ std::set<MVertex*>& corrected) {
+
+ bool check = true;
+
+ typename std::vector<MElementType*>::iterator sIter = srcList.begin();
+ typename std::vector<MElementType*>::iterator tIter = tgtList.begin();
+
+ for (;sIter!=srcList.end();++sIter,++tIter) {
+
+ MElementType* src = *sIter;
+ MElementType* tgt = *sIter;
+
+ for (int iVtx=0;iVtx<tgt->getNumVertices();iVtx++) {
+
+ MVertex* tgtVtx = tgt->getVertex(iVtx);
+ MVertex* srcVtx = src->getVertex(iVtx);
+
+ if (corrected.find(tgtVtx) != corrected.end()) {
+
+ check = check && correctPeriodicity(tgtVtx,srcVtx,tfo);
+ if (check) corrected.insert(tgtVtx);
+ }
+ }
+ }
+ return check;
+}
+
+bool correctPeriodicity(GFace* tgt) {
+
+ if (!correctVertexPeriodicity(tgt)) return false;
+
+
+ GFace* src = dynamic_cast<GFace*> (tgt->meshMaster());
+ if (!src) throw;
+
+ const std::vector<double>& tfo = src->affineTransform;
+
+ std::set<MVertex*> corr(tgt->mesh_vertices.begin(),tgt->mesh_vertices.end());
+
+ if (!correctPeriodicity(tgt->triangles ,src->triangles ,tfo,corr)) return false;
+ if (!correctPeriodicity(tgt->quadrangles,src->quadrangles,tfo,corr)) return false;
+ if (!correctPeriodicity(tgt->polygons ,src->polygons ,tfo,corr)) return false;
+
+ return true;
+}
+
void fourthPoint(double *p1, double *p2, double *p3, double *p4)
{
double c[3];
--
GitLab