diff --git a/Mesh/meshGFaceRecombine.cpp b/Mesh/meshGFaceRecombine.cpp
index b1b38c7347321f0e76c6f0814611e2c25a1f2600..e344ac5d3a77e4f7083f3ed8f5e496bb760e4899 100644
--- a/Mesh/meshGFaceRecombine.cpp
+++ b/Mesh/meshGFaceRecombine.cpp
@@ -1,4 +1,4 @@
-// Gmsh - Copyright (C) 1997-2011 C. Geuzaine, J.-F. Remacle
+// Gmsh - Copyright (C) 1997-2012 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to <gmsh@geuz.org>.
@@ -7,28 +7,21 @@
// Amaury Johnen (a.johnen@ulg.ac.be)
//
-#define REC2D_EDGE_BASE 1
-#define REC2D_EDGE_QUAD 1
-#define REC2D_ALIGNMENT .5
-#define REC2D_WAIT_TIME .01
-#define REC2D_NUM_ACTIO 550
-
#include "meshGFaceRecombine.h"
-#include "MElement.h"
+//#include "MElement.h"
#include "MTriangle.h"
#include "MQuadrangle.h"
-#include "meshGFaceOptimize.h"
-
-#include "drawContext.h"
-#include "FlGui.h"
-#include "Context.h"
-#include "OS.h"
+//#include "meshGFaceOptimize.h"
-// #define REC2D_SMOOTH
- #define REC2D_DRAW
+#ifdef REC2D_DRAW
+ #include "drawContext.h"
+ #include "FlGui.h"
+ #include "Context.h"
+ #include "OS.h"
+#endif
Recombine2D *Recombine2D::_current = NULL;
-bool Recombine2D::bo = 0;
+Rec2DData *Rec2DData::_current = NULL;
double **Rec2DVertex::_qualVSnum = NULL;
double **Rec2DVertex::_gains = NULL;
@@ -38,661 +31,455 @@ int otherParity(int a) {
return a + 1;
}
-/*
- TODO :
- - map -> myset<pair> (mapCornerVert, mapVert, mapEdge)
- - set -> myset
- (_vertices, _edges, _elements) -> set
- (_vertices, _edges, _elements) -> set
-*/
-
-void Recombine2D::associateParity(int pOld, int pNew)
-{
- std::set<Rec2DVertex*>::iterator it = _current->_parities[pOld].begin();
- for (; it != _current->_parities[pOld].end(); ++it)
- (*it)->setParity(pNew);
- _current->_parities.erase(pOld);
- pOld = otherParity(pOld);
- pNew = otherParity(pNew);
- it = _current->_parities[pOld].begin();
- for (; it != _current->_parities[pOld].end(); ++it)
- (*it)->setParity(pNew);
- _current->_parities.erase(pOld);
-}
-
/** Recombine2D **/
/*******************/
Recombine2D::Recombine2D(GFace *gf) : _gf(gf)
{
if (Recombine2D::_current != NULL) {
- Msg::Warning("An instance of recombination is already in execution");
+ Msg::Warning("[Recombine2D] An instance already in execution");
return;
}
Recombine2D::_current = this;
- backgroundMesh::set(gf);
+ backgroundMesh::set(_gf);
_bgm = backgroundMesh::current();
- _numChange = 0;
- _numEdge = _numVert = 0;
- _valEdge = _valVert = .0;
+ _data = new Rec2DData(gf->triangles.size(), gf->quadrangles.size());
Rec2DVertex::initStaticTable();
+ _numChange = 0;
- _tri = _gf->triangles;
+#ifdef REC2D_DRAW
+ //_data->_tri = _gf->triangles;
+ //_data->_quad = _gf->quadrangles;
+#endif
// Be able to compute geometrical angle at corners
- std::map<MVertex*, std::set<GEdge*> > mapCornerVert; // TOBE myset
+ std::map<MVertex*, AngleData> mapCornerVert;
{
std::list<GEdge*> listge = gf->edges();
std::list<GEdge*>::iterator itge = listge.begin();
for (; itge != listge.end(); ++itge) {
mapCornerVert[(*itge)->getBeginVertex()->getMeshElement(0)->getVertex(0)]
- .insert(*itge);
+ ._gEdges.push_back(*itge);
mapCornerVert[(*itge)->getEndVertex()->getMeshElement(0)->getVertex(0)]
- .insert(*itge);
- }
- }
-
- // Find all Vertices and edges and create the 'Rec2DElement'
- std::map<MVertex*, std::vector<MTriangle*> > mapVert; // TOBE myset
- std::map<MEdge, std::vector<MTriangle*>, Less_Edge> mapEdge; // TOBE myset
- for (unsigned int i = 0; i < gf->triangles.size(); ++i) {
- MTriangle *t = gf->triangles[i];
- Rec2DElement *el = new Rec2DElement(t);
- _elements[t] = el; // ISITTOKEEP keep only list of relem, map local
- for (int j = 0; j < 3; ++j) {
- mapVert[t->getVertex(j)].push_back(t);
- mapEdge[t->getEdge(j)].push_back(t);
+ ._gEdges.push_back(*itge);
}
}
- // Create the 'Rec2DVertex' and store iterator to vertices which have degree 4
- std::list<std::map<MVertex*, std::vector<MTriangle*> >::iterator> fourTri; // TOBE vector
+ // Create the 'Rec2DEdge', the 'Rec2DVertex' and the 'Rec2DElement'
{
- mapofVertices::iterator itV2N = _vertices.begin();
- std::map<MVertex*, std::vector<MTriangle*> >::iterator itvert = mapVert.begin();
- for (; itvert != mapVert.end(); ++itvert) {
- MVertex *v = itvert->first;
- int onWhat = 1;
- if (mapCornerVert.find(v) != mapCornerVert.end())
- onWhat = -1;
- Rec2DVertex *rv = new Rec2DVertex(v, itvert->second, onWhat, mapCornerVert);
- itV2N = _vertices.insert(itV2N, std::make_pair(v,rv)); // ISITTOKEEP keep only list of rvert, map local
- if (itvert->second.size() == 4)
- fourTri.push_back(itvert);
+ std::map<MVertex*, Rec2DVertex*> mapVert;
+ std::map<MVertex*, Rec2DVertex*>::iterator itV;
+ std::map<MVertex*, AngleData>::iterator itCorner;
+
+ // triangles
+ for (unsigned int i = 0; i < gf->triangles.size(); ++i) {
+ MTriangle *t = gf->triangles[i];
+ Rec2DElement *rel = new Rec2DElement(t);
+ Rec2DVertex *rv[3];
+ for (int j = 0; j < 3; ++j) {
+ MVertex *v = t->getVertex(j);
+ if ( (itCorner = mapCornerVert.find(v)) != mapCornerVert.end() ) {
+ if (!itCorner->second._rv)
+ itCorner->second._rv = new Rec2DVertex(v, false);
+ rv[j] = itCorner->second._rv;
+ itCorner->second._mElements.push_back((MElement*)t);
+ }
+ else if ( (itV = mapVert.find(v)) != mapVert.end() ) {
+ rv[j] = itV->second;
+ }
+ else {
+ rv[j] = new Rec2DVertex(v);
+ mapVert[v] = rv[j];
+ }
+ rv[j]->add(rel); //up data
+ }
+ for (int j = 0; j < 3; ++j) {
+ Rec2DEdge *re;
+ if ( (re = Rec2DVertex::getCommonEdge(rv[j], rv[(j+1)%3])) == NULL)
+ re = new Rec2DEdge(rv[j], rv[(j+1)%3]);
+ rel->add(re); //up data
+ }
}
- }
-
- // Create the 'Rec2DEdge' and store boundary edges
- std::map<MVertex*, std::list<MEdge> > boundV2E; // TOBE myset
- {
- mapofEdges::iterator itE2E = _edges.begin();
- std::map<MEdge, std::vector<MTriangle*> >::iterator itedge = mapEdge.begin();
- for (; itedge != mapEdge.end(); ++itedge) {
- MEdge e = itedge->first;
- Rec2DEdge *re = new Rec2DEdge(e, _vertices, itedge->second);
- itE2E = _edges.insert(itE2E, std::make_pair(e,re));
- _numEdge += REC2D_EDGE_BASE;
- _valEdge += REC2D_EDGE_BASE * re->getQual();
- if (itedge->second.size() == 1) {
- _numEdge += REC2D_EDGE_QUAD;
- _valEdge += REC2D_EDGE_QUAD * re->getQual();
- _vertices[e.getVertex(0)]->setOnBoundary();
- _vertices[e.getVertex(1)]->setOnBoundary();
- boundV2E[e.getVertex(0)].push_back(e);
- boundV2E[e.getVertex(1)].push_back(e);
+
+ // quadrangles
+ for (unsigned int i = 0; i < gf->quadrangles.size(); ++i) {
+ MQuadrangle *q = gf->quadrangles[i];
+ Rec2DElement *rel = new Rec2DElement(q);
+ Rec2DVertex *rv[4];
+ for (int j = 0; j < 4; ++j) {
+ MVertex *v = q->getVertex(j);
+ if ( (itCorner = mapCornerVert.find(v)) != mapCornerVert.end() ) {
+ if (!itCorner->second._rv)
+ itCorner->second._rv = new Rec2DVertex(v, false);
+ rv[j] = itCorner->second._rv;
+ itCorner->second._mElements.push_back((MElement*)q);
+ }
+ else if ( (itV = mapVert.find(v)) == mapVert.end() ) {
+ rv[j] = itV->second;
+ }
+ else {
+ rv[j] = new Rec2DVertex(v);
+ mapVert[v] = rv[j];
+ }
+ rv[j]->add(rel);
+ }
+ for (int j = 0; j < 4; ++j) {
+ Rec2DEdge *re;
+ if ( (re = Rec2DVertex::getCommonEdge(rv[i], rv[(i+1)%3])) == NULL) {
+ re = new Rec2DEdge(rv[i], rv[(i+1)%3]);
+ rv[i]->add(re);
+ rv[(i+1)%3]->add(re);
+ }
+ rel->add(re);
}
}
}
-
- // We know now on what are vertices, compute reward
+ // update corner
{
- mapofVertices::iterator itV2N = _vertices.begin();
- for (; itV2N != _vertices.end(); ++itV2N) {
- ++_numVert;
- _valVert += itV2N->second->getQual();
+ std::map<MVertex*, AngleData>::iterator it = mapCornerVert.begin();
+ for (; it != mapCornerVert.end(); ++it) {
+ double angle = _geomAngle(it->first,
+ it->second._gEdges,
+ it->second._mElements);
+ //new Rec2DVertex(it->second._rv, angle);
}
}
-
- // Be dealing with "parity" on boundaries
+ mapCornerVert.clear();
+ // update neighbouring, create the 'Rec2DTwoTri2Quad'
{
- std::map<MVertex*, std::list<MEdge> >::iterator itBound = boundV2E.begin();
- int actualParity = 0;
- bool p = 0;
- while (itBound != boundV2E.end()) {
- MVertex *firstV = itBound->first, *actualV = firstV;
- MEdge actualEdge = *itBound->second.begin();
- Rec2DVertex *rv = _vertices[firstV];
- rv->setParity(actualParity+p);
- _parities[actualParity+p].insert(rv);
- p = !p;
- if (actualEdge.getVertex(0) == actualV)
- actualV = actualEdge.getVertex(1);
- else
- actualV = actualEdge.getVertex(0);
- while (actualV != firstV) {
- rv = _vertices[actualV];
- rv->setParity(actualParity+p);
- _parities[actualParity+p].insert(rv);
- p = !p;
- itBound = boundV2E.find(actualV);
- std::list<MEdge>::iterator it = itBound->second.begin();
- if (*it == actualEdge)
- actualEdge = *(++it);
- else
- actualEdge = *it;
- if (actualEdge.getVertex(0) == actualV)
- actualV = actualEdge.getVertex(1);
- else
- actualV = actualEdge.getVertex(0);
- boundV2E.erase(itBound);
+ Rec2DData::iter_re it = Rec2DData::firstEdge();
+ for (; it != Rec2DData::lastEdge(); ++it) {
+ Rec2DVertex *rv0 = (*it)->getVertex(0), *rv1 = (*it)->getVertex(1);
+ if ((*it)->isOnBoundary()) {
+ rv0->setOnBoundary();
+ rv1->setOnBoundary();
+ }
+ else {
+ std::vector<Rec2DElement*> elem;
+ Rec2DVertex::getCommonElements(rv0, rv1, elem);
+ if (elem.size() != 2)
+ Msg::Error("[Recombine2D] %d elements instead of 2 for neighbour link",
+ elem.size() );
+ else {
+ elem[0]->addNeighbour(*it, elem[1]);
+ elem[1]->addNeighbour(*it, elem[0]);
+ new Rec2DTwoTri2Quad(elem[0], elem[1]);
+ }
}
- if (_vertices[actualV]->getParity() != actualParity+p)
- Msg::Error("[Recombine2D] Not even number of mesh edges on boundary !");
-
- boundV2E.erase(boundV2E.begin());
- itBound = boundV2E.begin();
- actualParity += 2;
}
}
- // Create the 'Rec2DTwoTri2Quad' and 'Rec2DCollapseTri'
+ // set parity on boundary, create 'Rec2DFourTri2Quad'
{
- mapofEdges::iterator itE2E = _edges.begin();
- std::map<MEdge, std::vector<MTriangle*> >::iterator ite = mapEdge.begin();
- for (; ite != mapEdge.end(); ++ite)
- if (ite->second.size() == 2) {
- Rec2DTwoTri2Quad *q = new Rec2DTwoTri2Quad(_edges[ite->first], ite->second);
- //Rec2DCollapseTri *c = new Rec2DCollapseTri(_edges[ite->first], ite->second);
- _actions.push_back(q);
- //_actions.insert(c);
- std::vector<MTriangle*>::iterator it = ite->second.begin();
- _mea[*it].insert(q);
- //_mea[*it].insert(c);
- ++it;
- _mea[*it].insert(q);
- //_mea[*it].insert(c);
+ Rec2DData::iter_rv it = Rec2DData::firstVertex();
+ for (; it != Rec2DData::lastVertex(); ++it) {
+ Rec2DVertex *rv = *it;
+ if (rv->getOnBoundary()) {
+ if (rv->getParity() == -1) {
+ int base = Rec2DData::getNewParity();
+ rv->setBoundaryParity(base, base+1);
+ }
+ }
+ else if (rv->getNumElements() == 4) {
+ ;
+ }
}
}
- // Create the 'Rec2DFourTri2Quad' and 'Rec2DCollapseTri'
- //{
- // std::list<std::map<MVertex*, std::vector<MTriangle*> >::iterator>::iterator it4;
- // for (it4 = fourTri.begin(); it4 != fourTri.end(); it4++) {
- // Rec2DVertex *rv = _vertices[(*it4)->first];
- // if (!rv->getIsOnBoundary()) {
- // Rec2DFourTri2Quad *q = new Rec2DFourTri2Quad(rv, (*it4)->second);
- // //Rec2DCollapseTri *c = new Rec2DCollapseTri(rv, (*it4)->second);
- // _actions.push_back(q);
- // //_actions.insert(c);
- // std::vector<MTriangle*>::iterator it = (*it4)->second.begin();
- // for (; it != (*it4)->second.end(); ++it) {
- // _mea[*it].insert(q);
- // //_mea[*it].insert(c);
- // }
- // }
- // }
- //}
-
- Msg::Info("State");
- Msg::Info("-----");
- Msg::Info("numEdge %d (%d), valEdge %g", _numEdge, _edges.size(), _valEdge);
- Msg::Info("numVert %d (%d), valVert %g", _numVert, _vertices.size(), _valVert);
- Msg::Info("global Value %g", Recombine2D::getGlobalValue(0,0,0,0));
- Msg::Info("num action %d", _actions.size());
- std::map<int, std::set<Rec2DVertex*> >::iterator it = _parities.begin();
- for (; it != _parities.end(); ++it) {
- Msg::Info("par %d, #%d", it->first, it->second.size());
- }
- mapofElementActions::iterator it2 = _mea.begin();
- int a = 0;
- for (; it2 != _mea.end(); ++it2) {
- a += it2->second.size();
- }
- if (a > 0)
- Msg::Info("%d elements in mapof elem to action, with average action %g", _mea.size(), (double)a/_mea.size());
- Msg::Info(" ");
+ _data->printState();
}
Recombine2D::~Recombine2D()
{
- Recombine2D::_current = NULL;
- //delete vertices, edges, elements, actions;
+ delete _data;
+ if (Recombine2D::_current == this)
+ Recombine2D::_current = NULL;
}
bool Recombine2D::recombine()
{
Rec2DAction *nextAction;
- int i = -1;
- while (_actions.size() > 0 && ++i < REC2D_NUM_ACTIO) {
- double t = Cpu();
- nextAction = *std::max_element(_actions.begin(), _actions.end(),
- lessRec2DAction() );
+ while (nextAction = Rec2DData::getBestAction()) {
+ FlGui::instance()->check();
+ double time = Cpu();
+#ifdef REC2D_DRAW
nextAction->color(0, 0, 200);
-
CTX::instance()->mesh.changed = (ENT_ALL);
- FlGui::instance()->check();
drawContext::global()->draw();
- //if (nextAction->getReward() > 0) {
- if (_isAllQuad(nextAction))
- _choose(nextAction);
- else {
- nextAction->color(190, 0, 0);
- _removeReference(nextAction);
- }
- /*}
- else {
- nextAction->color(0, 100, 0);
- _removeReference(nextAction);
- }*/
+#endif
+
+ ++_numChange;
+ nextAction->apply();
- _gf->triangles = _tri;
- _gf->quadrangles = _quad;
#ifdef REC2D_DRAW
- //Msg::Info("%d", REC2D_NUM_ACTIO-i);
- //Msg::GetAnswer("Continue ?", 1, "no", "yes");
+ _gf->triangles = _data->_tri;
+ _gf->quadrangles = _data->_quad;
CTX::instance()->mesh.changed = (ENT_ALL);
- FlGui::instance()->check();
- if (i%5 == 0)
drawContext::global()->draw();
- while (Cpu()-t < REC2D_WAIT_TIME)
+ while (Cpu()-time < REC2D_WAIT_TIME)
FlGui::instance()->check();
#endif
-#ifdef REC2D_SMOOTH
- laplaceSmoothing(_gf,1);
-#endif
-
- std::map<int, std::set<Rec2DVertex*> >::iterator it = _parities.begin();
- for (; it != _parities.end(); ++it) {
- //Msg::Info("par %d, #%d", it->first, it->second.size());
- }
+
+ delete nextAction;
}
-//
-// Msg::Info("State");
-// Msg::Info("-----");
-// Msg::Info("numEdge %d (%d), valEdge %g", _numEdge, _edges.size(), _valEdge);
-// Msg::Info("numVert %d (%d), valVert %g", _numVert, _vertices.size(), _valVert);
-// Msg::Info("global Value %g", Recombine2D::getGlobalValue(0,0,0,0));
-// Msg::Info("num action %d", _actions.size());
-// std::map<int, std::set<Rec2DVertex*> >::iterator it = _parities.begin();
-// for (; it != _parities.end(); ++it) {
-// Msg::Info("par %d, #%d", it->first, it->second.size());
-// }
-// mapofElementActions::iterator it2 = _mea.begin();
-// int a = 0;
-// for (; it2 != _mea.end(); ++it2) {
-// a += it2->second.size();
-// }
-// if (a > 0)
-// Msg::Info("%d elements in mapof elem to action, with average action %g", _mea.size(), (double)a/_mea.size());
-// Msg::Info(" ");
-//
return 1;
}
-bool Recombine2D::_isAllQuad(Rec2DAction *action)
+double Recombine2D::_geomAngle(MVertex *v,
+ std::vector<GEdge*> &gEdge,
+ std::vector<MElement*> &elem) //*
{
- //Msg::Info(" ---new isAllQuad---");
-
- int p[4];
- action->getParities(p);
-
- if (p[0] < 0 && p[1] < 0 && p[2] < 0 && p[3] < 0) {
- Msg::Info("is isolated");
- return true;
- }
-
- if (action->isObsolete()) {
- Msg::Error("[Recombine2D] obsolete action (allQuad)");
- return false;
+ if (gEdge.size() != 2) {
+ Msg::Error("[Recombine2D] Wrong number of edge : %d, returning pi/2",
+ gEdge.size() );
+ return M_PI / 2.;
}
+ static const double prec = 100.;
- std::set<Rec2DVertex*> setRv;
- std::set<Rec2DElement*> neighbours, checked;
- for (int i = 0; i < Recombine2D::getNumVertAllQuad(); ++i) {
- if (Recombine2D::getVertAllQuad(i)->hasTriangle())
- Recombine2D::getVertAllQuad(i)->getElements(neighbours);
- else
- setRv.insert(Recombine2D::getVertAllQuad(i));
- }
- std::set<Rec2DVertex*>::iterator itit = setRv.begin();
- for (; itit != setRv.end(); ++itit)
- Recombine2D::removeVertAllQuad(*itit);
- setRv.clear();
+ SVector3 vectv = SVector3(v->x(), v->y(), v->z());
+ SVector3 firstDer0, firstDer1;
- int min = 100, index;
- for (int i = 0; i < 4; ++i) {
- if (p[i] > -1 && min > p[i]) {
- min = p[i];
- index = i;
+ for (unsigned int k = 0; k < 2; ++k) {
+ GEdge *ge = gEdge[k];
+ SVector3 vectlb = ge->position(ge->getLowerBound());
+ SVector3 vectub = ge->position(ge->getUpperBound());
+ vectlb -= vectv;
+ vectub -= vectv;
+ double normlb, normub;
+ if ((normlb = norm(vectlb)) > prec * (normub = norm(vectub)))
+ firstDer1 = -1 * ge->firstDer(ge->getUpperBound());
+ else if (normub > prec * normlb)
+ firstDer1 = ge->firstDer(ge->getLowerBound());
+ else {
+ Msg::Error("[Recombine2D] Bad precision, returning pi/2");
+ return M_PI / 2.;
}
+ if (k == 0)
+ firstDer0 = firstDer1;
}
- //Msg::Info("Passsed through there, [%d %d %d %d] -> min %d", p[0], p[1], p[2], p[3], min);
- std::map<Rec2DVertex*, int> newParVert;
- std::map<Rec2DVertex*, int>::iterator it2;
- std::map<int, int> parAssoc; // first parity value become second parity value
- std::map<int, int>::iterator itassoc;
+ firstDer0 = firstDer0 * (1 / norm(firstDer0));
+ firstDer1 = firstDer1 * (1 / norm(firstDer1));
- for (int i = 0; i < 4; i += 2) {
- int par;
- if ((index/2) * 2 == i)
- par = min;
- else
- par = otherParity(min);
- for (int j = 0; j < 2; ++j) {
- if (p[i+j] == -1) {
- newParVert.insert(std::make_pair(action->getRVertex(i+j), par));
- action->getRVertex(i+j)->getElements(neighbours);
- }
- else if (p[i+j] != par) {
- if (parAssoc.find(p[i+j]) == parAssoc.end()) {
- parAssoc[p[i+j]] = par;
- parAssoc[otherParity(p[i+j])] = otherParity(par);
- _getElemToLook(p[i+j], neighbours);
- _getElemToLook(otherParity(p[i+j]), neighbours);
- }
- }
- }
- }
+ double angle1 = acos(dot(firstDer0, firstDer1));
+ double angle2 = 2. * M_PI - angle1;
- while (neighbours.size() > 0) {
- //Msg::Info("num vert %d, num neigh %d", newParVert.size(), neighbours.size());
-
- std::set<Rec2DElement*>::iterator itel = neighbours.begin();
- Rec2DElement *element = *itel;
- if (element->isQuad()) {
- neighbours.erase(itel);
- continue;
- }
-
- Rec2DVertex *v[3];
- element->getVertices(v);
- int p[3];
- for (int i = 0; i < 3; ++i) {
- if ((it2 = newParVert.find(v[i])) != newParVert.end()) {
- p[i] = it2->second;
- while ((itassoc = parAssoc.find(p[i])) != parAssoc.end()) {
- p[i] = itassoc->second;
- }
- it2->second = p[i];
- }
- else {
- p[i] = v[i]->getParity();
- if ((itassoc = parAssoc.find(p[i])) != parAssoc.end()) {
- do
- p[i] = itassoc->second;
- while ((itassoc = parAssoc.find(p[i])) != parAssoc.end());
- newParVert[v[i]] = p[i];
- }
- }
- }
- //Msg::Info("tri %d [%d %d %d]", (*itel)->getNum(), p[0], p[1], p[2]);
-
- bool hasIdentical = false;
- for (int i = 0; i < 3; ++i) {
- if (p[i] > -1 && p[i] == p[(i+1)%3]) hasIdentical = true;
- }
- if (!hasIdentical) {
- neighbours.erase(itel);
- continue;
- }
- if (p[0] == p[1] && p[0] == p[2]) {
- Msg::Info("3 identical par");
- return false;
- }
- //Msg::Info("has identical");
-
- bool hasAction = false;
- std::map<Rec2DVertex*, std::vector<int> > suggestions;
- for (int i = 0; i < element->getNumActions(); ++i) {
- if (element->getAction(i)->whatImply(v, p, suggestions))
- hasAction = true;
- }
- //Msg::Info("suggest %d", suggestions.size());
-
- if (!hasAction) {
- Msg::Info("No action %d", (*itel)->getNum());
- return false;
- }
-
- std::map<Rec2DVertex*, std::vector<int> >::iterator itsug;
- itsug = suggestions.begin();
- for (; itsug != suggestions.end(); ++itsug) {
- int par = itsug->second.front();
- bool onlyOnePar = true;
- for (unsigned int i = 1; i < itsug->second.size(); ++i) {
- if (itsug->second[i] != par) {
- onlyOnePar = false;
- break;
- }
- }
- if (onlyOnePar) {
- if ((itassoc = parAssoc.find(par)) != parAssoc.end()) {
- do
- par = itassoc->second;
- while ((itassoc = parAssoc.find(par)) != parAssoc.end());
- }
-
- int oldPar;
- Rec2DVertex *v = itsug->first;
-
- if ((it2 = newParVert.find(v)) != newParVert.end()) {
- oldPar = it2->second;
- while ((itassoc = parAssoc.find(oldPar)) != parAssoc.end()) {
- oldPar = itassoc->second;
- }
- it2->second = oldPar;
- }
- else {
- oldPar = v->getParity();
- if ((itassoc = parAssoc.find(oldPar)) != parAssoc.end()) {
- do
- oldPar = itassoc->second;
- while ((itassoc = parAssoc.find(oldPar)) != parAssoc.end());
- newParVert[v] = oldPar;
- }
- }
-
- //Msg::Info("a %d, %d", par, oldPar);
- if (oldPar == -1) {
- //Msg::Info("b");
- newParVert[v] = par;
- v->getElements(neighbours);
- }
- else if ((par/2)*2 != (oldPar/2)*2) {
- //Msg::Info("c");
- if (oldPar < par) {
- int a = oldPar;
- oldPar = par;
- par = a;
- }
- parAssoc[oldPar] = par;
- parAssoc[otherParity(oldPar)] = otherParity(par);
- _getElemToLook(oldPar, neighbours);
- _getElemToLook(otherParity(oldPar), neighbours);
- }
- else if (par%2 != oldPar%2) {
- //Msg::Info("d");
- Msg::Info("not all quad");
- return false;
- }
- }
+ double angleMesh = .0;
+ std::vector<MElement*>::iterator it2 = elem.begin();
+ for (; it2 != elem.end(); ++it2) {
+ MElement *el = *it2;
+ int k = 0, numV = el->getNumVertices();
+ while (el->getVertex(k) != v && k < numV) ++k;
+ if (k == numV) {
+ Msg::Error("[Recombine2D] Wrong element, returning pi/2");
+ return M_PI / 2.;
}
- neighbours.erase(itel);
+ angleMesh += angle3Vertices(el->getVertex((k+1) % numV), v,
+ el->getVertex((k+2) % numV) );
}
- Msg::Info("all quad");
- return true;
+
+ if (angleMesh < M_PI)
+ return angle1;
+ return angle2;
}
-void Recombine2D::_choose(Rec2DAction *action)
-{
- action->apply(); //to check
- _removeIncompatible(action); //to check
- ++_numChange; //to check
-}
-void Recombine2D::_removeIncompatible(Rec2DAction *action) //to check
+/** Rec2DData **/
+/*****************/
+Rec2DData::Rec2DData(int numTri, int numQuad)
{
- std::set<Rec2DElement*> touchedTri;
- action->getRTriangles(touchedTri);
- std::set<Rec2DElement*>::iterator itTri;
-
- std::set<Rec2DAction*> incomp;
- {
- itTri = touchedTri.begin();
- for (; itTri != touchedTri.end(); ++itTri)
- (*itTri)->getActions(incomp);
+ if (Rec2DData::_current != NULL) {
+ Msg::Error("[Rec2DData] An instance in execution");
+ return;
}
+ Rec2DData::_current = this;
+ _numEdge = _numVert = 0;
+ _valEdge = _valVert = .0;
- //Msg::Warning("size %d", incomp.size());
-
- std::set<Rec2DAction*>::iterator itAction = incomp.begin();
- for (; itAction != incomp.end(); ++itAction) {
- delete *itAction;
- }
+ //_elements.reserve((int) (numTri + numQuad) * 1.1);
+ //_edges.reserve((int) (numTri * 1.8 + numQuad * 2.4));
+ //_vertices.reserve((int) (numTri * .6 + numQuad * 1.2));
}
-void Recombine2D::_removeReference(Rec2DAction *action)
-{
- delete action;
-}
-
-//void Recombine2D::_getNeighbours(Rec2DAction *ra, std::set<Rec2DAction*> &actions)
-//{
-// std::set<Rec2DElement*> touchedTri, neighbourTri;
-// ra->getRTriangles(touchedTri);
-// std::set<Rec2DElement*>::iterator itTri;
-//
-// {
-// itTri = touchedTri.begin();
-// for (; itTri != touchedTri.end(); ++itTri)
-// (*itTri)->getNeighbours(neighbourTri);
-// //itTri = touchedTri.begin();
-// //for (; itTri != touchedTri.end(); ++itTri)
-// // neighbourTri.erase(*itTri);
-// }
-//
-// std::set<Rec2DAction*> incomp, neighbours;
-// {
-// itTri = touchedTri.begin();
-// for (; itTri != touchedTri.end(); ++itTri)
-// (*itTri)->getActions(neighbours);
-// itTri = neighbourTri.begin();
-// for (; itTri != neighbourTri.end(); ++itTri)
-// (*itTri)->getActions(neighbours);
-// }
-//
-// std::set<Rec2DAction*>::iterator itAction = incomp.begin();
-// for (; itAction != incomp.end(); ++itAction) {
-// bool b = neighbours.erase(*itAction);
-// if (!b)
-// Msg::Error("[Recombine2D] didn't found incomp action");
-// }
-//
-// for (itAction = neighbours.begin(); itAction != neighbours.end(); ++itAction) {
-// actions.insert(*itAction);
-// }
-//}
-
-void Recombine2D::_getElemToLook(int par, std::set<Rec2DElement*> &elem)
-{
- std::set<Rec2DVertex*>::iterator it = _parities[par].begin();
- std::set<Rec2DVertex*>::iterator itend = _parities[par].end();
- for (; it != itend; ++it) {
- (*it)->getElements(elem);
- }
+Rec2DData::~Rec2DData()
+{
+ if (Rec2DData::_current == this)
+ Rec2DData::_current = NULL;
}
-
-double Recombine2D::getGlobalValue(int numEdge, double valEdge,
- int numVert, double valVert )
+#ifdef REC2D_DRAW
+void Rec2DData::add(Rec2DElement *rel)
{
- double a = (_current->_valVert + valVert) / (_current->_numVert + numVert);
- return (_current->_valEdge + valEdge) / (_current->_numEdge + numEdge) * a * a;
+ _current->_elements.insert(rel);
+
+ MTriangle *t = rel->getMTriangle();
+ if (t)
+ _current->_tri.push_back(t);
+ MQuadrangle *q = rel->getMQuadrangle();
+ if (q)
+ _current->_quad.push_back(q);
+
}
+#endif
-//Rec2DVertex* Recombine2D::getRVertex(MVertex *v)
-//{
-// mapofVertices::iterator it = _current->_vertices.find(v);
-// if (it != _current->_vertices.end())
-// return it->second;
-// Msg::Error("[Recombine2D::getRVertex] vertex not found");
-// return NULL;
-//}
-
-Rec2DEdge* Recombine2D::getREdge(MEdge e)
+void Rec2DData::remove(Rec2DEdge *re)
{
- mapofEdges::iterator it = _current->_edges.find(e);
- if (it != _current->_edges.end())
- return it->second;
- Msg::Error("[Recombine2D::getREdge] edge not found");
- return NULL;
+ /*bool b = false;
+ unsigned int i = 0;
+ while (i < _current->_edges.size()) {
+ if (_current->_edges[i] == re) {
+ _current->_edges[i] = _current->_edges.back();
+ _current->_edges.pop_back();
+ if (b)
+ Msg::Error("[Rec2DData] Two or more times same edge");
+ b = true;
+ }
+ else
+ ++i;
+ }*/
+ _current->_edges.erase(re);
}
-Rec2DElement* Recombine2D::getRElement(MElement *el)
+void Rec2DData::remove(Rec2DVertex *rv)
{
- mapofElements::iterator it = _current->_elements.find(el);
- if (it != _current->_elements.end())
- return it->second;
- Msg::Error("[Recombine2D::getRElement] element not found");
- return NULL;
+ _current->_vertices.erase(rv);
}
-void Recombine2D::remove(MEdge e)
+void Rec2DData::remove(Rec2DElement *rel)
{
- mapofEdges::iterator it = _current->_edges.find(e);
- if (it != _current->_edges.end())
- _current->_edges.erase(it);
- else
- Msg::Error("[Recombine2D::getREdge] edge not found");
+ /*bool b = false;
+ unsigned int i = 0;
+ while (i < _current->_elements.size()) {
+ if (_current->_elements[i] == rel) {
+ _current->_elements[i] = _current->_elements.back();
+ _current->_elements.pop_back();
+ if (b)
+ Msg::Error("[Rec2DData] Two or more times same element");
+ b = true;
+ }
+ else
+ ++i;
+ }*/
+ _current->_elements.erase(rel);
+#ifdef REC2D_DRAW
+ MTriangle *t = rel->getMTriangle();
+ if (t) {
+ for (unsigned int i = 0; i < _current->_tri.size(); ++i) {
+ if (_current->_tri[i] == t) {
+ _current->_tri[i] = _current->_tri.back();
+ _current->_tri.pop_back();
+ return;
+ }
+ }
+ Msg::Info("[Rec2DData] Didn't erased mtriangle :(");
+ }
+ MQuadrangle *q = rel->getMQuadrangle();
+ if (q) {
+ for (unsigned int i = 0; i < _current->_quad.size(); ++i) {
+ if (_current->_quad[i] == q) {
+ _current->_quad[i] = _current->_quad.back();
+ _current->_quad.pop_back();
+ return;
+ }
+ }
+ Msg::Info("[Rec2DData] Didn't erased mquadrangle :(");
+ }
+#endif
}
-void Recombine2D::remove(MVertex *v)
+void Rec2DData::remove(Rec2DAction *ra)
{
- mapofVertices::iterator it = _current->_vertices.find(v);
- if (it != _current->_vertices.end())
- _current->_vertices.erase(it);
- else
- Msg::Error("[Recombine2D::getRVertex] vertex not found");
+ std::list<Rec2DAction*>::iterator it = _current->_actions.begin();
+ while (it != _current->_actions.end()) {
+ if (*it == ra)
+ it = _current->_actions.erase(it);
+ else
+ ++it;
+ }
}
-void Recombine2D::remove(Rec2DElement *tri)
+void Rec2DData::printState()
{
- Rec2DVertex *v[3];
- tri->getVertices(v);
- for (int i = 0; i < 3; ++i) {
- v[i]->remove(tri);
+ Msg::Info("State");
+ Msg::Info("-----");
+ Msg::Info("numEdge %d (%d), valEdge %g => %g", _numEdge, _edges.size(), _valEdge, _valEdge/_numEdge);
+ Msg::Info("numVert %d (%d), valVert %g => %g", _numVert, _vertices.size(), _valVert, _valVert/_numVert);
+ Msg::Info("Element (%d)", _elements.size());
+ Msg::Info("global Value %g", Rec2DData::getGlobalValue());
+ Msg::Info("num action %d", _actions.size());
+ std::map<int, std::vector<Rec2DVertex*> >::iterator it = _parities.begin();
+ for (; it != _parities.end(); ++it) {
+ Msg::Info("par %d, #%d", it->first, it->second.size());
}
- MTriangle *t = tri->getMTriangle();
- std::vector<MTriangle*>::iterator it = _current->_tri.begin();
- for (; it != _current->_tri.end(); ++it) {
- if (*it == t) {
- _current->_tri.erase(it);
- return;
- }
+ Msg::Info(" ");
+ iter_re ite;
+ long double valEdge = .0;
+ for (ite = firstEdge(); ite != lastEdge(); ++ite) {
+ valEdge += (*ite)->getQual();
+ }
+ Msg::Info("valEdge : %g >< %g", valEdge, _valEdge);
+ iter_rv itv;
+ long double valVert = .0;
+ for (itv = firstVertex(); itv != lastVertex(); ++itv) {
+ valVert += (*itv)->getQual();
}
- Msg::Error("[Recombine2D] triangle was not there");
+ Msg::Info("valVert : %g >< %g", valVert, _valVert);
+ Msg::Info(" : %g", valVert);
+ Msg::Info(" : %g", _valVert);
+}
+
+int Rec2DData::getNewParity()
+{
+ if (!_current->_parities.size())
+ return 0;
+ std::map<int, std::vector<Rec2DVertex*> >::iterator it;
+ it = _current->_parities.end();
+ --it;
+ return (it->first/2)*2 + 2;
}
-void Recombine2D::remove(Rec2DAction *ra)
+void Rec2DData::removeParity(Rec2DVertex *rv, int p)
{
- setofRec2DAction::iterator itAction = _current->_actions.begin();
+ std::map<int, std::vector<Rec2DVertex*> >::iterator it;
+ if ( (it = _current->_parities.find(p)) == _current->_parities.end() ) {
+ Msg::Error("[Rec2DData] Don't have parity %d", p);
+ return;
+ }
bool b = false;
- while (itAction != _current->_actions.end()) {
- if (*itAction == ra) {
- _current->_actions.erase(itAction);
+ std::vector<Rec2DVertex*> *vect = &it->second;
+ unsigned int i = 0;
+ while (i < vect->size()) {
+ if (vect->at(i) == rv) {
+ vect->at(i) = vect->back();
+ vect->pop_back();
+ if (b)
+ Msg::Error("[Rec2DData] Two or more times same vertex");
b = true;
- break;
}
- ++itAction;
+ else
+ ++i;
}
- if (!b)
- Msg::Error("[Recombine2D] action to delete not found");
}
-void Recombine2D::add(MQuadrangle *quad)
+double Rec2DData::getGlobalValue()
+{
+ double a = (_current->_valVert) / (_current->_numVert);
+ return a * (_current->_valEdge) / (_current->_numEdge);
+}
+
+double Rec2DData::getGlobalValue(int numEdge, double valEdge,
+ int numVert, double valVert)
+{
+ double a = (_current->_valVert + valVert) / (_current->_numVert + numVert);
+ return a * (_current->_valEdge + valEdge) / (_current->_numEdge + numEdge);
+}
+
+Rec2DAction* Rec2DData::getBestAction()
{
- _current->_quad.push_back(quad);
+ if (_current->_actions.size() == 0)
+ return NULL;
+ return *std::max_element(_current->_actions.begin(),
+ _current->_actions.end(), lessRec2DAction());
}
@@ -709,112 +496,76 @@ bool Rec2DAction::operator<(Rec2DAction &other)
return getReward() < other.getReward();
}
-MQuadrangle* Rec2DAction::_createQuad(std::vector<Rec2DEdge*> &boundEdge) const
+double Rec2DAction::getReward()
{
- //Msg::Info("Imhere");
- if (boundEdge.size() != 4) {
- Msg::Error("[Rec2DAction] Not 4 edges for quad creation");
- return new MQuadrangle(NULL, NULL, NULL, NULL);
- }
- MVertex *v1, *v2, *v3, *v4;
- v1 = boundEdge[0]->getRVertex(0)->getMVertex();
- v2 = boundEdge[0]->getRVertex(1)->getMVertex();
- int I;
- for (unsigned int i = 1; i < 4; ++i) {
- if (v2 == boundEdge[i]->getRVertex(0)->getMVertex()) {
- v3 = boundEdge[i]->getRVertex(1)->getMVertex();
- I = i;
- break;
- }
- if (v2 == boundEdge[i]->getRVertex(1)->getMVertex()) {
- v3 = boundEdge[i]->getRVertex(0)->getMVertex();
- I = i;
- break;
- }
- }
- for (unsigned int i = 1; i < 4; ++i) {
- if (i == I) continue;
- if (v3 == boundEdge[i]->getRVertex(0)->getMVertex()) {
- v4 = boundEdge[i]->getRVertex(1)->getMVertex();
- break;
- }
- if (v3 == boundEdge[i]->getRVertex(1)->getMVertex()) {
- v4 = boundEdge[i]->getRVertex(0)->getMVertex();
- break;
- }
- }
- return new MQuadrangle(v1, v2, v3, v4);
+ if (_lastUpdate < Recombine2D::getNumChange())
+ _computeGlobVal();
+ return _globValIfExecuted - Rec2DData::getGlobalValue();
}
/** Rec2DTwoTri2Quad **/
/************************/
-Rec2DTwoTri2Quad::Rec2DTwoTri2Quad(Rec2DEdge *re, std::vector<MTriangle*> &tri)
+Rec2DTwoTri2Quad::Rec2DTwoTri2Quad(Rec2DElement *el0, Rec2DElement *el1)
{
- int i;
- std::set<MEdge, Less_Edge> extEdges;
- {
- std::vector<MTriangle*>::iterator it = tri.begin();
- for (i = 0; it != tri.end() && i < 2; ++it, ++i) {
- _rtriangles[i] = Recombine2D::getRElement(*it);
- _rtriangles[i]->add(this);
- for (int j = 0; j < 3; ++j)
- extEdges.insert((*it)->getEdge(j));
- }
- if (it != tri.end() || i < 2)
- Msg::Error("[Rec2DTwoTri2Quad] Wrong number of triangles");
- extEdges.erase(re->getMEdge());
- }
+ _triangles[0] = el0;
+ _triangles[1] = el1;
+ _edges[4] = Rec2DElement::getCommonEdge(el0, el1);
- _edges[0] = re;
- std::set<MEdge>::iterator ite = extEdges.begin();
- for (i = 1; ite != extEdges.end() && i < 5; ++ite, ++i)
- _edges[i] = Recombine2D::getREdge(*ite);
- if (ite != extEdges.end() || i < 5)
- Msg::Error("[Rec2DTwoTri2Quad] Wrong number of edges");
+ std::vector<Rec2DEdge*> edges;
+ el0->getMoreEdges(edges);
+ el1->getMoreEdges(edges);
+ int k = -1;
+ for (unsigned int i = 0; i < edges.size(); ++i) {
+ if (edges[i] != _edges[4])
+ _edges[++k] = edges[i];
+ }
+ if (k > 3)
+ Msg::Error("[Rec2DTwoTri2Quad] too much edges");
- _vertices[0] = re->getRVertex(0);
- _vertices[1] = re->getRVertex(1);
- _vertices[2] = _rtriangles[0]->getOtherVertex(_vertices[0], _vertices[1]);
- _vertices[3] = _rtriangles[1]->getOtherVertex(_vertices[0], _vertices[1]);
+ _vertices[0] = _edges[4]->getVertex(0);
+ _vertices[1] = _edges[4]->getVertex(1);
+ _vertices[2] = _triangles[0]->getOtherVertex(_vertices[0], _vertices[1]);
+ _vertices[3] = _triangles[1]->getOtherVertex(_vertices[0], _vertices[1]);
- _computeGlobVal();
+ _triangles[0]->add(this);
+ _triangles[1]->add(this);
+ Rec2DData::add(this);
}
Rec2DTwoTri2Quad::~Rec2DTwoTri2Quad()
{
- _rtriangles[0]->remove(this);
- _rtriangles[1]->remove(this);
+ if (_triangles[0])
+ _triangles[0]->removeT(this);
+ if (_triangles[1])
+ _triangles[1]->removeT(this);
+ Rec2DData::remove(this);
}
-bool Rec2DTwoTri2Quad::isObsolete()
+void Rec2DTwoTri2Quad::_computeGlobVal()
{
- int p[4];
- p[0] = _vertices[0]->getParity();
- p[1] = _vertices[1]->getParity();
- p[2] = _vertices[2]->getParity();
- p[3] = _vertices[3]->getParity();
- return Rec2DTwoTri2Quad::isObsolete(p);
+ double valEdge = - REC2D_EDGE_BASE * _edges[4]->getQual();
+ for (int i = 0; i < 4; ++i)
+ valEdge += REC2D_EDGE_QUAD * _edges[i]->getQual();
+
+ double valVert = _vertices[0]->getGain(-1) + _vertices[1]->getGain(-1);
+
+ _globValIfExecuted =
+ Rec2DData::getGlobalValue(4*REC2D_EDGE_QUAD - REC2D_EDGE_BASE,
+ valEdge, 0, valVert );
+ _lastUpdate = Recombine2D::getNumChange();
}
-bool Rec2DTwoTri2Quad::isObsolete(int *p)
+void Rec2DTwoTri2Quad::color(int a, int b, int c)
{
- if (p[0]>-1 && p[1]>-1 && p[0]/2 == p[1]/2 && p[0]%2 != p[1]%2 ||
- p[2]>-1 && p[3]>-1 && p[2]/2 == p[3]/2 && p[2]%2 != p[3]%2 ||
- p[0]>-1 && (p[0] == p[2] || p[0] == p[3]) ||
- p[1]>-1 && (p[1] == p[2] || p[1] == p[3]) )
- return true;
- return false;
+ unsigned int col = CTX::instance()->packColor(a, b, c, 255);
+ _triangles[0]->getMElement()->setCol(col);
+ _triangles[1]->getMElement()->setCol(col);
}
-void Rec2DTwoTri2Quad::apply() //to check
+void Rec2DTwoTri2Quad::apply()
{
- if (isObsolete()) {
- Msg::Error("[Rec2DTwoTri2Quad] Applying obsolet action");
- return;
- }
-
- int min = 100, index = -1;
+ /*int min = 100, index = -1;
for (int i = 0; i < 4; ++i) {
if (_vertices[i]->getParity() > -1 && min > _vertices[i]->getParity()) {
min = _vertices[i]->getParity();
@@ -848,836 +599,635 @@ void Rec2DTwoTri2Quad::apply() //to check
Recombine2D::associateParity(_vertices[0]->getParity(), par);
}
}
- }
+ }*/
- double valEdge = - REC2D_EDGE_BASE * _edges[0]->getQual();
- for (int i = 1; i < 5; ++i)
- valEdge += REC2D_EDGE_QUAD * _edges[i]->getQual();
+ _triangles[0]->removeT(this);
+ _triangles[1]->removeT(this);
- Recombine2D::addNumEdge(4*REC2D_EDGE_QUAD - REC2D_EDGE_BASE);
- Recombine2D::addValEdge(valEdge);
- //Recombine2D::addValVert(_vertices[0]->getGain(-1) + _vertices[1]->getGain(-1));
+ std::vector<Rec2DAction*> actions;
+ _triangles[0]->getUniqueActions(actions);
+ _triangles[1]->getUniqueActions(actions);
+ for (unsigned int i = 0; i < actions.size(); ++i)
+ delete actions[i];
- Recombine2D::remove(_rtriangles[0]); // FIXME not performant at all !!!
- Recombine2D::remove(_rtriangles[1]);
- std::vector<Rec2DEdge*> edges;
- for (int i = 1; i < 5; ++i)
- edges.push_back(_edges[i]);
- Recombine2D::add(_createQuad(edges));
+ delete _triangles[0];
+ delete _triangles[1];
+ _triangles[0] = NULL;
+ _triangles[1] = NULL;
+ delete _edges[4];
+
+ /*new Rec2DCollapse(*/new Rec2DElement(_edges)/*)*/;
- _vertices[0]->remove(_edges[0]);
- _vertices[1]->remove(_edges[0]);
- Recombine2D::remove(_edges[0]->getMEdge()); // could also delete REdge
- delete _edges[0];
}
-void Rec2DTwoTri2Quad::getParities(int *par)
+
+/** Rec2DEdge **/
+/*****************/
+Rec2DEdge::Rec2DEdge(Rec2DVertex *rv0, Rec2DVertex *rv1)
+: _rv0(rv0), _rv1(rv1), _weight(REC2D_EDGE_BASE), _boundary(-1), _lastUpdate(-2)
{
- par[0] = _vertices[0]->getParity();
- par[1] = _vertices[1]->getParity();
- par[2] = _vertices[2]->getParity();
- par[3] = _vertices[3]->getParity();
+ _rv0->add(this);
+ _rv1->add(this);
+ Rec2DData::add(this);
+ _computeQual();
+ Rec2DData::addEdge(_weight, _weight*getQual());
}
-void Rec2DTwoTri2Quad::_computeGlobVal()
+Rec2DEdge::~Rec2DEdge()
{
- double valEdge = - REC2D_EDGE_BASE * _edges[0]->getQual();
- for (int i = 1; i < 5; ++i)
- valEdge += REC2D_EDGE_QUAD * _edges[i]->getQual();
-
- double valVert = _vertices[0]->getGain(-1) + _vertices[1]->getGain(-1);
-
- _globValIfExecuted = Recombine2D::getGlobalValue(4*REC2D_EDGE_QUAD - REC2D_EDGE_BASE,
- valEdge, 0, valVert );
- _lastUpdate = Recombine2D::getNumChange();
+ _rv0->remove(this);
+ _rv1->remove(this);
+ Rec2DData::remove(this);
+ Rec2DData::addEdge(-_weight, -_weight*getQual());
}
-void Rec2DTwoTri2Quad::color(int a, int b, int c)
+void Rec2DEdge::_computeQual() //*
{
- unsigned int col = CTX::instance()->packColor(a, b, c, 255);
- _rtriangles[0]->getMElement()->setCol(col);
- _rtriangles[1]->getMElement()->setCol(col);
+ double adimLength = _straightAdimLength();
+ double alignment = _straightAlignment();
+ if (adimLength > 1)
+ adimLength = 1/adimLength;
+ _qual = adimLength * ((1-REC2D_ALIGNMENT) + REC2D_ALIGNMENT * alignment);
+ _lastUpdate = Recombine2D::getNumChange();
+}
+
+double Rec2DEdge::getQual()
+{
+ if (_lastUpdate < Recombine2D::getNumChange() &&
+ _rv0->getLastMove() > _lastUpdate ||
+ _rv1->getLastMove() > _lastUpdate ) {
+ _computeQual();
+ }
+ return _qual;
}
+double Rec2DEdge::addNeighbour()
+{
+ if (++_boundary > 1)
+ Msg::Error("[Rec2DEdge] Too much boundary element");
+}
-bool Rec2DTwoTri2Quad::whatImply(
- Rec2DVertex **rv, int *par,
- std::map<Rec2DVertex*, std::vector<int> > &suggestions)
+double Rec2DEdge::addQuadNeighbour()
{
- int num = 0, p[4];
- for (int i = 0; i < 4; ++i) {
- bool b = false;
- for (int j = 0; j < 3; ++j) {
- if (_vertices[i] == rv[j]) {
- ++num;
- p[i] = par[j];
- b = true;
- }
- }
- if (!b)
- p[i] = _vertices[i]->getParity();
+ _weight += REC2D_EDGE_QUAD;
+ if (_weight > REC2D_EDGE_BASE + 2*REC2D_EDGE_QUAD)
+ Msg::Error("[Rec2DEdge] Weight too higher : %d (%d max)",
+ _weight, REC2D_EDGE_BASE + 2*REC2D_EDGE_QUAD );
+ Rec2DData::addEdge(REC2D_EDGE_QUAD, REC2D_EDGE_QUAD*getQual());
+}
+
+
+Rec2DElement* Rec2DEdge::getSingleElement(Rec2DEdge *re)
+{
+ std::vector<Rec2DElement*> elem;
+ Rec2DVertex::getCommonElements(re->getVertex(0), re->getVertex(1), elem);
+ if (elem.size() == 1)
+ return elem[0];
+ if (elem.size() != 0)
+ Msg::Error("[Rec2DEdge] Edge bound %d elements, returning NULL", elem.size());
+ return NULL;
+}
+
+void Rec2DEdge::swap(Rec2DVertex *oldRV, Rec2DVertex *newRV)
+{
+ if (_rv0 == oldRV) {
+ _rv0 = newRV;
+ return;
}
- if (num != 3)
- Msg::Error("[Rec2DTwoTri2Quad] %d corresponding vertices instead of 3", num);
-
- if (Rec2DTwoTri2Quad::isObsolete(p)) {
- //Msg::Info(" obsolete action");
- return false;
+ if (_rv1 == oldRV) {
+ _rv1 = newRV;
+ return;
}
+ Msg::Error("[Rec2DEdge] oldRV not found");
+}
+
+double Rec2DEdge::_straightAdimLength() const
+{
+ double dx, dy, dz, *xyz0 = new double[3], *xyz1 = new double[3];
+ _rv0->getxyz(xyz0);
+ _rv1->getxyz(xyz1);
+ dx = xyz0[0] - xyz1[0];
+ dy = xyz0[1] - xyz1[1];
+ dz = xyz0[2] - xyz1[2];
+ double length = sqrt(dx * dx + dy * dy + dz * dz);
+ delete[] xyz0;
+ delete[] xyz1;
- //Msg::Info(" [%d %d %d %d]", p[0], p[1], p[2], p[3]);
+ double lc0 = (*Recombine2D::bgm())(_rv0->u(), _rv0->v(), .0);
+ double lc1 = (*Recombine2D::bgm())(_rv1->u(), _rv1->v(), .0);
- int min = 100, index = -1;
- for (int i = 0; i < 4; ++i) {
- if (p[i] > -1 && min > p[i]) {
- min = p[i];
- index = i;
- }
- }
- if (index == -1) {
- Msg::Error("[Rec2DTwoTri2Quad] no parities");
- return false;
- }
+ return length * (1/lc0 + 1/lc1) / 2;
+}
+
+double Rec2DEdge::_straightAlignment() const
+{
+ double angle0 = Recombine2D::bgm()->getAngle(_rv0->u(), _rv0->v(), .0);
+ double angle1 = Recombine2D::bgm()->getAngle(_rv1->u(), _rv1->v(), .0);
+ double angleEdge = atan2(_rv0->u()-_rv1->u(), _rv0->v()-_rv1->v());
- for (int i = 0; i < 4; i += 2) {
- int par;
- if ((index/2)*2 == i)
- par = min;
- else
- par = otherParity(min);
- for (int j = 0; j < 2; ++j) {
- if (p[i+j] != par)
- suggestions[_vertices[i+j]].push_back(par);
- }
- }
- //Msg::Info(" ok");
- return true;
-}
-
-
-///** Rec2DFourTri2Quad **/
-///*************************/
-//Rec2DFourTri2Quad::Rec2DFourTri2Quad(Rec2DVertex *rv, std::vector<MTriangle*> &tri)
-//{
-// int i, j;
-// std::set<MEdge, Less_Edge> edges;
-// {
-// std::vector<MTriangle*>::iterator it = tri.begin();
-// for (i = 0; it != tri.end() && i < 4; ++it, ++i) {
-// for (j = 0; j < 3; ++j)
-// edges.insert((*it)->getEdge(j));
-// _rtriangles[i] = Recombine2D::getRElement(*it);
-// _rtriangles[i]->add(this);
-// }
-// if (it != tri.end() || i < 4)
-// Msg::Error("[Rec2DFourTri2Quad] Wrong number of triangles");
-// }
-//
-// std::set<MEdge>::iterator ite = edges.begin();
-// for (i = 0, j = 4; ite != edges.end() && (i < 4 || j < 8); ++ite) {
-// if ((*ite).getVertex(0) == rv->getMVertex() ||
-// (*ite).getVertex(1) == rv->getMVertex() )
-// _edges[i++] = Recombine2D::getREdge(*ite);
-// else if (j < 8)
-// _edges[j++] = Recombine2D::getREdge(*ite);
-// else
-// Msg::Error("[Rec2DFourTri2Quad] Too much exterior edges");
-// }
-// if (edges.size() > 8 || ite != edges.end() || i < 4 || i > 4 || j < 8)
-// Msg::Error("[Rec2DFourTri2Quad] Wrong number of edges");
-//
-// _vertices[4] = rv;
-// // the 4 other must be in order : 2 non adjacent + last 2
-// _vertices[1] = _edges[4]->getRVertex(0);
-// _vertices[3] = _edges[4]->getRVertex(1);
-// for (int i = 5; i < 8; ++i) {
-// if (_edges[i]->getRVertex(0) == _vertices[1])
-// _vertices[0] = _edges[i]->getRVertex(1);
-// if (_edges[i]->getRVertex(1) == _vertices[1])
-// _vertices[0] = _edges[i]->getRVertex(0);
-// if (_edges[i]->getRVertex(0) == _vertices[3])
-// _vertices[2] = _edges[i]->getRVertex(1);
-// if (_edges[i]->getRVertex(1) == _vertices[3])
-// _vertices[2] = _edges[i]->getRVertex(0);
-// }
-//
-// _computeGlobVal();
-//}
-//
-//Rec2DFourTri2Quad::~Rec2DFourTri2Quad()
-//{
-// _rtriangles[0]->remove(this);
-// _rtriangles[1]->remove(this);
-// _rtriangles[2]->remove(this);
-// _rtriangles[3]->remove(this);
-//}
-//
-//bool Rec2DFourTri2Quad::isObsolete()
-//{
-// int p[4];
-// p[0] = _vertices[0]->getParity();
-// p[1] = _vertices[1]->getParity();
-// p[2] = _vertices[2]->getParity();
-// p[3] = _vertices[3]->getParity();
-// return Rec2DTwoTri2Quad::isObsolete(p);
-//}
-//
-//void Rec2DFourTri2Quad::apply()
-//{
-// double valEdge = 0;
-// for (int i = 0; i < 4; ++i)
-// valEdge -= REC2D_EDGE_BASE * _edges[i]->getQual();
-// for (int i = 4; i < 8; ++i)
-// valEdge += REC2D_EDGE_QUAD * _edges[i]->getQual();
-//
-// double valVert = - _vertices[4]->getQual();
-// for (int i = 0; i < 4; ++i)
-// valEdge += _vertices[i]->getGain(-1);
-//
-// Recombine2D::addNumEdge(4*REC2D_EDGE_QUAD - 4*REC2D_EDGE_BASE);
-// Recombine2D::addNumVert(-1);
-// Recombine2D::addValEdge(valEdge);
-// Recombine2D::addValVert(valVert);
-//
-// Recombine2D::remove(_rtriangles[0]); // FIXME not performant at all !!!
-// Recombine2D::remove(_rtriangles[1]);
-// Recombine2D::remove(_rtriangles[2]);
-// Recombine2D::remove(_rtriangles[3]);
-// std::vector<Rec2DEdge*> edges;
-// for (int i = 4; i < 8; ++i)
-// edges.push_back(_edges[i]);
-// Recombine2D::add(_createQuad(edges));
-//
-// for (int i = 0; i < 4; ++i) {
-// _edges[i]->getRVertex(0)->remove(_edges[i]);
-// _edges[i]->getRVertex(1)->remove(_edges[i]);
-// Recombine2D::remove(_edges[i]->getMEdge());
-// delete _edges[i];
-// }
-// Recombine2D::remove(_vertices[4]->getMVertex());
-// delete _vertices[4];
-//}
-//
-//void Rec2DFourTri2Quad::getParities(int *par)
-//{
-// par[0] = _vertices[0]->getParity();
-// par[1] = _vertices[1]->getParity();
-// par[2] = _vertices[2]->getParity();
-// par[3] = _vertices[3]->getParity();
-//}
-//
-//void Rec2DFourTri2Quad::_computeGlobVal()
-//{
-// double valEdge = 0;
-// for (int i = 0; i < 4; ++i)
-// valEdge -= REC2D_EDGE_BASE * _edges[i]->getQual();
-// for (int i = 4; i < 8; ++i)
-// valEdge += REC2D_EDGE_QUAD * _edges[i]->getQual();
-//
-// double valVert = - _vertices[4]->getQual();
-// for (int i = 0; i < 4; ++i)
-// valEdge += _vertices[i]->getGain(-1);
-//
-// _globValIfExecuted = Recombine2D::getGlobalValue(4*REC2D_EDGE_QUAD - 4*REC2D_EDGE_BASE,
-// valEdge, -1, valVert );
-// _lastUpdate = Recombine2D::getNumChange();
-//}
-//
-//void Rec2DFourTri2Quad::color(int a, int b, int c)
-//{
-// unsigned int col = CTX::instance()->packColor(a, b, c, 255);
-// _rtriangles[0]->getMElement()->setCol(col);
-// _rtriangles[1]->getMElement()->setCol(col);
-// _rtriangles[2]->getMElement()->setCol(col);
-// _rtriangles[3]->getMElement()->setCol(col);
-//}
-//
-//
-///** Rec2DCollapseTri **/
-///*********************/
-//Rec2DCollapseTri::Rec2DCollapseTri(Rec2DVertex *rv, std::vector<MTriangle*> &tri)
-//{
-// int i;
-// std::set<MEdge, Less_Edge> extEdges;
-// {
-// std::vector<MTriangle*>::iterator it = tri.begin();
-// for (i = 0; it != tri.end() && i < 4; ++it, ++i) {
-// for (int j = 0; j < 3; ++j) {
-// if ((*it)->getEdge(j).getVertex(0) != rv->getMVertex() &&
-// (*it)->getEdge(j).getVertex(1) != rv->getMVertex() )
-// extEdges.insert((*it)->getEdge(j));
-// }
-// _rtriangles[i] = Recombine2D::getRElement(*it);
-// _rtriangles[i]->add(this);
-// }
-// if (it != tri.end() || i < 4)
-// Msg::Error("[Rec2DFourTri2Quad] Wrong number of triangles");
-// }
-//
-// _vertices[4] = rv;
-// // the 4 other must be in order : 2 non adjacent + last 2
-// std::set<MEdge>::iterator ite = extEdges.begin();
-// Rec2DEdge *re = Recombine2D::getREdge(*ite);
-// _vertices[1] = re->getRVertex(0);
-// _vertices[3] = re->getRVertex(1);
-// ++ite;
-// for (; ite != extEdges.end(); ++ite) {
-// re = Recombine2D::getREdge(*ite);
-// if (re->getRVertex(0) == _vertices[1])
-// _vertices[0] = re->getRVertex(1);
-// if (re->getRVertex(1) == _vertices[1])
-// _vertices[0] = re->getRVertex(0);
-// if (re->getRVertex(0) == _vertices[3])
-// _vertices[2] = re->getRVertex(1);
-// if (re->getRVertex(1) == _vertices[3])
-// _vertices[2] = re->getRVertex(0);
-// }
-//
-// _computeGlobVal();
-//}
-//
-//Rec2DCollapseTri::Rec2DCollapseTri(Rec2DEdge *re, std::vector<MTriangle*> &tri)
-//: _edge(re)
-//{
-// int i;
-// std::set<MEdge, Less_Edge> extEdges;
-// {
-// std::vector<MTriangle*>::iterator it = tri.begin();
-// for (i = 0; it != tri.end() && i < 2; ++it, ++i) {
-// for (int j = 0; j < 3; ++j)
-// extEdges.insert((*it)->getEdge(j));
-// _rtriangles[i] = Recombine2D::getRElement(*it);
-// _rtriangles[i]->add(this);
-// }
-// if (it != tri.end() || i < 2)
-// Msg::Error("[Rec2DTwoTri2Quad] Wrong number of triangles");
-// extEdges.erase(re->getMEdge());
-// }
-// _rtriangles[2] = _rtriangles[3] = NULL;
-//
-// _vertices[0] = re->getRVertex(0);
-// _vertices[1] = re->getRVertex(1);
-// _vertices[2] = _rtriangles[0]->getOtherVertex(_vertices[0], _vertices[1]);
-// _vertices[3] = _rtriangles[1]->getOtherVertex(_vertices[0], _vertices[1]);
-// _vertices[4] = NULL;
-//
-// _computeGlobVal();
-//}
-//
-//Rec2DCollapseTri::~Rec2DCollapseTri()
-//{
-// //_rtriangles[0]->remove(this);
-// //_rtriangles[1]->remove(this);
-//}
-//
-//bool Rec2DCollapseTri::isObsolete()
-//{
-// int p0 = _vertices[0]->getParity();
-// int p1 = _vertices[1]->getParity();
-// int p2 = _vertices[2]->getParity();
-// int p3 = _vertices[3]->getParity();
-// int b0 = _vertices[0]->getIsOnBoundary();
-// int b1 = _vertices[1]->getIsOnBoundary();
-// int b2 = _vertices[2]->getIsOnBoundary();
-// int b3 = _vertices[3]->getIsOnBoundary();
-// if ((b0 && b1) || (p0>-1 && p1>-1 && p0/2 == p1/2 && p0 % 2 != p1 % 2) &&
-// (b2 && b3) || (p2>-1 && p3>-1 && p2/2 == p3/2 && p2 % 2 != p3 % 2) )
-// return true;
-// return false;
-//}
-//
-//void Rec2DCollapseTri::apply()
-//{
-// return;
-//}
-//
-//void Rec2DCollapseTri::getParities(int *par)
-//{
-// par[0] = _vertices[0]->getParity();
-// par[1] = _vertices[1]->getParity();
-// par[2] = _vertices[2]->getParity();
-// par[3] = _vertices[3]->getParity();
-//}
-//
-//void Rec2DCollapseTri::_computeGlobVal()
-//{
-// int b0 = _vertices[0]->getIsOnBoundary();
-// int b1 = _vertices[1]->getIsOnBoundary();
-// int b2 = _vertices[2]->getIsOnBoundary();
-// int b3 = _vertices[3]->getIsOnBoundary();
-//
-// std::map<Rec2DVertex*, int> neighbourVert;
-// std::set<Rec2DEdge*> neighbourEdge;
-// if (!b0 || !b1) {
-// _vertices[0]->getNeighbours(neighbourVert, neighbourEdge);
-// _vertices[1]->getNeighbours(neighbourVert, neighbourEdge);
-// if (_vertices[4]) {
-// neighbourVert.erase(_vertices[4]);
-// _vertices[4]->getNeighbourEdges(neighbourEdge);
-// }
-// //else
-// // neighbourEdge.insert(_edge);
-//
-// int numEdge = 0, numVert = 0;
-// double valEdge = .0, valVert = .0;
-// {
-// std::set<Rec2DEdge*>::iterator it = neighbourEdge.begin();
-// for (; it != neighbourEdge.end(); ++it) {
-// valEdge -= (*it)->getWeightedQual();
-// numEdge -= (*it)->getWeight();
-// }
-// valVert -= _vertices[0]->getQual();
-// valVert -= _vertices[1]->getQual();
-// if (_vertices[4]) {
-// valVert -= _vertices[4]->getQual();
-// numVert -= 1;
-// valVert += _vertices[2]->getGain(-2);
-// valVert += _vertices[3]->getGain(-2);
-// }
-// else {
-// valVert += _vertices[2]->getGain(-1);
-// valVert += _vertices[3]->getGain(-1);
-// }
-// }
-// if (b0)
-// _qualCavity(valVert, valEdge, numEdge, neighbourVert, _vertices[0]);
-// else if (b1)
-// _qualCavity(valVert, valEdge, numEdge, neighbourVert, _vertices[1]);
-// else
-// _qualCavity(valVert, valEdge, numEdge, neighbourVert);
-// numVert -= 1;
-//
-// _globValIfExecuted =
-// Recombine2D::getGlobalValue(numEdge, valEdge, numVert, valVert);
-// }
-// else
-// _globValIfExecuted = -1.;
-//
-//
-// neighbourVert.clear();
-// neighbourEdge.clear();
-// if (!b2 || !b3) {
-// _vertices[2]->getNeighbours(neighbourVert, neighbourEdge);
-// _vertices[3]->getNeighbours(neighbourVert, neighbourEdge);
-// if (_vertices[4]) {
-// neighbourVert.erase(_vertices[4]);
-// _vertices[4]->getNeighbourEdges(neighbourEdge);
-// }
-// else
-// neighbourEdge.insert(_edge);
-//
-// int numEdge = 0, numVert = 0;
-// double valEdge = .0, valVert = .0;
-// {
-// std::set<Rec2DEdge*>::iterator it = neighbourEdge.begin();
-// for (; it != neighbourEdge.end(); ++it) {
-// valEdge -= (*it)->getWeightedQual();
-// numEdge -= (*it)->getWeight();
-// }
-// valVert -= _vertices[2]->getQual();
-// valVert -= _vertices[3]->getQual();
-// if (_vertices[4]) {
-// valVert -= _vertices[4]->getQual();
-// numVert -= 1;
-// valVert += _vertices[0]->getGain(-2);
-// valVert += _vertices[1]->getGain(-2);
-// }
-// else {
-// valVert += _vertices[0]->getGain(-2);
-// valVert += _vertices[1]->getGain(-2);
-// }
-// }
-// if (b2)
-// _qualCavity(valVert, valEdge, numEdge, neighbourVert, _vertices[2]);
-// else if (b3)
-// _qualCavity(valVert, valEdge, numEdge, neighbourVert, _vertices[3]);
-// else
-// _qualCavity(valVert, valEdge, numEdge, neighbourVert);
-// numVert -= 1;
-//
-// _globValIfExecuted2 =
-// Recombine2D::getGlobalValue(numEdge, valEdge, numVert, valVert);
-// }
-// else
-// _globValIfExecuted2 = -1.;
-//
-// _lastUpdate = Recombine2D::getNumChange();
-//}
-//
-//void Rec2DCollapseTri::_qualCavity(double &valVert, double &valEdge,
-// int &numEdge,
-// std::map<Rec2DVertex*, int> &vert,
-// Rec2DVertex *imposedVert )
-//{
-// Recombine2D::bo = true;
-// std::map<Rec2DVertex*, int>::iterator it;
-// Rec2DVertex *rv = imposedVert;
-// if (rv == NULL) {
-// double u, v = u = .0;
-// it = vert.begin();
-// for (; it != vert.end(); ++it) {
-// u += it->first->u();
-// v += it->first->v();
-// }
-// u /= vert.size();
-// v /= vert.size();
-//
-// rv = new Rec2DVertex(u, v);
-// }
-// valVert = rv->getQual(vert.size());
-// it = vert.begin();
-// for (; it != vert.end(); ++it) {
-// Rec2DEdge re(rv, it->first);
-// valEdge += it->second * re.getQual();
-// numEdge += it->second;
-// }
-// if (imposedVert == NULL) {
-// rv->deleteMVertex();
-// delete rv;
-// }
-//}
-//
-//void Rec2DCollapseTri::color(int a, int b, int c)
-//{
-// unsigned int col = CTX::instance()->packColor(a, b, c, 255);
-// _rtriangles[0]->getMElement()->setCol(col);
-// _rtriangles[1]->getMElement()->setCol(col);
-// if (_rtriangles[2]) {
-// _rtriangles[2]->getMElement()->setCol(col);
-// _rtriangles[3]->getMElement()->setCol(col);
-// }
-//}
-//
-//
+ double alpha0 = angleEdge - angle0;
+ double alpha1 = angleEdge - angle1;
+ crossField2d::normalizeAngle(alpha0);
+ crossField2d::normalizeAngle(alpha1);
+ alpha0 = 1 - 4 * std::min(alpha0, .5 * M_PI - alpha0) / M_PI;
+ alpha1 = 1 - 4 * std::min(alpha1, .5 * M_PI - alpha1) / M_PI;
+
+ double lc0 = (*Recombine2D::bgm())(_rv0->u(), _rv0->v(), .0);
+ double lc1 = (*Recombine2D::bgm())(_rv1->u(), _rv1->v(), .0);
+
+ return (alpha0/lc0 + alpha1/lc1) / (1/lc0 + 1/lc1);
+}
+
+Rec2DVertex* Rec2DEdge::getOtherVertex(Rec2DVertex *rv) const
+{
+ if (rv == _rv0)
+ return _rv1;
+ if (rv == _rv1)
+ return _rv0;
+ Msg::Error("[Rec2DEdge] You are wrong, I don't have this vertex !");
+ return NULL;
+}
+
+
/** Rec2DVertex **/
/*******************/
-Rec2DVertex::Rec2DVertex(MVertex *v, std::vector<MTriangle*> tri, int onWhat,
- std::map<MVertex*, std::set<GEdge*> > gEdges)
-: _v(v), _onWhat(onWhat), _parity(-1), _lastChange(-1), _numEl(tri.size())
+Rec2DVertex::Rec2DVertex(MVertex *v, bool toSave)
+: _v(v), _onWhat(1), _parity(-1), _lastMove(-1), _angle(4*M_PI), _isMesh(toSave)
{
- for (unsigned int i = 0; i < tri.size(); ++i)
- _elements.push_back(Recombine2D::getRElement(tri[i]));
-
- if (_onWhat == -1) {
- std::map<MVertex*, std::set<GEdge*> >::iterator it = gEdges.find(_v);
- if (it != gEdges.end())
- _angle = _computeAngle(_v, tri, it->second);
- else {
- Msg::Warning("[meshRecombine2D] Can't compute corner angle, setting angle to pi/2");
- _angle = M_PI / 2.;
- }
- }
reparamMeshVertexOnFace(_v, Recombine2D::getGFace(), _param);
- _computeQual();
+ if (_isMesh) {
+ Rec2DData::add(this);
+ Rec2DData::addVert(1, getQual());
+ }
}
-Rec2DVertex::Rec2DVertex(double u, double v)
-: _onWhat(1), _parity(-1), _lastChange(-1)
+Rec2DVertex::Rec2DVertex(Rec2DVertex *rv, double ang)
+: _v(rv->_v), _onWhat(-1), _parity(-1), _lastMove(-1), _angle(ang), _isMesh(true)
{
- _param[0] = u;
- _param[1] = v;
-
- GPoint p = Recombine2D::getGFace()->point(_param);
- _v = new MVertex(p.x(), p.y(), p.z(), Recombine2D::getGFace());
- static int i = 0;
- if (++i < 20)
- Msg::Info("resulting point = [%g, %g, %g]", p.x(), p.y(), p.z());
+ _edges = rv->_edges;
+ _elements = rv->_elements;
+ _param = rv->_param;
+ for (int i = 0; i < _edges.size(); ++i) {
+ _edges[i]->swap(rv, this);
+ }
+ Rec2DData::add(this);
+ Rec2DData::addVert(1, getQual());
+ delete rv;
}
-double Rec2DVertex::_computeAngle(MVertex *v,
- std::vector<MTriangle*> &tri,
- std::set<GEdge*> &setEdge)
+Rec2DVertex::~Rec2DVertex()
{
- if (setEdge.size() != 2) {
- Msg::Warning("[meshGFaceRecombine] Wrong number of edge : %d, returning pi/2", setEdge.size());
- return M_PI / 2.;
+ Rec2DData::remove(this);
+ Rec2DData::addVert(-1, -getQual());
+}
+
+void Rec2DVertex::initStaticTable()
+{
+ // _qualVSnum[onWhat][numEl]; onWhat={0:edge,1:face}
+ // _gains[onWhat][numEl]; onWhat={0:edge,1:face} (earning of adding an element)
+ if (_qualVSnum == NULL) {
+ int nE = 10, nF = 20; //edge, face
+
+ _qualVSnum = new double*[2];
+ _qualVSnum[0] = new double[nE];
+ _qualVSnum[1] = new double[nF];
+ _qualVSnum[0][0] = -10.;
+ _qualVSnum[1][0] = -10.;
+ for (int i = 1; i < nE; ++i)
+ _qualVSnum[0][i] = 1. - fabs(2./i - 1.);
+ for (int i = 1; i < nF; ++i)
+ _qualVSnum[1][i] = 1. - fabs(4./i - 1.);
+
+ _gains = new double*[2];
+ _gains[0] = new double[nE-1];
+ for (int i = 0; i < nE-1; ++i)
+ _gains[0][i] = _qualVSnum[0][i+1] - _qualVSnum[0][i];
+ _gains[1] = new double[nF-1];
+ for (int i = 0; i < nF-1; ++i)
+ _gains[1][i] = _qualVSnum[1][i+1] - _qualVSnum[1][i];
}
- static const double prec = 100.;
-
- SVector3 vectv = SVector3(v->x(), v->y(), v->z());
- SVector3 firstDer0, firstDer1;
-
- std::set<GEdge*>::iterator it = setEdge.begin();
- for (int k = 0; k < 2; ++it, ++k) {
- GEdge *ge = *it;
- SVector3 vectlb = ge->position(ge->getLowerBound());
- SVector3 vectub = ge->position(ge->getUpperBound());
- vectlb -= vectv;
- vectub -= vectv;
- double normlb, normub;
- if ((normlb = norm(vectlb)) > prec * (normub = norm(vectub)))
- firstDer1 = ge->firstDer(ge->getUpperBound());
- else if (normub > prec * normlb)
- firstDer1 = ge->firstDer(ge->getLowerBound());
- else {
- Msg::Warning("[meshGFaceRecombine] Bad precision, returning pi/2");
- return M_PI / 2.;
- }
- if (k == 0)
- firstDer0 = firstDer1;
+}
+
+Rec2DEdge* Rec2DVertex::getCommonEdge(Rec2DVertex *rv0, Rec2DVertex *rv1)
+{
+ for (unsigned int i = 0; i < rv0->_edges.size(); ++i) {
+ if (rv1->has(rv0->_edges[i]))
+ return rv0->_edges[i];
}
-
- firstDer0 = firstDer0 * (1 / norm(firstDer0));
- firstDer1 = firstDer1 * (1 / norm(firstDer1));
-
- double angle1 = acos(dot(firstDer0, firstDer1));
- double angle2 = 2. * M_PI - angle1;
-
- double angleMesh = .0;
- std::vector<MTriangle*>::iterator it2 = tri.begin();
- for (; it2 != tri.end(); ++it2) {
- MTriangle *t = *it2;
- int k = 0;
- while (t->getVertex(k) != v && k < 3)
- ++k;
- if (k == 3) {
- Msg::Warning("[meshGFaceRecombine] Wrong triangle, returning pi/2");
- return M_PI / 2.;
+ //Msg::Warning("[Rec2DVertex] didn't find edge, returning NULL");
+ return NULL;
+}
+
+void Rec2DVertex::getCommonElements(Rec2DVertex *rv0, Rec2DVertex *rv1,
+ std::vector<Rec2DElement*> &elem )
+{
+ for (unsigned int i = 0; i < rv0->_elements.size(); ++i) {
+ if (rv1->has(rv0->_elements[i]))
+ elem.push_back(rv0->_elements[i]);
+ }
+}
+
+bool Rec2DVertex::setBoundaryParity(int p0, int p1)
+{
+ if (_parity > -1) {
+ Msg::Error("[Rec2DVertex] Are you kidding me ? Already have a parity !");
+ return false;
+ }
+ setParity(p0);
+ int num = 0;
+ Rec2DVertex *nextRV = NULL;
+ for (unsigned int i = 0; i < _edges.size(); ++i) {
+ if (_edges[i]->isOnBoundary()) {
+ nextRV = _edges[i]->getOtherVertex(this);
+ ++num;
}
- angleMesh += angle3Vertices(t->getVertex((k+1)%3), v, t->getVertex((k+2)%3));
}
-
- if (angleMesh < M_PI)
- return angle1;
- return angle2;
+ if (num != 2)
+ Msg::Error("[Rec2DVertex] What's happening ? Am I on boundary or not ? TELL ME !");
+ if (nextRV)
+ return nextRV->_recursiveBoundParity(this, p1, p0); // alternate parity
+ Msg::Error("[Rec2DVertex] Have I really to say that I didn't find neighbouring vertex ?");
+ return false;
}
-double Rec2DVertex::getQual(int numEdge)
+bool Rec2DVertex::_recursiveBoundParity(Rec2DVertex *prevRV, int p0, int p1)
{
- if (numEdge > 1) {
- switch (_onWhat) {
- case 1 :
- return _qualVSnum[_onWhat][numEdge-1];
- case 0 :
- return _qualVSnum[_onWhat][numEdge-2];
- case -1 :
- return 1. - fabs(2./M_PI * _angle/(numEdge-1) - 1.);
- default :
- Msg::Error("[Rec2DVertex] Don't know on what is the vertex");
- return -1.;
+ if (_parity == p0)
+ return true;
+ if (_parity > -1) {
+ Msg::Error("[Rec2DVertex] Sorry, have parity %d, can't set it to %d... "
+ "You failed ! Try again !", _parity, p0);
+ return false;
+ }
+ setParity(p0);
+ int num = 0;
+ Rec2DVertex *nextRV = NULL;
+ for (unsigned int i = 0; i < _edges.size(); ++i) {
+ if (_edges[i]->isOnBoundary() && _edges[i]->getOtherVertex(this) != prevRV) {
+ nextRV = _edges[i]->getOtherVertex(this);
+ ++num;
}
}
- if (_lastChange < Recombine2D::getNumChange())
- _computeQual();
- return _qual;
+ if (num != 1)
+ Msg::Error("[Rec2DVertex] Holy shit !");
+ if (nextRV)
+ return nextRV->_recursiveBoundParity(this, p1, p0); // alternate parity
+ Msg::Error("[Rec2DVertex] Have I really to say that I didn't find next vertex ?");
+ return false;
}
-bool Rec2DVertex::hasTriangle()
+void Rec2DVertex::setParity(int p)
{
- bool b = false;
- for (int i = 0; i < _elements.size(); ++i)
- if (!_elements[i]->isQuad()) b = true;
- return b;
+ if (_parity > -1) {
+ Rec2DData::removeParity(this, _parity);
+ }
+ _parity = p;
+ Rec2DData::addParity(this, _parity);
}
-void Rec2DVertex::_computeQual()
+double Rec2DVertex::getQual(int numEl) const
{
+ int nEl = numEl > -1 ? numEl : _elements.size();
if (_onWhat > -1)
- _qual = _qualVSnum[_onWhat][_numEl-1];
- else
- _qual = 1. - fabs(2./M_PI * _angle/_numEl - 1.);
- _lastChange = Recombine2D::getNumChange();
+ return _qualVSnum[_onWhat][nEl];
+ if (nEl == 0)
+ return -10.;
+ return 1. - fabs(2./M_PI * _angle/nEl - 1.);
}
-double Rec2DVertex::getGain(int n)
+double Rec2DVertex::getGain(int n) const
{
+ if (!n)
+ return .0;
+ if (_elements.size() + n < 0) {
+ Msg::Error("[Rec2DVertex] gain for %d elements unavailable",
+ _elements.size() + n );
+ return .0;
+ }
+
if (_onWhat > -1) {
switch (n) {
case 1 :
- return _gains[_onWhat][_numEl-1];
+ return _gains[_onWhat][_elements.size()];
case -1 :
- return - _gains[_onWhat][_numEl-2];
+ return - _gains[_onWhat][_elements.size()-1];
default :
- return _qualVSnum[_onWhat][_numEl-1+n] - _qualVSnum[_onWhat][_numEl-1];
+ return _qualVSnum[_onWhat][_elements.size()+n]
+ - _qualVSnum[_onWhat][_elements.size()-1];
}
}
- else
- return fabs(2./M_PI * _angle/_numEl - 1.)
- - fabs(2./M_PI * _angle/(_numEl + n) - 1.);
+
+ if (_elements.size() == 0)
+ return 11 - fabs(2./M_PI * _angle/(_elements.size() + n) - 1.);
+
+ if (_elements.size() + n == 0)
+ return fabs(2./M_PI * _angle/_elements.size() - 1.) - 11;
+
+ return fabs(2./M_PI * _angle/_elements.size() - 1.)
+ - fabs(2./M_PI * _angle/(_elements.size() + n) - 1.);
}
-void Rec2DVertex::getNeighbours(std::map<Rec2DVertex*, int> &vert,
- std::set<Rec2DEdge*> &edge )
+void Rec2DVertex::add(Rec2DEdge *re)
{
- Rec2DVertex *rv;
- std::map<Rec2DVertex*, int>::iterator it;
- for (unsigned int i = 0; i < _edges.size(); ++i){
- edge.insert(_edges[i]);
- if (_edges[i]->getRVertex(0) != this)
- rv = _edges[i]->getRVertex(0);
- else
- rv = _edges[i]->getRVertex(1);
- if ((it = vert.find(rv)) != vert.end())
- it->second += _edges[i]->getWeight() - REC2D_EDGE_BASE;
- else
- vert[rv] = _edges[i]->getWeight();
+ for (unsigned int i = 0; i < _edges.size(); ++i) {
+ if (_edges[i] == re) {
+ Msg::Warning("[Rec2DVertex] Edge was already there");
+ return;
+ }
}
+ _edges.push_back(re);
}
-void Rec2DVertex::getNeighbourEdges(std::set<Rec2DEdge*> &edge)
+bool Rec2DVertex::has(Rec2DEdge *re) const
{
- for (unsigned int i = 0; i < _edges.size(); ++i)
- edge.insert(_edges[i]);
+ for (unsigned int i = 0; i < _edges.size(); ++i) {
+ if (_edges[i] == re)
+ return true;
+ }
+ return false;
}
-void Rec2DVertex::initStaticTable()
+void Rec2DVertex::remove(Rec2DEdge *re)
{
- if (_qualVSnum == NULL) {
- // _qualVSnum[onWhat][nEl]; onWhat={0:edge,1:face} / nEl=_numEl-1
- // _gains[onWhat][nEl]; onWhat={0:edge,1:face} / nEl=_numEl-1 (earning of adding an element)
- int nE = 10, nF = 20;
-
- _qualVSnum = new double*[2];
- _qualVSnum[0] = new double[nE];
- for (int i = 0; i < nE; ++i)
- _qualVSnum[0][i] = 1. - fabs(2./(i+1) - 1.);
- _qualVSnum[1] = new double[nF];
- for (int i = 0; i < nF; ++i)
- _qualVSnum[1][i] = 1. - fabs(4./(i+1) - 1.);
-
- _gains = new double*[2];
- _gains[0] = new double[nE-1];
- for (int i = 0; i < nE-1; ++i)
- _gains[0][i] = _qualVSnum[0][i+1] - _qualVSnum[0][i];
- _gains[1] = new double[nF-1];
- for (int i = 0; i < nF-1; ++i)
- _gains[1][i] = _qualVSnum[1][i+1] - _qualVSnum[1][i];
+ unsigned int i = 0;
+ while (i < _edges.size()) {
+ if (_edges[i] == re) {
+ _edges[i] = _edges.back();
+ _edges.pop_back();
+ return;
+ }
+ ++i;
}
+ Msg::Warning("[Rec2DVertex] Didn't removed edge, didn't have it");
}
-void Rec2DVertex::remove(Rec2DEdge *re)
+void Rec2DVertex::add(Rec2DElement *rel)
{
- std::vector<Rec2DEdge*>::iterator it = _edges.begin();
- bool b = false;
- for (; it != _edges.end() ; ++it)
- if (*it == re) {
- _edges.erase(it--);
- b = true;
- break;
+ for (unsigned int i = 0; i < _elements.size(); ++i) {
+ if (_elements[i] == rel) {
+ Msg::Warning("[Rec2DVertex] Element was already there");
+ return;
}
- if (b) {
- Recombine2D::addValVert(getGain(-1));
- --_numEl;
}
- else
- Msg::Error("[Rec2DVertex] Edge not found");
+ if (_isMesh)
+ Rec2DData::addValVert(getGain(1));
+ _elements.push_back(rel);
}
-void Rec2DVertex::getElements(std::set<Rec2DElement*> &elem)
+bool Rec2DVertex::has(Rec2DElement *rel) const
{
- for (unsigned int i = 0; i < _elements.size(); ++i)
- elem.insert(_elements[i]);
+ for (unsigned int i = 0; i < _elements.size(); ++i) {
+ if (_elements[i] == rel)
+ return true;
+ }
+ return false;
}
+void Rec2DVertex::remove(Rec2DElement *rel)
+{
+ unsigned int i = 0;
+ while (i < _elements.size()) {
+ if (_elements[i] == rel) {
+ _elements[i] = _elements.back();
+ _elements.pop_back();
+ return;
+ }
+ ++i;
+ }
+ Msg::Warning("[Rec2DVertex] Didn't removed element, didn't have it");
+}
-/** Rec2DEdge **/
-/*****************/
-Rec2DEdge::Rec2DEdge(MEdge e, mapofVertices &vert, std::vector<MTriangle*> &tri)
-: _lastChange(-1), _qual(.0), _weight(REC2D_EDGE_BASE),
- _rv1(vert[e.getVertex(0)]), _rv2(vert[e.getVertex(1)])
+
+/** Rec2DElement **/
+/********************/
+Rec2DElement::Rec2DElement(MTriangle *t) : _mEl((MElement *)t), _numEdge(3)
{
- _rv1->add(this);
- _rv2->add(this);
- if (tri.size() == 2) {
- Rec2DElement *t0 = Recombine2D::getRElement(tri[0]);
- Rec2DElement *t1 = Recombine2D::getRElement(tri[1]);
- t0->add(this, t1);
- t1->add(this, t0);
+ for (int i = 0; i < 4; ++i) {
+ _edges[i] = NULL;
+ _elements[i] = NULL;
}
- else if (tri.size() == 1)
- Recombine2D::getRElement(tri[0])->add(this);
- else
- Msg::Error("[Rec2DEdge] Wrong number of triangles");
+ Rec2DData::add(this);
}
-Rec2DEdge::Rec2DEdge(Rec2DVertex *rv1, Rec2DVertex *rv2)
-: _lastChange(-1), _qual(.0), _weight(REC2D_EDGE_BASE), _rv1(rv1), _rv2(rv2)
+Rec2DElement::Rec2DElement(MQuadrangle *q) : _mEl((MElement *)q), _numEdge(4)
{
+ for (int i = 0; i < 4; ++i) {
+ _edges[i] = NULL;
+ _elements[i] = NULL;
+ }
+ Rec2DData::add(this);
}
-double Rec2DEdge::getQual()
+Rec2DElement::Rec2DElement(Rec2DEdge **edges) : _mEl(NULL), _numEdge(4)
{
- if (_lastChange < Recombine2D::getNumChange())
- _computeQual();
- return _qual;
+ for (int i = 0; i < 4; ++i) {
+ _edges[i] = edges[i];
+ _elements[i] = Rec2DEdge::getSingleElement(edges[i]);
+ if (_elements[i])
+ _elements[i]->addNeighbour(_edges[i], this);
+ }
+ std::vector<Rec2DVertex*> vertices(4);
+ getVertices(vertices);
+ for (unsigned int i = 0; i < 4; ++i) {
+ vertices[i]->add(this);
+ }
+ Rec2DData::add(this);
}
-void Rec2DEdge::_computeQual() //*
+Rec2DElement::~Rec2DElement()
{
- //if (!_rv1->update() && !_rv2->update())
- // return;
- static int i = 0; if (++i < 2) Msg::Warning("FIXME compute qual edge and update of vertex");
- //_rv1->update();
- //_rv2->update();
- double adimLength = _straightAdimLength();
- double alignment = _straightAlignment();
- if (adimLength > 1)
- adimLength = 1/adimLength;
- _qual = adimLength * ((1-REC2D_ALIGNMENT) + REC2D_ALIGNMENT * alignment);
- _lastChange = Recombine2D::getNumChange();
+ if (_actions.size())
+ Msg::Error("[Rec2DElement] I didn't want to be deleted :'(");
+ for (int i = 0; i < _numEdge; ++i) {
+ if (_elements[i])
+ _elements[i]->removeNeighbour(_edges[i], this);
+ }
+ std::vector<Rec2DVertex*> vertices(_numEdge);
+ getVertices(vertices);
+ for (unsigned int i = 0; i < _numEdge; ++i) {
+ vertices[i]->remove(this);
+ }
+ Rec2DData::remove(this);
}
-double Rec2DEdge::_straightAdimLength()
+void Rec2DElement::add(Rec2DEdge *re)
{
- double dx, dy, dz, **xyz;
- xyz = new double*[2];
- xyz[0] = new double[3];
- xyz[1] = new double[3];
- _rv1->getxyz(xyz[0]);
- _rv2->getxyz(xyz[1]);
- dx = xyz[0][0] - xyz[1][0];
- dy = xyz[0][1] - xyz[1][1];
- dz = xyz[0][2] - xyz[1][2];
- double length = sqrt(dx * dx + dy * dy + dz * dz);
- delete[] xyz[0];
- delete[] xyz[1];
- delete[] xyz;
-
- double lc1 = (*Recombine2D::bgm())(_rv1->u(), _rv1->v(), .0);
- double lc2 = (*Recombine2D::bgm())(_rv2->u(), _rv2->v(), .0);
+ int i;
+ for (i = 0; i < _numEdge; ++i) {
+ if (_edges[i] == re) {
+ Msg::Warning("[Rec2DElement] Edge was already there");
+ return;
+ }
+ if (_edges[i] == NULL) {
+ _edges[i] = re;
+ break;
+ }
+ }
+ if (i == _numEdge)
+ Msg::Error("[Rec2DElement] Already %d edges, can't add anymore", _numEdge);
- return length * (1/lc1 + 1/lc2) / 2;
+ if (_numEdge == 4)
+ re->addQuadNeighbour();
+ re->addNeighbour();
}
-double Rec2DEdge::_straightAlignment()
+bool Rec2DElement::has(Rec2DEdge *re) const
{
- double angle1 = Recombine2D::bgm()->getAngle(_rv1->u(), _rv1->v(), .0);
- double angle2 = Recombine2D::bgm()->getAngle(_rv1->u(), _rv1->v(), .0);
- double angleEdge = atan2(_rv1->u()-_rv2->u(), _rv1->v()-_rv2->v());
-
- double alpha1 = angleEdge - angle1;
- double alpha2 = angleEdge - angle2;
- crossField2d::normalizeAngle(alpha1);
- crossField2d::normalizeAngle(alpha2);
- alpha1 = 1 - 4 * std::min(alpha1, .5 * M_PI - alpha1) / M_PI;
- alpha2 = 1 - 4 * std::min(alpha2, .5 * M_PI - alpha2) / M_PI;
-
- double lc1 = (*Recombine2D::bgm())(_rv1->u(), _rv1->v(), .0);
- double lc2 = (*Recombine2D::bgm())(_rv1->u(), _rv1->v(), .0);
-
- return (alpha1/lc1 + alpha2/lc2) / (1/lc1 + 1/lc2);
+ for (unsigned int i = 0; i < _numEdge; ++i) {
+ if (_edges[i] == re)
+ return true;
+ }
+ return false;
}
+void Rec2DElement::add(Rec2DAction *ra)
+{
+ for (unsigned int i = 0; i < _actions.size(); ++i) {
+ if (_actions[i] == ra) {
+ Msg::Warning("[Rec2DElement] Action was already there");
+ return;
+ }
+ }
+ _actions.push_back(ra);
+}
-/** Rec2DElement **/
-/********************/
-void Rec2DElement::add(Rec2DEdge *re)
+void Rec2DElement::removeT(Rec2DAction *ra)
+{
+ unsigned int i = 0;
+ while (i < _actions.size()) {
+ if (_actions[i] == ra) {
+ _actions[i] = _actions.back();
+ _actions.pop_back();
+ return;
+ }
+ ++i;
+ }
+ Msg::Warning("[Rec2DElement] Didn't removed action, didn't have it");
+}
+
+void Rec2DElement::addNeighbour(Rec2DEdge *re, Rec2DElement *rel)
+{
+ for (int i = 0; i < _numEdge; ++i) {
+ if (_edges[i] == re) {
+ if (_elements[i] != NULL && _elements[i] != rel)
+ Msg::Error("[Rec2DElement] Have already a neighbour element");
+ _elements[i] = rel;
+ return;
+ }
+ }
+ Msg::Error("[Rec2DElement] Edge not found");
+}
+
+
+void Rec2DElement::removeNeighbour(Rec2DEdge *re, Rec2DElement *rel)
+{
+ for (int i = 0; i < _numEdge; ++i) {
+ if (_edges[i] == re) {
+ if (_elements[i] == rel)
+ _elements[i] = NULL;
+ else
+ Msg::Error("[Rec2DElement] I didn't know this element was my neighbour...");
+ return;
+ }
+ }
+ Msg::Error("[Rec2DElement] Edge not found");
+}
+
+
+void Rec2DElement::getMoreEdges(std::vector<Rec2DEdge*> &edges) const
+{
+ for (int i = 0; i < _numEdge; ++i)
+ edges.push_back(_edges[i]);
+}
+
+void Rec2DElement::getVertices(std::vector<Rec2DVertex*> &verts) const
+{
+ verts.resize(_numEdge);
+ int i = 0, k = 0;
+ while (k < _numEdge) {
+ verts[k] = _edges[i/2]->getVertex(i%2);
+ bool b = true;
+ for (unsigned int j = 0; j < k; ++j) {
+ if (verts[k] == verts[j])
+ b = false;
+ }
+ if (b) ++k;
+ ++i;
+ }
+}
+
+void Rec2DElement::getUniqueActions(std::vector<Rec2DAction*> &vectRA) const
+{
+ for (unsigned int i = 0; i < _actions.size(); ++i) {
+ unsigned int j = -1;
+ while (++j < vectRA.size() && vectRA[j] != _actions[i]);
+ if (j == vectRA.size())
+ vectRA.push_back(_actions[i]);
+ }
+}
+
+Rec2DEdge* Rec2DElement::getCommonEdge(Rec2DElement *rel0, Rec2DElement *rel1)
+{
+ for (unsigned int i = 0; i < rel0->_numEdge; ++i) {
+ if (rel1->has(rel0->_edges[i]))
+ return rel0->_edges[i];
+ }
+ Msg::Error("[Rec2DElement] didn't find edge, returning NULL");
+ return NULL;
+}
+
+Rec2DVertex* Rec2DElement::getOtherVertex(Rec2DVertex *rv1, Rec2DVertex *rv2) const
{
- int i = -1;
- while (_redges[++i] && i < 4);
- if (_redges[i])
- Msg::Error("[Rec2DElement] too much edges");
- _redges[i] = re;
+ if (_numEdge == 4) {
+ Msg::Error("[Rec2DElement] I'm not a triangle");
+ return NULL;
+ }
+ for (int i = 0; i < 2; ++i) {
+ Rec2DVertex *rv = _edges[i]->getVertex(0);
+ if (rv != rv1 && rv != rv2)
+ return rv;
+ rv = _edges[i]->getVertex(1);
+ if (rv != rv1 && rv != rv2)
+ return rv;
+ }
+ Msg::Error("[Rec2DElement] I should not be here... Why this happen to me ?");
+ return NULL;
}
-void Rec2DElement::add(Rec2DEdge *re, Rec2DElement *el)
+MQuadrangle* Rec2DElement::_createQuad() const
{
- int i = -1;
- while (_redges[++i] && i < 4);
- if (_redges[i])
- Msg::Error("[Rec2DElement] too much edges");
- _redges[i] = re;
- _relements[i] = el;
+ if (_numEdge != 4) {
+ Msg::Error("[Rec2DElement] Why do you ask me to do this ? You know I can't do it ! COULD YOU LEAVE ME KNOW ?");
+ return new MQuadrangle(NULL, NULL, NULL, NULL);
+ }
+ MVertex *v1, *v2, *v3, *v4;
+ v1 = _edges[0]->getVertex(0)->getMVertex();
+ v2 = _edges[0]->getVertex(1)->getMVertex();
+ int I;
+ for (unsigned int i = 1; i < 4; ++i) {
+ if (v2 == _edges[i]->getVertex(0)->getMVertex()) {
+ v3 = _edges[i]->getVertex(1)->getMVertex();
+ I = i;
+ break;
+ }
+ if (v2 == _edges[i]->getVertex(1)->getMVertex()) {
+ v3 = _edges[i]->getVertex(0)->getMVertex();
+ I = i;
+ break;
+ }
+ }
+ for (unsigned int i = 1; i < 4; ++i) {
+ if (i == I) continue;
+ if (v3 == _edges[i]->getVertex(0)->getMVertex()) {
+ v4 = _edges[i]->getVertex(1)->getMVertex();
+ break;
+ }
+ if (v3 == _edges[i]->getVertex(1)->getMVertex()) {
+ v4 = _edges[i]->getVertex(0)->getMVertex();
+ break;
+ }
+ }
+ return new MQuadrangle(v1, v2, v3, v4);
}
diff --git a/Mesh/meshGFaceRecombine.h b/Mesh/meshGFaceRecombine.h
index ebd9f872adcd0d1f8aa077ac9f9fece409bba74a..f786283249f6fd2d452b13f1c501144e97cce1bd 100644
--- a/Mesh/meshGFaceRecombine.h
+++ b/Mesh/meshGFaceRecombine.h
@@ -1,4 +1,4 @@
-// Gmsh - Copyright (C) 1997-2011 C. Geuzaine, J.-F. Remacle
+// Gmsh - Copyright (C) 1997-2012 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to <gmsh@geuz.org>.
@@ -11,414 +11,306 @@
#define _MESH_GFACE_RECOMBINE_H_
#include "GFace.h"
-#include "GModel.h"
-#include "MEdge.h"
#include "BackgroundMesh.h"
+//#include "GModel.h"
+//#include "MEdge.h"
+
+#define REC2D_EDGE_BASE 1
+#define REC2D_EDGE_QUAD 1
+#define REC2D_ALIGNMENT .5
+#define REC2D_WAIT_TIME .05
+#define REC2D_NUM_ACTIO 1000
+// #define REC2D_SMOOTH
+ #define REC2D_DRAW
class Rec2DVertex;
class Rec2DEdge;
class Rec2DElement;
class Rec2DAction;
+class Rec2DData;
-typedef std::list<Rec2DAction*> setofRec2DAction;
-typedef std::map<MVertex*, Rec2DVertex*> mapofVertices;
-typedef std::map<MEdge, Rec2DEdge*, Less_Edge> mapofEdges;
-typedef std::map<MElement*, Rec2DElement*> mapofElements;
-typedef std::map<MElement*, std::set<Rec2DAction*> > mapofElementActions;
-typedef std::map<MQuadrangle*, std::set<MElement*> > mapofAdjacencies;
-
-struct lessRec2DAction {
- bool operator()(Rec2DAction*, Rec2DAction*) const;
-};
+//typedef std::list<Rec2DAction*> setofRec2DAction;
+//typedef std::map<MVertex*, Rec2DVertex*> mapofVertices;
+//typedef std::map<MEdge, Rec2DEdge*, Less_Edge> mapofEdges;
+//typedef std::map<MElement*, Rec2DElement*> mapofElements;
+//typedef std::map<MElement*, std::set<Rec2DAction*> > mapofElementActions;
+//typedef std::map<MQuadrangle*, std::set<MElement*> > mapofAdjacencies;
class Recombine2D {
private :
- int _numEdge, _numVert, _numChange;
- double _valEdge, _valVert;
GFace *_gf;
backgroundMesh *_bgm;
- mapofVertices _vertices;
- mapofEdges _edges;
- mapofElements _elements;
- setofRec2DAction _actions;
- mapofElementActions _mea;
- mapofAdjacencies _ma;
- std::map<int, std::set<Rec2DVertex*> > _parities;
- std::vector<Rec2DAction*> _obsoleteAction;
+ Rec2DData *_data;
static Recombine2D *_current;
- std::vector<MTriangle*> _tri;
- std::vector<MQuadrangle*> _quad;
- std::vector<Rec2DVertex*> _checkAllQuad;
+ int _numChange;
public :
- static bool bo;
-
Recombine2D(GFace*);
- Recombine2D(GModel*);
~Recombine2D();
bool recombine();
- void apply(bool){return;}
-
- static int getNumVertAllQuad() {return _current->_checkAllQuad.size();}
- static Rec2DVertex* getVertAllQuad(int i) {return _current->_checkAllQuad[i];}
- static void removeVertAllQuad(Rec2DVertex *rv) {
- std::vector<Rec2DVertex*>::iterator it = _current->_checkAllQuad.begin();
- while (it != _current->_checkAllQuad.end()) {
- if (*it == rv)
- it = _current->_checkAllQuad.erase(it);
- else
- ++it;
- }
- }
- static void addVertAllQuad(Rec2DVertex *rv) {
- _current->_checkAllQuad.push_back(rv);
- }
- //static inline double edgeReward(double adimLength, double alignment)
- // {return adimLength * ((1-REC2D_ALIGNMENT) + REC2D_ALIGNMENT * alignment);}
+ static inline GFace* getGFace() {return _current->_gf;}
static inline int getNumChange() {return _current->_numChange;}
static inline backgroundMesh* bgm() {return _current->_bgm;}
+ //recombine : if _current==this ok !
+
+ private :
+ double _geomAngle(MVertex*,
+ std::vector<GEdge*>&,
+ std::vector<MElement*>&);
+};
+
+class Rec2DData {
+ private :
+ int _numEdge, _numVert;
+ double _valEdge, _valVert;
+ static Rec2DData *_current;
+
+ std::set<Rec2DEdge*> _edges;
+ std::set<Rec2DVertex*> _vertices;
+ std::set<Rec2DElement*> _elements;
+
+ std::list<Rec2DAction*> _actions;
+ std::map<int, std::vector<Rec2DVertex*> > _parities;
+
+ public :
+ Rec2DData(int numTri, int numQuad);
+ ~Rec2DData();
+
+ void printState();
+#ifdef REC2D_DRAW
+ std::vector<MTriangle*> _tri;
+ std::vector<MQuadrangle*> _quad;
+#endif
+
+ static double getGlobalValue();
static double getGlobalValue(int numEdge, double valEdge,
int numVert, double valVert );
- //static Rec2DVertex* getRVertex(MVertex*);
- static Rec2DEdge* getREdge(MEdge);
- static Rec2DElement* getRElement(MElement*);
- static void remove(MEdge);
- static void remove(MVertex*);
- static inline GFace* getGFace() {return _current->_gf;}
- static inline void addObsleteAction(Rec2DAction *a) {
- _current->_obsoleteAction.push_back(a);
- } // correct ?
- static inline void addNumEdge(int a) {_current->_numEdge += a;}
- static inline void addNumVert(int a) {_current->_numVert += a;}
- static inline void addValEdge(double a) {_current->_valEdge += a;}
- static inline void addValVert(double a) {_current->_valVert += a;}
+ static Rec2DAction* getBestAction();
+
+ typedef std::set<Rec2DEdge*>::iterator iter_re;
+ typedef std::set<Rec2DVertex*>::iterator iter_rv;
+ typedef std::set<Rec2DElement*>::iterator iter_rel;
+ static inline iter_re firstEdge() {return _current->_edges.begin();}
+ static inline iter_rv firstVertex() {return _current->_vertices.begin();}
+ static inline iter_rel firstElement() {return _current->_elements.begin();}
+ static inline iter_re lastEdge() {return _current->_edges.end();}
+ static inline iter_rv lastVertex() {return _current->_vertices.end();}
+ static inline iter_rel lastElement() {return _current->_elements.end();}
+
+ static inline void add(Rec2DEdge *re) {_current->_edges.insert(re);}
+ static inline void add(Rec2DVertex *rv) {_current->_vertices.insert(rv);}
+#ifndef REC2D_DRAW
+ static inline void add(Rec2DElement *rel) {_current->_elements.insert(rel);}
+#else
+ static void add(Rec2DElement *rel);
+#endif
+ static inline void add(Rec2DAction *ra) {_current->_actions.push_back(ra);}
+ static void remove(Rec2DEdge*);
+ static void remove(Rec2DVertex*);
static void remove(Rec2DElement*);
static void remove(Rec2DAction*);
- static void add(MQuadrangle*);
-
- static inline void addVert2Par(Rec2DVertex *rv, int a)
- {_current->_parities[a].insert(rv);}
- static void associateParity(int pOld, int pNew);
- static inline int highPar() {
- std::map<int, std::set<Rec2DVertex*> >::iterator it;
- it = _current->_parities.end();
- --it;
- return it->first + 1;
- }
- private :
- bool _isAllQuad(Rec2DAction *action);
- void _choose(Rec2DAction*);
- void _removeIncompatible(Rec2DAction*);
- void _removeReference(Rec2DAction*);
- void _getNeighbours(Rec2DAction*, std::set<Rec2DAction*>&);
- void _getElemToLook(int par, std::set<Rec2DElement*>&);
+ static int getNewParity();
+ static void removeParity(Rec2DVertex*, int);
+ static inline void addParity(Rec2DVertex *rv, int p)
+ {_current->_parities[p].push_back(rv);}
+
+ static inline void addVert(int num, double val) {
+ _current->_numVert += num;
+ _current->_valVert += val;
+ }
+ static inline void addValVert(double val) {_current->_valVert += val;}
+ static inline void addEdge(int num, double val) {
+ _current->_numEdge += num;
+ _current->_valEdge += val;
+ }
+};
+
+struct lessRec2DAction {
+ bool operator()(Rec2DAction*, Rec2DAction*) const;
};
class Rec2DAction {
- private :
- int _position;
-
protected :
- int _lastUpdate; // + _lastCheck
double _globValIfExecuted;
+ int _lastUpdate;
public :
- virtual ~Rec2DAction() {Recombine2D::remove(this);}
+ virtual inline ~Rec2DAction() {Rec2DData::remove(this);}
bool operator<(Rec2DAction&);
- inline int getVectPosition() const {return _position;}
- inline void setVectPosition(int p) {_position = p;}
-
- virtual int getNum() = 0;
- virtual double getReward() {
- if (_lastUpdate < Recombine2D::getNumChange())
- _computeGlobVal();
- return _globValIfExecuted - Recombine2D::getGlobalValue(0, .0, 0, .0);
- }
- virtual bool isObsolete() = 0;
- virtual void getRTriangles(std::set<Rec2DElement*>&) = 0;
- virtual Rec2DElement* getRTriangle(int) = 0;
- virtual void apply() = 0;
- virtual void getParities(int*) = 0;
- virtual Rec2DVertex* getRVertex(int) = 0;
+ virtual double getReward();
virtual void color(int, int, int) = 0;
- virtual bool whatImply(Rec2DVertex**, int*,
- std::map<Rec2DVertex*, std::vector<int> >&) = 0;
-
- protected :
- MQuadrangle* _createQuad(std::vector<Rec2DEdge*>&) const;
-
+ virtual void apply() = 0;
+
private :
virtual void _computeGlobVal() = 0;
};
class Rec2DTwoTri2Quad : public Rec2DAction {
private :
- Rec2DElement *_rtriangles[2];
- Rec2DEdge *_edges[5]; // 1 embedded, 4 boundary
+ Rec2DElement *_triangles[2];
+ Rec2DEdge *_edges[5]; // 4 boundary, 1 embedded
Rec2DVertex *_vertices[4]; // 4 boundary (2 on embedded edge + 2)
public :
- Rec2DTwoTri2Quad(Rec2DEdge*, std::vector<MTriangle*>&);
+ Rec2DTwoTri2Quad(Rec2DElement*, Rec2DElement*);
~Rec2DTwoTri2Quad();
- bool isObsolete();
- void getRTriangles(std::set<Rec2DElement*> &tri) {
- tri.insert(_rtriangles[0]);
- tri.insert(_rtriangles[1]);
- }
- Rec2DElement* getRTriangle(int a) {
- return _rtriangles[a];
- }
- void apply();
- void getParities(int*);
- Rec2DVertex* getRVertex(int i) {return _vertices[i];}
- void color(int, int, int);
- virtual bool whatImply(Rec2DVertex**, int*,
- std::map<Rec2DVertex*, std::vector<int> >&);
- static bool isObsolete(int*);
-
- int getNum() {return 1;}
- private :
- void _computeGlobVal();
-};
-
-/*class Rec2DFourTri2Quad : public Rec2DAction {
- private :
- Rec2DElement *_rtriangles[4];
- Rec2DEdge *_edges[8]; // 4 embedded, 4 boundary
- Rec2DVertex *_vertices[5]; // 4 boundary (2 non adjacent + 2), 1 embedded
- public :
- Rec2DFourTri2Quad(Rec2DVertex*, std::vector<MTriangle*>&);
- ~Rec2DFourTri2Quad();
- bool isObsolete();
- void getRTriangles(std::set<Rec2DElement*> &tri) {
- tri.insert(_rtriangles[0]);
- tri.insert(_rtriangles[1]);
- tri.insert(_rtriangles[2]);
- tri.insert(_rtriangles[3]);
- }
- Rec2DElement* getRTriangle(int a) {
- return _rtriangles[a];
- }
- void apply();
- void getParities(int*);
- Rec2DVertex* getRVertex(int i) {return _vertices[i];}
- void color(int, int, int);
- virtual void whatImply(Rec2DVertex**, int*,
- std::map<Rec2DVertex*, std::vector<int> >&)
- {return;}
-
- int getNum() {return 2;}
+ virtual void color(int, int, int);
+ virtual void apply();
+
private :
- void _computeGlobVal();
+ virtual void _computeGlobVal();
};
-class Rec2DCollapseTri : public Rec2DAction {
+class Rec2DEdge {
private :
- Rec2DElement *_rtriangles[4]; // either 2 or 4 triangles
- Rec2DEdge *_edge; // 1 embedded or NULL
- Rec2DVertex *_vertices[5]; // 4 boundary (2 on embedded edge + 2), 1 embedded or NULL
- double _globValIfExecuted2;
+ Rec2DVertex *_rv0, *_rv1;
+ double _qual;
+ int _lastUpdate, _weight;
+ int _boundary;
public :
- Rec2DCollapseTri(Rec2DVertex*, std::vector<MTriangle*>&);
- Rec2DCollapseTri(Rec2DEdge*, std::vector<MTriangle*>&);
- ~Rec2DCollapseTri();
- double getReward() {
- return std::max(_globValIfExecuted, _globValIfExecuted2)
- - Recombine2D::getGlobalValue(0, .0, 0, .0);
- }
- bool isObsolete();
- void getRTriangles(std::set<Rec2DElement*> &tri) {
- tri.insert(_rtriangles[0]);
- tri.insert(_rtriangles[1]);
- if (_rtriangles[2]) {
- tri.insert(_rtriangles[2]);
- tri.insert(_rtriangles[3]);
- }
- }
- Rec2DElement* getRTriangle(int a) {
- return _rtriangles[a];
- }
- void apply();
- void getParities(int*);
- Rec2DVertex* getRVertex(int i) {return _vertices[i];}
- void color(int, int, int);
- virtual void whatImply(Rec2DVertex**, int*,
- std::map<Rec2DVertex*, std::vector<int> >&)
- {return;}
-
- int getNum() {return 3;}
- private :
- void _computeGlobVal();
- void _qualCavity(double &valVert, double &valEdge, int &numEdge,
- std::map<Rec2DVertex*, int> &vertices,
- Rec2DVertex *imposedVert = NULL );
-};*/
-//
-
-class Rec2DListAction {
- private :
+ Rec2DEdge(Rec2DVertex*, Rec2DVertex*);
+ ~Rec2DEdge();
+ double getQual();
+ double addNeighbour();
+ double addQuadNeighbour();
- public :
+ inline bool isOnBoundary() const {return !_boundary;}
+ inline Rec2DVertex* getVertex(int i) const {if (i) return _rv1; return _rv0;}
+
+ static Rec2DElement* getSingleElement(Rec2DEdge*);
+
+ void swap(Rec2DVertex *oldRV, Rec2DVertex *newRV);
+ Rec2DVertex* getOtherVertex(Rec2DVertex*) const;
+
+ private :
+ void _computeQual();
+ double _straightAdimLength() const;
+ double _straightAlignment() const;
+};
+
+struct AngleData {
+ std::vector<GEdge*> _gEdges;
+ std::vector<MElement*> _mElements;
+ Rec2DVertex *_rv;
+
+ AngleData() : _rv(NULL) {}
};
class Rec2DVertex {
private :
MVertex *_v;
+ const double _angle;
std::vector<Rec2DEdge*> _edges;
std::vector<Rec2DElement*> _elements;
- int _onWhat, _numEl, _parity, _lastChange; //lastMove
+ int _onWhat, _parity, _lastMove;
// _onWhat={-1:corner,0:edge,1:face,(2:volume)}
- double _angle, _qual;
SPoint2 _param;
+ bool _isMesh;
static double **_qualVSnum;
static double **_gains;
public :
- Rec2DVertex(MVertex*, std::vector<MTriangle*>, int onWhat,
- std::map<MVertex*, std::set<GEdge*> >);
- Rec2DVertex(double u, double v);
-
- bool hasTriangle();
- inline void add(Rec2DEdge *re) {_edges.push_back(re);}
- void remove(Rec2DEdge*);
- inline void setOnBoundary() {if (_onWhat > 0) _onWhat = 0;
- static bool a=1; if(a){Msg::Warning("FIXME boundary");a=0;}}
- inline bool getIsOnBoundary() {return _onWhat < 1;}
- double getQual(int numEdge = -1);
- inline int getParity() {return _parity;}
- inline void setParity(int p) {
- _parity = p;
- Recombine2D::addVert2Par(this, p);
- }
- inline MVertex* getMVertex() {return _v;}
- void getxyz(double *xyz) {
- xyz[0] = _v->x(); xyz[1] = _v->y(); xyz[2] = _v->z();}
- inline double u() {return _param[0];}
- inline double v() {return _param[1];}
- inline void deleteMVertex() {delete _v;}
-
- double getGain(int);
- void getNeighbours(std::map<Rec2DVertex*, int>&, std::set<Rec2DEdge*>&);
- void getNeighbourEdges(std::set<Rec2DEdge*>&);
- void getElements(std::set<Rec2DElement*>&);
-
- void remove(Rec2DElement *rel) {
- std::vector<Rec2DElement*>::iterator it = _elements.begin();
- while (it != _elements.end()) {
- if (*it == rel)
- it = _elements.erase(it);
- else
- ++it;
- }
- }
+ Rec2DVertex(MVertex*, bool toSave = true);
+ Rec2DVertex(Rec2DVertex*, double angle);
+ ~Rec2DVertex();
static void initStaticTable();
+ static Rec2DEdge* getCommonEdge(Rec2DVertex*, Rec2DVertex*);
+ static void getCommonElements(Rec2DVertex*, Rec2DVertex*,
+ std::vector<Rec2DElement*>&);
+
+ inline void setOnBoundary() {if (_onWhat > 0) _onWhat = 0;}
+ inline bool getOnBoundary() const {return _onWhat < 1;}
+ bool setBoundaryParity(int p0, int p1);
+ inline int getParity() const {return _parity;}
+ void setParity(int);
+ inline int getNumElements() const {return _elements.size();}
+ inline MVertex* getMVertex() const {return _v;}
+
+ inline int getLastMove() const {return _lastMove;}
+ inline void getxyz(double *xyz) const {
+ xyz[0] = _v->x();
+ xyz[1] = _v->y();
+ xyz[2] = _v->z();
+ }
+ inline double u() const {return _param[0];}
+ inline double v() const {return _param[1];}
- private :
- void _computeQual();
- double _computeAngle(MVertex*, std::vector<MTriangle*>&, std::set<GEdge*>&);
-};
-
-class Rec2DEdge {
- private :
- Rec2DVertex *_rv1, *_rv2;
- double _qual;
- int _lastChange, _weight; //lastUpdate
+ double getQual(int numEl = -1) const;
+ double getGain(int) const;
- public :
- Rec2DEdge(MEdge, mapofVertices&, std::vector<MTriangle*>&);
- Rec2DEdge(Rec2DVertex*, Rec2DVertex*);
+ void add(Rec2DEdge*);
+ bool has(Rec2DEdge*) const;
+ void remove(Rec2DEdge*);
- double getQual();
- inline MEdge getMEdge() {return MEdge(_rv1->getMVertex(), _rv2->getMVertex());}
- inline Rec2DVertex* getRVertex(int i) {if (i) return _rv2; return _rv1;}
- int getWeight() {return _weight;}
- double getWeightedQual() {return _weight * _qual;}
+ void add(Rec2DElement*);
+ bool has(Rec2DElement*) const;
+ void remove(Rec2DElement*);
private :
- void _computeQual();
- double _straightAdimLength();
- double _straightAlignment();
+ bool _recursiveBoundParity(Rec2DVertex *prev, int p0, int p1);
};
class Rec2DElement {
private :
- Rec2DEdge *_redges[4]; // if triangle, _redges[3] == NULL
- Rec2DElement *_relements[4]; // NULL if no neighbour
- MElement *_mEl; // can be NULL
+ int _numEdge;
+ Rec2DEdge *_edges[4];
+ Rec2DElement *_elements[4]; // NULL if no neighbour
+ MElement *_mEl; // can be NULL
std::vector<Rec2DAction*> _actions;
public :
- Rec2DElement(MTriangle *t) : _mEl((MElement *)t) {
- for (int i = 0; i < 4; ++i) {
- _redges[i] = NULL;
- }
- }
+ Rec2DElement(MTriangle*);
+ Rec2DElement(MQuadrangle*);
+ Rec2DElement(Rec2DEdge**);
+ ~Rec2DElement();
void add(Rec2DEdge*);
- void add(Rec2DEdge*, Rec2DElement*);
- inline MElement* getMElement() {return _mEl;}
- inline MTriangle* getMTriangle()
- {if (_redges[3])return NULL; return (MTriangle*) _mEl;}
- inline bool isQuad() {return _redges[3];}
- inline void add(Rec2DAction *ra) {_actions.push_back(ra);}
- void remove(Rec2DAction *ra) {
- std::vector<Rec2DAction*>::iterator it = _actions.begin();
- while (it != _actions.end()) {
- if (*it == ra) it = _actions.erase(it);
- else ++it;
- }
- }
- inline int getNum() {return _mEl->getNum();}
- inline int getNumActions() {return _actions.size();}
- inline Rec2DAction* getAction(int i) {return _actions[i];}
- void getVertices(Rec2DVertex **v) {
- v[0] = _redges[0]->getRVertex(0);
- v[1] = _redges[0]->getRVertex(1);
- if (_redges[1]->getRVertex(0) == v[0] ||
- _redges[1]->getRVertex(0) == v[1] )
- v[2] = _redges[1]->getRVertex(1);
- else
- v[2] = _redges[1]->getRVertex(0);
- }
- void getActions(std::set<Rec2DAction*> &actions) {
- for (unsigned int i = 0; i < _actions.size(); ++i)
- actions.insert(_actions[i]);
- }
- void getNeighbours(std::set<Rec2DElement*> &elem) {
- for (int i = 0; i < 4; ++i) {
- if (_relements[i])
- elem.insert(_relements[i]);
+ bool has(Rec2DEdge*) const;
+ void add(Rec2DAction*);
+ void removeT(Rec2DAction*);
+ void addNeighbour(Rec2DEdge*, Rec2DElement*);
+ void removeNeighbour(Rec2DEdge*, Rec2DElement*);
+ inline MElement* getMElement() const {return _mEl;}
+#ifdef REC2D_DRAW
+ MTriangle* getMTriangle() {
+ if (_numEdge == 3) {
+ if (_mEl)
+ return (MTriangle*) _mEl;
+ else
+ Msg::Error("[Rec2DElement] Do you thing I'll create a triangle for you ?");
}
+ return NULL;
}
- Rec2DVertex* getOtherVertex(Rec2DVertex *rv1, Rec2DVertex *rv2) {
- if (_redges[3]) {
- Msg::Error("[Rec2DElement] I'm not a triangle");
- return NULL;
+ MQuadrangle* getMQuadrangle() {
+ if (_numEdge == 4) {
+ if (!_mEl)
+ _mEl = (MElement*) _createQuad();
+ return (MQuadrangle*) _mEl;
}
- std::set<Rec2DVertex*> vert;
- for (int i = 0; i < 3; ++i) {
- vert.insert(_redges[i]->getRVertex(0));
- vert.insert(_redges[i]->getRVertex(1));
- }
- std::set<Rec2DVertex*>::iterator it = vert.begin();
- for (; it != vert.end(); ++it) {
- if (*it != rv1 && *it != rv2)
- return *it;
- }
- Msg::Error("[Rec2DElement] I should not be here... Why this happen to me ?");
return NULL;
}
+#endif
+
+ void getMoreEdges(std::vector<Rec2DEdge*>&) const;
+ void getVertices(std::vector<Rec2DVertex*>&) const;
+ void getUniqueActions(std::vector<Rec2DAction*>&) const;
+ static Rec2DEdge* getCommonEdge(Rec2DElement*, Rec2DElement*);
+
+ Rec2DVertex* getOtherVertex(Rec2DVertex*, Rec2DVertex*) const;
+
+ private :
+ MQuadrangle* _createQuad() const;
};
#endif
\ No newline at end of file