From 70194ee52a26227de52334d2f47b4569c8ce2c4a Mon Sep 17 00:00:00 2001
From: Jean-Francois Remacle <jean-francois.remacle@uclouvain.be>
Date: Wed, 18 May 2016 10:21:36 +0000
Subject: [PATCH] discrete face for PA
---
Geo/discreteFace.cpp | 100 ++++++++++++++++++++++++++++-
Mesh/Generator.cpp | 7 +-
Mesh/meshGRegionRelocateVertex.cpp | 25 +++++---
Mesh/meshGRegionRelocateVertex.h | 1 +
Mesh/yamakawa.cpp | 19 +++---
Parser/Gmsh.y | 1 +
6 files changed, 129 insertions(+), 24 deletions(-)
diff --git a/Geo/discreteFace.cpp b/Geo/discreteFace.cpp
index 84f5d987e7..0e2b65093e 100644
--- a/Geo/discreteFace.cpp
+++ b/Geo/discreteFace.cpp
@@ -301,6 +301,105 @@ void discreteFace::checkAndFixOrientation(){
} // end while
}
+// should be put in discreteEdge
+void splitDiscreteEdge ( discreteEdge *de , MVertex *v , discreteEdge *replacements[2]) {
+ if (v->onWhat() != de) {
+ Msg::Fatal ("cannot split a discreteEdge on a point that is NOT classified on it");
+ }
+
+ MVertex *vstart = de->getBeginVertex()->mesh_vertices[0];
+ MVertex *vend = de->getEndVertex()->mesh_vertices[0];
+
+ // We create a new Model vertex and we classify
+ discreteVertex *newV = new discreteVertex (de->model(),NEWPOINT());
+ newV->mesh_vertices.push_back (v);
+ v->setEntity (newV);
+ de->model()->add(newV);
+
+ // We create 2 new Model edges and we classify vertices
+ // The creation of discrete edge geometries have already been done (or not)
+ // FIXME !!!!!
+ discreteEdge *newE[2];
+ newE[0] = new discreteEdge (de->model(),NEWLINE(),de->getBeginVertex(),newV);
+ newE[1] = new discreteEdge (de->model(),NEWLINE(),newV, de->getEndVertex());
+ de->model()->add(newE[0]);
+ de->model()->add(newE[1]);
+
+ int current = 0;
+
+ for (unsigned int i=0;i<de->lines.size();i++){
+ MLine *l = de->lines[i];
+ newE[current]->lines.push_back(l);
+ if (l->getVertex(0) != vstart && l->getVertex(0) != v) {
+ l->getVertex(0)->setEntity(newE[current]);
+ newE[current]->mesh_vertices.push_back(l->getVertex(0));
+ }
+ if (l->getVertex(1) == vend) current ++;
+ }
+ de->mesh_vertices.clear();
+ de->lines.clear();
+
+ // We replace de by its 2 parts in every face that is adjacent to de
+ std::list<GFace*> faces = de->faces();
+ for (std::list<GFace*>::iterator it = faces.begin(); it != faces.end(); ++it){
+ GFace *gf = *it;
+ if (gf->geomType() == GEntity::DiscreteSurface){
+ discreteFace *df = static_cast<discreteFace*> (gf);
+ if (df){
+ std::vector<int> tagEdges;
+ std::list<GEdge*> edges = df->edges();
+ for (std::list<GEdge*>::iterator it2 = edges.begin(); it2 != edges.end(); ++it2){
+ if (*it2 == de){
+ tagEdges.push_back (newE[0]->tag());
+ tagEdges.push_back (newE[1]->tag());
+ }
+ else tagEdges.push_back (de->tag());
+ }
+
+ df->setBoundEdges (df->model(), tagEdges);
+ }
+ else {
+ Msg::Fatal("this only applies to discrete geometries");
+ }
+ }
+ }
+}
+
+void addInternalBoundaries ( discreteFace *df) {
+ // create new discreteEdges that split the domain
+
+}
+
+/*
+void splitDiscreteFace ( discreteFace *df , std::vector<triangulation*> &partition, std::vector<discreteEdge*> &internalBoundaries) {
+ std::map<MEdge, GEdge*> _dictionnary;
+
+ std::vector<discreteEdge*> allBoundaries = internalBoundaries;
+ std::list<GEdge*> edges = df->edges();
+ allBoundaries.insert(allBoundaries.begin(), edges.begin(), edges.end());
+
+ for (unsigned int i=0;i<allBoundaries.size();i++){
+ for (unsigned int j=0;i<allBoundaries[i]->lines.size();j++){
+ MLine *l = allBoundaries[i]->lines[j];
+ MEdge e (l->getVertex(0),l->getVertex(1));
+ _dictionnary[e] = allBoundaries[i];
+ }
+ }
+ for (unsigned int i=0;i<partition.size();i++){
+ std::set<int> tags;
+ for (int j=0;j<partition[i]->tri.size();j++){
+ for (int k=0;k<3;k++){
+ MEdge e = partition[i]->tri[j]->getEdge(k);
+ std::map<MEdge, GEdge*> :: iterator it = _dictionnary.find(e);
+ if (it != _dictionnary.end()){
+ tags.insert(it->second->tag());
+ }
+ }
+ }
+ }
+}
+
+*/
void discreteFace::split(triangulation* trian,std::vector<triangulation*> &partition,int nPartitions)
{
#if defined(HAVE_METIS)
@@ -353,7 +452,6 @@ void discreteFace::split(triangulation* trian,std::vector<triangulation*> &parti
partition[i] = new triangulation(elem[i],this);
-
for(int i=0; i<nPartitions; i++){// setting GEdge's
//"update" old GEdge's
diff --git a/Mesh/Generator.cpp b/Mesh/Generator.cpp
index 8daad2869e..b81bf497ce 100644
--- a/Mesh/Generator.cpp
+++ b/Mesh/Generator.cpp
@@ -35,6 +35,7 @@
#include "Options.h"
#include "simple3D.h"
#include "yamakawa.h"
+#include "meshGRegionRelocateVertex.h"
#include "pointInsertion.h"
@@ -785,7 +786,7 @@ static void Mesh3D(GModel *m)
}
double a = Cpu();
- CTX::instance()->mesh.recombine3DLevel = 2;
+ // CTX::instance()->mesh.recombine3DLevel = 2;
if (CTX::instance()->mesh.recombine3DLevel >= 0){
Recombinator rec;
@@ -796,10 +797,12 @@ static void Mesh3D(GModel *m)
sup.execute(gr);
}
PostOp post;
- post.execute(gr,CTX::instance()->mesh.recombine3DLevel,3);
+ printf("-----------> %d %d\n",CTX::instance()->mesh.recombine3DLevel,CTX::instance()->mesh.recombine3DConformity);
+ post.execute(gr,CTX::instance()->mesh.recombine3DLevel,CTX::instance()->mesh.recombine3DConformity);
// CTX::instance()->mesh.recombine3DConformity);
// 0: no pyramid, 1: single-step, 2: two-steps (conforming),
// true: fill non-conformities with trihedra
+ RelocateVertices(gr);
// while(LaplaceSmoothing (gr)){
// }
nb_elements_recombination += post.get_nb_elements();
diff --git a/Mesh/meshGRegionRelocateVertex.cpp b/Mesh/meshGRegionRelocateVertex.cpp
index ac2c899880..3a4ac88484 100644
--- a/Mesh/meshGRegionRelocateVertex.cpp
+++ b/Mesh/meshGRegionRelocateVertex.cpp
@@ -9,6 +9,7 @@
#include "MHexahedron.h"
#include "Context.h"
#include "meshGFaceOptimize.h"
+#include "meshGRegionRelocateVertex.h"
static double objective_function (double xi, MVertex *ver,
double xTarget, double yTarget, double zTarget,
@@ -238,18 +239,22 @@ void RelocateVertices (GFace* gf, int niter, double tol) {
}
+void RelocateVertices (GRegion* region, double tol) {
+ v2t_cont adj;
+ buildVertexToElement(region->tetrahedra, adj);
+ buildVertexToElement(region->pyramids, adj);
+ buildVertexToElement(region->prisms, adj);
+ buildVertexToElement(region->hexahedra, adj);
+ v2t_cont::iterator it = adj.begin();
+ while (it != adj.end()){
+ _relocateVertex( it->first, it->second, tol);
+ ++it;
+ }
+}
+
void RelocateVertices (std::vector<GRegion*> ®ions, double tol) {
for(unsigned int k = 0; k < regions.size(); k++){
- v2t_cont adj;
- buildVertexToElement(regions[k]->tetrahedra, adj);
- buildVertexToElement(regions[k]->pyramids, adj);
- buildVertexToElement(regions[k]->prisms, adj);
- buildVertexToElement(regions[k]->hexahedra, adj);
- v2t_cont::iterator it = adj.begin();
- while (it != adj.end()){
- _relocateVertex( it->first, it->second, tol);
- ++it;
- }
+ RelocateVertices (regions[k], tol);
}
}
diff --git a/Mesh/meshGRegionRelocateVertex.h b/Mesh/meshGRegionRelocateVertex.h
index 2b94d81867..56c8ac45f7 100644
--- a/Mesh/meshGRegionRelocateVertex.h
+++ b/Mesh/meshGRegionRelocateVertex.h
@@ -3,6 +3,7 @@
#include <vector>
class GRegion;
class GFace;
+void RelocateVertices (GRegion* region, double tol = 1.e-2);
void RelocateVertices (std::vector<GRegion*> ®ions, double tol = 1.e-2);
void RelocateVertices (GFace*, int niter, double tol = 1.e-3);
#endif
diff --git a/Mesh/yamakawa.cpp b/Mesh/yamakawa.cpp
index e3bf3407fe..6eabdaf276 100644
--- a/Mesh/yamakawa.cpp
+++ b/Mesh/yamakawa.cpp
@@ -5014,14 +5014,14 @@ void PostOp::execute(GRegion* gr,int level, int conformity){
rearrange(gr);
}
-if(conformity == 2 || conformity == 3){
- init_markings(gr);
- build_vertex_to_tetrahedra(gr);
- build_vertex_to_pyramids(gr);
- pyramids2(gr);
- rearrange(gr);
-}
-
+ if(conformity == 2 || conformity == 3){
+ init_markings(gr);
+ build_vertex_to_tetrahedra(gr);
+ build_vertex_to_pyramids(gr);
+ pyramids2(gr);
+ rearrange(gr);
+ }
+
if (conformity == 3 || conformity == 4){
init_markings_hex(gr);
build_vertex_to_tetrahedra(gr);
@@ -5694,9 +5694,6 @@ void PostOp::trihedra(GRegion* gr){
std::vector<MPyramid*> opt2;
std::map<MElement*,bool>::iterator it;
- hexahedra.clear();
- prisms.clear();
-
for(i=0;i<gr->getNumMeshElements();i++){
element = gr->getMeshElement(i);
if(eight(element)){
diff --git a/Parser/Gmsh.y b/Parser/Gmsh.y
index 0a1ad9c387..3abe92aefc 100644
--- a/Parser/Gmsh.y
+++ b/Parser/Gmsh.y
@@ -4627,6 +4627,7 @@ Constraints :
"slave (%d) edges", List_Nbr($10), List_Nbr($5));
}
else{
+ Msg::Warning("This syntax is deprecated\n Use \"Periodic Surface { List } = { List } Translate { dx , dy , dz } ;\"");
int j_master = (int)$8;
int j_slave = (int)$3;
std::map<int,int> edgeCounterParts;
--
GitLab