diff --git a/Geo/Cell.cpp b/Geo/Cell.cpp
index 60341b93d21bac14cf7ac144555ff6d96acf07b1..0b30b209c2ef9203dca71c53462a2d3122fa30de 100755
--- a/Geo/Cell.cpp
+++ b/Geo/Cell.cpp
@@ -32,7 +32,7 @@ bool Less_Cell::operator()(const Cell* c1, const Cell* c2) const
Cell::Cell(MElement* image) :
_combined(false), _index(0), _immune(false), _image(NULL),
- _deleteImage(false), _inSubdomain(false)
+ _delimage(false), _subdomain(false)
{
_dim = image->getDim();
_image = image;
@@ -43,11 +43,15 @@ Cell::Cell(MElement* image) :
Cell::~Cell()
{
- if(_deleteImage) delete _image;
+ if(_delimage) delete _image;
}
int Cell::getNumFacets() const
-{
+{
+ if(_image == NULL){
+ printf("ERROR: No image mesh element for cell. \n");
+ return 0;
+ }
if(getDim() == 0) return 0;
else if(getDim() == 1) return 2;
else if(getDim() == 2) return _image->getNumEdges();
@@ -57,6 +61,10 @@ int Cell::getNumFacets() const
void Cell::getFacetVertices(const int num, std::vector<MVertex*> &v) const
{
+ if(_image == NULL){
+ printf("ERROR: No image mesh element for cell. \n");
+ return;
+ }
if(getDim() == 0) return;
else if(getDim() == 1) { v.resize(1); v[0] = getVertex(num); }
else if(getDim() == 2) _image->getEdgeVertices(num, v);
@@ -66,6 +74,10 @@ void Cell::getFacetVertices(const int num, std::vector<MVertex*> &v) const
int Cell::getFacetOri(Cell* cell)
{
+ if(_image == NULL){
+ printf("ERROR: No image mesh element for cell. \n");
+ return 0;
+ }
std::vector<MVertex*> v;
for(int i = 0; i < cell->getNumVertices(); i++) {
v.push_back(cell->getVertex(i));
@@ -112,57 +124,61 @@ void Cell::printCell()
for(int i = 0; i < this->getNumVertices(); i++){
printf("%d ", this->getSortedVertex(i));
}
- printf(", in subdomain: %d\n", _inSubdomain);
+ printf(", in subdomain: %d\n", inSubdomain());
printf("Combined: %d \n" , isCombined() );
};
void Cell::restoreCell(){
- _boundary = _obd;
- _coboundary = _ocbd;
+ _bd = _obd;
+ _cbd = _ocbd;
_combined = false;
_index = 0;
_immune = false;
}
-bool Cell::addBoundaryCell(int orientation, Cell* cell, bool org)
+bool Cell::addBoundaryCell(int orientation, Cell* cell,
+ bool orig, bool other)
{
- if(org) _obd.insert( std::make_pair(cell, orientation ) );
- biter it = _boundary.find(cell);
- if(it != _boundary.end()){
+ if(orig) _obd.insert( std::make_pair(cell, orientation ) );
+ if(other) cell->addCoboundaryCell(orientation, this);
+ biter it = _bd.find(cell);
+ if(it != _bd.end()){
(*it).second = (*it).second + orientation;
if((*it).second == 0) {
- _boundary.erase(it);
+ _bd.erase(it);
(*it).first->removeCoboundaryCell(this,false);
return false;
}
return true;
}
- _boundary.insert( std::make_pair(cell, orientation ) );
+ _bd.insert( std::make_pair(cell, orientation ) );
return true;
}
-bool Cell::addCoboundaryCell(int orientation, Cell* cell, bool org)
+bool Cell::addCoboundaryCell(int orientation, Cell* cell,
+ bool orig, bool other)
{
- if(org) _ocbd.insert( std::make_pair(cell, orientation ) );
- biter it = _coboundary.find(cell);
- if(it != _coboundary.end()){
+ if(orig) _ocbd.insert( std::make_pair(cell, orientation ) );
+ if(other) cell->addBoundaryCell(orientation, this);
+ biter it = _cbd.find(cell);
+ if(it != _cbd.end()){
(*it).second = (*it).second + orientation;
if((*it).second == 0) {
- _coboundary.erase(it);
+ _cbd.erase(it);
(*it).first->removeBoundaryCell(this,false);
return false;
}
return true;
}
- _coboundary.insert( std::make_pair(cell, orientation ) );
+ _cbd.insert( std::make_pair(cell, orientation ) );
return true;
}
int Cell::removeBoundaryCell(Cell* cell, bool other)
{
- biter it = _boundary.find(cell);
- if(it != _boundary.end()){
- _boundary.erase(it);
+ biter it = _bd.find(cell);
+ if(it != _bd.end()){
+ _bd.erase(it);
if(other) (*it).first->removeCoboundaryCell(this, false);
return (*it).second;
}
@@ -172,20 +188,20 @@ int Cell::removeBoundaryCell(Cell* cell, bool other)
int Cell::removeCoboundaryCell(Cell* cell, bool other)
{
- biter it = _coboundary.find(cell);
- if(it != _coboundary.end()){
- _coboundary.erase(it);
+ biter it = _cbd.find(cell);
+ if(it != _cbd.end()){
+ _cbd.erase(it);
if(other) (*it).first->removeBoundaryCell(this, false);
return (*it).second;
}
return 0;
}
-bool Cell::hasBoundary(Cell* cell, bool org)
+bool Cell::hasBoundary(Cell* cell, bool orig)
{
- if(!org){
- biter it = _boundary.find(cell);
- if(it != _boundary.end()) return true;
+ if(!orig){
+ biter it = _bd.find(cell);
+ if(it != _bd.end()) return true;
return false;
}
else{
@@ -195,11 +211,11 @@ bool Cell::hasBoundary(Cell* cell, bool org)
}
}
-bool Cell::hasCoboundary(Cell* cell, bool org)
+bool Cell::hasCoboundary(Cell* cell, bool orig)
{
- if(!org){
- biter it = _coboundary.find(cell);
- if(it != _coboundary.end()) return true;
+ if(!orig){
+ biter it = _cbd.find(cell);
+ if(it != _cbd.end()) return true;
return false;
}
else{
@@ -209,25 +225,25 @@ bool Cell::hasCoboundary(Cell* cell, bool org)
}
}
-void Cell::printBoundary(bool org)
+void Cell::printBoundary(bool orig)
{
- for(biter it = firstBoundary(org); it != lastBoundary(org); it++){
+ for(biter it = firstBoundary(orig); it != lastBoundary(orig); it++){
printf("Boundary cell orientation: %d ", (*it).second);
Cell* cell2 = (*it).first;
cell2->printCell();
}
- if(firstBoundary(org) == lastBoundary(org)){
+ if(firstBoundary(orig) == lastBoundary(orig)){
printf("Cell boundary is empty. \n");
}
}
-void Cell::printCoboundary(bool org)
+void Cell::printCoboundary(bool orig)
{
- for(biter it = firstCoboundary(org); it != lastCoboundary(org); it++){
+ for(biter it = firstCoboundary(orig); it != lastCoboundary(orig); it++){
printf("Coboundary cell orientation: %d, ", (*it).second);
Cell* cell2 = (*it).first;
cell2->printCell();
- if(firstCoboundary(org) == lastCoboundary(org)){
+ if(firstCoboundary(orig) == lastCoboundary(orig)){
printf("Cell coboundary is empty. \n");
}
}
@@ -244,7 +260,7 @@ CombinedCell::CombinedCell(Cell* c1, Cell* c2, bool orMatch, bool co) : Cell()
_index = c1->getIndex();
_dim = c1->getDim();
- _inSubdomain = c1->inSubdomain();
+ _subdomain = c1->inSubdomain();
_combined = true;
_image = NULL;
diff --git a/Geo/Cell.h b/Geo/Cell.h
index e497d154dfa5edf30654581bf03304c7319969da..101de646bcb6cdda0299090d1e8e61ac432079ca 100644
--- a/Geo/Cell.h
+++ b/Geo/Cell.h
@@ -41,7 +41,7 @@ class Cell
// whether this cell belongs to a subdomain
// used in relative homology computation
- bool _inSubdomain;
+ bool _subdomain;
// whether this cell a combinded cell of elemetary cells
bool _combined;
@@ -53,57 +53,62 @@ class Cell
bool _immune;
// mutable list of cells on the boundary and on the coboundary of this cell
- std::map< Cell*, int, Less_Cell > _boundary;
- std::map< Cell*, int, Less_Cell > _coboundary;
-
+ std::map< Cell*, int, Less_Cell > _bd;
+ std::map< Cell*, int, Less_Cell > _cbd;
// immutable original boundary and coboundary before the reduction of the
// cell complex
std::map<Cell*, int, Less_Cell > _obd;
std::map<Cell*, int, Less_Cell > _ocbd;
- // The mesh element that is the image of this cell
+ // the mesh element that is the image of this cell
MElement* _image;
-
- // Whether to delete the mesh element when done
+ // whether to delete the mesh element when done
// (created for homology computation only)
- bool _deleteImage;
+ bool _delimage;
// sorted vertices of this cell (used for ordering of the cells)
std::vector<int> _vs;
virtual MVertex* getVertex(int vertex) const {
+ if(_image == NULL) printf("ERROR: No image mesh element for cell. \n");
return _image->getVertex(vertex); }
public:
Cell() : _combined(false), _index(0), _immune(false), _image(NULL),
- _deleteImage(false), _inSubdomain(false) {}
+ _delimage(false), _subdomain(false) {}
Cell(MElement* image);
-
virtual ~Cell();
- virtual MElement* getImageMElement() const { return _image; };
-
-
+ // the mesh element this cell is represented by
+ virtual MElement* getImageMElement() const {
+ if(_image == NULL) printf("ERROR: No image mesh element for cell. \n");
+ return _image; }
// get the number of vertices this cell has
- virtual int getNumVertices() const { return _image->getNumVertices(); }
+ virtual int getNumVertices() const {
+ if(_image == NULL) return _vs.size();
+ else return _image->getNumVertices(); }
// get the number of facets of this cell
virtual int getNumFacets() const;
// get the vertices on a facet of this cell
virtual void getFacetVertices(const int num, std::vector<MVertex*> &v) const;
// get boundary cell orientation
virtual int getFacetOri(Cell* cell);
+ // get/set whether the image mesh element should be deleted when
+ // this cell gets deleted
+ // (was the mesh element in the original mesh,
+ // is it needed after homology computation for visualization?)
+ virtual void setDeleteImage(bool delimage) { _delimage = delimage; }
+ virtual bool getDeleteImage() const { return _delimage; }
+
virtual int getDim() const { return _dim; };
virtual int getIndex() const { return _index; };
virtual void setIndex(int index) { _index = index; };
virtual void setImmune(bool immune) { _immune = immune; };
virtual bool getImmune() const { return _immune; };
- virtual bool inSubdomain() const { return _inSubdomain; }
- virtual void setInSubdomain(bool subdomain) { _inSubdomain = subdomain; }
- virtual void setDeleteImage(bool deleteImage) {
- _deleteImage = deleteImage; };
- virtual bool getDeleteImage() const { return _deleteImage; };
+ virtual bool inSubdomain() const { return _subdomain; }
+ virtual void setInSubdomain(bool subdomain) { _subdomain = subdomain; }
virtual int getSortedVertex(int vertex) const { return _vs.at(vertex); }
// restores the cell information to its original state before reduction
@@ -114,45 +119,47 @@ class Cell
// (co)boundary cell iterator
typedef std::map<Cell*, int, Less_Cell>::iterator biter;
- biter firstBoundary(bool org=false){
- return org ? _obd.begin() : _boundary.begin(); }
- biter lastBoundary(bool org=false){
- return org ? _obd.end() : _boundary.end(); }
- biter firstCoboundary(bool org=false){
- return org ? _ocbd.begin() : _coboundary.begin(); }
- biter lastCoboundary(bool org=false){
- return org ? _ocbd.end() : _coboundary.end(); }
-
- virtual int getBoundarySize() { return _boundary.size(); }
- virtual int getCoboundarySize() { return _coboundary.size(); }
+ biter firstBoundary(bool orig=false){
+ return orig ? _obd.begin() : _bd.begin(); }
+ biter lastBoundary(bool orig=false){
+ return orig ? _obd.end() : _bd.end(); }
+ biter firstCoboundary(bool orig=false){
+ return orig ? _ocbd.begin() : _cbd.begin(); }
+ biter lastCoboundary(bool orig=false){
+ return orig ? _ocbd.end() : _cbd.end(); }
+
+ virtual int getBoundarySize() { return _bd.size(); }
+ virtual int getCoboundarySize() { return _cbd.size(); }
// get the cell boundary
virtual void getBoundary(std::map<Cell*, int, Less_Cell >& boundary,
- bool org=false){
- org ? boundary = _obd : boundary = _boundary; }
+ bool orig=false){
+ orig ? boundary = _obd : boundary = _bd; }
virtual void getCoboundary(std::map<Cell*, int, Less_Cell >& coboundary,
- bool org=false){
- org ? coboundary = _ocbd : coboundary = _coboundary; }
+ bool orig=false){
+ orig ? coboundary = _ocbd : coboundary = _cbd; }
// add (co)boundary cell
- virtual bool addBoundaryCell(int orientation, Cell* cell, bool org=false);
- virtual bool addCoboundaryCell(int orientation, Cell* cell, bool org=false);
+ virtual bool addBoundaryCell(int orientation, Cell* cell,
+ bool orig=false, bool other=false);
+ virtual bool addCoboundaryCell(int orientation, Cell* cell,
+ bool orig=false, bool other=false);
// remove (co)boundary cell
virtual int removeBoundaryCell(Cell* cell, bool other=true);
virtual int removeCoboundaryCell(Cell* cell, bool other=true);
// true if has given cell on (original) (co)boundary
- virtual bool hasBoundary(Cell* cell, bool org=false);
- virtual bool hasCoboundary(Cell* cell, bool org=false);
+ virtual bool hasBoundary(Cell* cell, bool orig=false);
+ virtual bool hasCoboundary(Cell* cell, bool orig=false);
- virtual void clearBoundary() { _boundary.clear(); }
- virtual void clearCoboundary() { _coboundary.clear(); }
+ virtual void clearBoundary() { _bd.clear(); }
+ virtual void clearCoboundary() { _cbd.clear(); }
// print cell debug info
virtual void printCell();
- virtual void printBoundary(bool org=false);
- virtual void printCoboundary(bool org=false);
+ virtual void printBoundary(bool orig=false);
+ virtual void printCoboundary(bool orig=false);
// tools for combined cells
virtual bool isCombined() { return _combined; }
@@ -194,21 +201,11 @@ class CombinedCell : public Cell{
// list of cells this cell is a combination of
std::list< std::pair<int, Cell*> > _cells;
- MVertex* getVertex(int vertex) const {
- printf("ERROR: No mesh vertices for combined cell."); }
-
public:
CombinedCell(Cell* c1, Cell* c2, bool orMatch, bool co=false);
~CombinedCell() {}
- MElement* getImageMElement() const {
- printf("ERROR: No image mesh element for combined cell."); }
- int getNumFacets() const { return 0; }
- void getFacetVertices(const int num, std::vector<MVertex*> &v) const {}
- int getFacetOri(std::vector<MVertex*> &v) { return 0; }
- int getFacetOri(Cell* cell) { return 0; }
-
int getNumVertices() const { return _vs.size(); }
int getSortedVertex(int vertex) const { return _vs.at(vertex); }
diff --git a/Geo/CellComplex.cpp b/Geo/CellComplex.cpp
index 7cd49b9f4f56a058ff62b91af8b4a764857861c0..a938889832d84a36396b3162d1d5dc47ae89d7f2 100644
--- a/Geo/CellComplex.cpp
+++ b/Geo/CellComplex.cpp
@@ -285,11 +285,15 @@ int CellComplex::coreduction(int dim, int omitted)
int CellComplex::reduceComplex(bool omit)
{
- //printf("Cell complex before reduction: %d volumes, %d faces, %d edges and %d vertices.\n", getSize(3), getSize(2), getSize(1), getSize(0));
-
+ printf("Cell Complex: \n %d volumes, %d faces, %d edges and %d vertices. \n",
+ getSize(3), getSize(2), getSize(1), getSize(0));
+
int count = 0;
for(int i = 3; i > 0; i--) count = count + reduction(i);
+ printf(" %d volumes, %d faces, %d edges and %d vertices. \n",
+ getSize(3), getSize(2), getSize(1), getSize(0));
+
if(omit){
int omitted = 0;
_store.clear();
@@ -309,8 +313,16 @@ int CellComplex::reduceComplex(bool omit)
for(int j = 3; j > 0; j--) reduction(j, omitted);
}
}
-
- //printf("Cell complex after reduction: %d volumes, %d faces, %d edges and %d vertices.\n", getSize(3), getSize(2), getSize(1), getSize(0));
+
+ combine(3);
+ reduction(2);
+ combine(2);
+ reduction(1);
+ combine(1);
+
+ printf(" %d volumes, %d faces, %d edges and %d vertices. \n",
+ getSize(3), getSize(2), getSize(1), getSize(0));
+
return 0;
}
@@ -326,12 +338,12 @@ void CellComplex::removeSubdomain()
for(unsigned int i = 0; i < toRemove.size(); i++) removeCell(toRemove[i]);
}
-int CellComplex::coreduceComplex()
+int CellComplex::coreduceComplex(bool omit)
{
- //printf("Cell complex before coreduction: %d volumes, %d faces, %d edges and %d vertices.\n", getSize(3), getSize(2), getSize(1), getSize(0));
+ printf("Cell Complex: \n %d volumes, %d faces, %d edges and %d vertices. \n",
+ getSize(3), getSize(2), getSize(1), getSize(0));
int count = 0;
-
removeSubdomain();
for(int dim = 0; dim < 4; dim++){
@@ -344,23 +356,38 @@ int CellComplex::coreduceComplex()
}
}
- int omitted = 0;
- _store.clear();
-
- while (getSize(0) != 0){
- citer cit = firstCell(0);
- Cell* cell = *cit;
+ printf(" %d volumes, %d faces, %d edges and %d vertices. \n",
+ getSize(3), getSize(2), getSize(1), getSize(0));
- removeCell(cell, false);
- std::set< Cell*, Less_Cell > omittedCells;
- omitted++;
- _store.push_back(omittedCells);
- _store.at(omitted-1).insert(cell);
- coreduction(cell, omitted);
- }
- //printf("Cell complex after coreduction: %d volumes, %d faces, %d edges and %d vertices.\n", getSize(3), getSize(2), getSize(1), getSize(0));
+ int omitted = 0;
+ if(omit){
+ _store.clear();
+
+ while (getSize(0) != 0){
+ citer cit = firstCell(0);
+ Cell* cell = *cit;
+
+ removeCell(cell, false);
+ std::set< Cell*, Less_Cell > omittedCells;
+ omitted++;
+ _store.push_back(omittedCells);
+ _store.at(omitted-1).insert(cell);
+ coreduction(cell, omitted);
+ }
+ }
+ cocombine(0);
+ coreduction(1);
+ cocombine(1);
+ coreduction(2);
+ cocombine(2);
+ coreduction(3);
+ coherent();
+
+ printf(" %d volumes, %d faces, %d edges and %d vertices. \n",
+ getSize(3), getSize(2), getSize(1), getSize(0));
+
return omitted;
}
@@ -485,18 +512,7 @@ int CellComplex::combine(int dim)
return count;
}
-void CellComplex::printComplex(int dim)
-{
- for (citer cit = firstCell(dim); cit != lastCell(dim); cit++){
- Cell* cell = *cit;
- cell->printCell();
- cell->printBoundary();
- cell->printCoboundary();
- //printf("--- \n");
- }
-}
-
-bool CellComplex::checkCoherence()
+bool CellComplex::coherent()
{
bool coherent = true;
for(int i = 0; i < 4; i++){
@@ -510,19 +526,13 @@ bool CellComplex::checkCoherence()
int ori = (*it).second;
citer cit = _cells[bdCell->getDim()].find(bdCell);
if(cit == lastCell(bdCell->getDim())){
- //printf("Warning! Boundary cell not in cell complex! Boundary removed. \n");
- //cell->printCell();
- //bdCell->printCell();
+ printf("Warning! Boundary cell not in cell complex! Boundary removed. \n");
cell->removeBoundaryCell(bdCell);
coherent = false;
}
if(!bdCell->hasCoboundary(cell)){
- //printf("Warning! Incoherent boundary/coboundary pair! Fixed. \n");
- /*cell->printCell();
- cell->printBoundary();
- bdCell->printCell();
- bdCell->printCoboundary();*/
- bdCell->addCoboundaryCell(ori, cell);
+ printf("Warning! Incoherent boundary/coboundary pair! Fixed. \n");
+ bdCell->addCoboundaryCell(ori, cell);
coherent = false;
}
@@ -535,19 +545,13 @@ bool CellComplex::checkCoherence()
int ori = (*it).second;
citer cit = _cells[cbdCell->getDim()].find(cbdCell);
if(cit == lastCell(cbdCell->getDim())){
- //printf("Warning! Coboundary cell not in cell complex! Coboundary removed. \n");
- cell->printCell();
- cbdCell->printCell();
- cell->removeCoboundaryCell(cbdCell);
+ printf("Warning! Coboundary cell not in cell complex! Coboundary removed. \n");
+ cell->removeCoboundaryCell(cbdCell);
coherent = false;
}
if(!cbdCell->hasBoundary(cell)){
- //printf("Warning! Incoherent coboundary/boundary pair! Fixed. \n");
- cell->printCell();
- cell->printCoboundary();
- cbdCell->printCell();
- cbdCell->printBoundary();
- cbdCell->addBoundaryCell(ori, cell);
+ printf("Warning! Incoherent coboundary/boundary pair! Fixed. \n");
+ cbdCell->addBoundaryCell(ori, cell);
coherent = false;
}
@@ -558,18 +562,13 @@ bool CellComplex::checkCoherence()
return coherent;
}
-bool CellComplex::hasCell(Cell* cell, bool org)
+bool CellComplex::hasCell(Cell* cell, bool orig)
{
- if(!org){
- citer cit = _cells[cell->getDim()].find(cell);
- if( cit == lastCell(cell->getDim()) ) return false;
- else return true;
- }
- else{
- citer cit = _ocells[cell->getDim()].find(cell);
- if( cit == lastOrgCell(cell->getDim()) ) return false;
- else return true;
- }
+ citer cit;
+ if(!orig) cit = _cells[cell->getDim()].find(cell);
+ else cit = _ocells[cell->getDim()].find(cell);
+ if( cit == lastCell(cell->getDim(), orig) ) return false;
+ else return true;
}
void CellComplex::getCells(std::set<Cell*, Less_Cell>& cells,
@@ -601,5 +600,15 @@ void CellComplex::restoreComplex()
_store.clear();
}
+void CellComplex::printComplex(int dim)
+{
+ for (citer cit = firstCell(dim); cit != lastCell(dim); cit++){
+ Cell* cell = *cit;
+ cell->printCell();
+ cell->printBoundary();
+ cell->printCoboundary();
+ //printf("--- \n");
+ }
+}
#endif
diff --git a/Geo/CellComplex.h b/Geo/CellComplex.h
index d533fb5fffab6afc155a084e6bb6ca13dba0cf94..246eff7bb4d20c77aa03b3b2964aacb82d456fd4 100644
--- a/Geo/CellComplex.h
+++ b/Geo/CellComplex.h
@@ -70,49 +70,45 @@ class CellComplex
~CellComplex();
// get the number of certain dimensional cells
- int getSize(int dim, bool org=false){
- if(!org) return _cells[dim].size();
+ int getSize(int dim, bool orig=false){
+ if(!orig) return _cells[dim].size();
else return _ocells[dim].size(); }
int getDim() {return _dim; }
-
bool simplicial() { return _simplicial; }
+ // get dim-dimensional cells
+ // domain = 0: cells in domain relative to subdomain
+ // domain = 1: cells in domain
+ // domain = 2: cells in subdomain
void getCells(std::set<Cell*, Less_Cell>& cells, int dim, int domain=0);
- std::set<Cell*, Less_Cell> getOrgCells(int dim){ return _ocells[dim]; }
+ std::set<Cell*, Less_Cell> getOrigCells(int dim){ return _ocells[dim]; }
// iterator for the cells of same dimension
typedef std::set<Cell*, Less_Cell>::iterator citer;
// iterators to the first and last cells of certain dimension
- citer firstCell(int dim) {return _cells[dim].begin(); }
- citer lastCell(int dim) {return _cells[dim].end(); }
- citer firstOrgCell(int dim) {return _ocells[dim].begin(); }
- citer lastOrgCell(int dim) {return _ocells[dim].end(); }
-
+ citer firstCell(int dim, bool orig=false) {
+ return orig ? _ocells[dim].begin() : _cells[dim].begin(); }
+ citer lastCell(int dim, bool orig=false) {
+ return orig ? _ocells[dim].end() : _cells[dim].end(); }
// true if cell complex has given cell
- bool hasCell(Cell* cell, bool org=false);
+ bool hasCell(Cell* cell, bool orig=false);
// check whether two cells both belong to subdomain or if neither one does
bool inSameDomain(Cell* c1, Cell* c2) const {
return ( (!c1->inSubdomain() && !c2->inSubdomain())
|| (c1->inSubdomain() && c2->inSubdomain() )); }
- // (co)reduction of this cell complex
- // removes (co)reduction pairs of cell of dimension dim and dim-1
- int reduction(int dim, int omitted=0);
- int coreduction(int dim, int omitted=0);
-
- // full (co)reduction of this cell complex (all dimensions)
- int reduceComplex(bool omit=true);
- int coreduceComplex();
// remove cells in subdomain from this cell complex
void removeSubdomain();
- // print the vertices of cells of certain dimension
- void printComplex(int dim);
-
+ // (co)reduction of this cell complex
+ // removes (co)reduction pairs of cell of dimension dim and dim-1
+ int reduction(int dim, int omitted=0);
+ int coreduction(int dim, int omitted=0);
+
// Cell combining for reduction and coreduction
int combine(int dim);
int cocombine(int dim);
@@ -121,18 +117,27 @@ class CellComplex
// as coboundary cell and vice versa
// also check whether all (co)boundary cells of a cell
// belong to this cell complex
- bool checkCoherence();
+ bool coherent();
+
+ // full (co)reduction of this cell complex (all dimensions, with combining)
+ // (with highest dimensional cell omitting?)
+ int reduceComplex(bool omit=true);
+ int coreduceComplex(bool imit=true);
int eulerCharacteristic(){
return getSize(0) - getSize(1) + getSize(2) - getSize(3);}
void printEuler(){
printf("Euler characteristic: %d. \n", eulerCharacteristic()); }
+ // get cells omitted by (co)reduction
int getNumOmitted() { return _store.size(); }
std::set<Cell*, Less_Cell> getOmitted(int i) { return _store.at(i); }
+ // restore the cell complex to its original state before (co)reduction
void restoreComplex();
-
+
+ // print the vertices of cells of certain dimension
+ void printComplex(int dim);
};
#endif
diff --git a/Geo/Homology.cpp b/Geo/Homology.cpp
index c256e5734f861cdd04abf990744ee23af0dd2319..eb88889d13c31017977685dd8251eb1044e4f204 100644
--- a/Geo/Homology.cpp
+++ b/Geo/Homology.cpp
@@ -113,37 +113,16 @@ void Homology::findGenerators(CellComplex* cellComplex)
bool ownComplex = false;
if(cellComplex==NULL){
cellComplex = createCellComplex(_domainEntities,
- _subdomainEntities);
+ _subdomainEntities);
ownComplex = true;
}
std::string domainString = getDomainString(_domain, _subdomain);
Msg::Info("Reducing the Cell Complex...");
Msg::StatusBar(1, false, "Reducing...");
- double t1 = Cpu();
-
- printf("Cell Complex: \n %d volumes, %d faces, %d edges and %d vertices. \n",
- cellComplex->getSize(3), cellComplex->getSize(2),
- cellComplex->getSize(1), cellComplex->getSize(0));
-
- int omitted = cellComplex->reduceComplex();
-
- printf(" %d volumes, %d faces, %d edges and %d vertices. \n",
- cellComplex->getSize(3), cellComplex->getSize(2),
- cellComplex->getSize(1), cellComplex->getSize(0));
- cellComplex->combine(3);
- cellComplex->reduction(2);
- cellComplex->combine(2);
- cellComplex->reduction(1);
- cellComplex->combine(1);
-
- printf(" %d volumes, %d faces, %d edges and %d vertices. \n",
- cellComplex->getSize(3), cellComplex->getSize(2),
- cellComplex->getSize(1), cellComplex->getSize(0));
-
- cellComplex->checkCoherence();
-
+ double t1 = Cpu();
+ int omitted = cellComplex->reduceComplex();
double t2 = Cpu();
Msg::Info("Cell Complex reduction complete (%g s).", t2 - t1);
Msg::Info("%d volumes, %d faces, %d edges and %d vertices.",
@@ -191,7 +170,7 @@ void Homology::findGenerators(CellComplex* cellComplex)
}
}
_generators.push_back(chain->createPGroup());
- delete chain;
+ delete chain;
}
if(j == cellComplex->getDim() && cellComplex->getNumOmitted() > 0){
for(int i = 0; i < cellComplex->getNumOmitted(); i++){
@@ -230,37 +209,17 @@ void Homology::findDualGenerators(CellComplex* cellComplex)
bool ownComplex = false;
if(cellComplex==NULL){
cellComplex = createCellComplex(_domainEntities,
- _subdomainEntities);
+ _subdomainEntities);
ownComplex = true;
}
Msg::Info("Reducing Cell Complex...");
Msg::StatusBar(1, false, "Reducing...");
- double t1 = Cpu();
- printf("Cell Complex: \n %d volumes, %d faces, %d edges and %d vertices. \n",
- cellComplex->getSize(3), cellComplex->getSize(2),
- cellComplex->getSize(1), cellComplex->getSize(0));
+ double t1 = Cpu();
int omitted = cellComplex->coreduceComplex();
-
- printf(" %d volumes, %d faces, %d edges and %d vertices. \n",
- cellComplex->getSize(3), cellComplex->getSize(2),
- cellComplex->getSize(1), cellComplex->getSize(0));
-
- cellComplex->cocombine(0);
- cellComplex->coreduction(1);
- cellComplex->cocombine(1);
- cellComplex->coreduction(2);
- cellComplex->cocombine(2);
- cellComplex->coreduction(3);
- cellComplex->checkCoherence();
-
- printf(" %d volumes, %d faces, %d edges and %d vertices. \n",
- cellComplex->getSize(3), cellComplex->getSize(2),
- cellComplex->getSize(1), cellComplex->getSize(0));
-
-
double t2 = Cpu();
+
Msg::Info("Cell Complex reduction complete (%g s).", t2 - t1);
Msg::Info("%d volumes, %d faces, %d edges and %d vertices.",
cellComplex->getSize(3), cellComplex->getSize(2),
@@ -349,31 +308,11 @@ void Homology::findHomSequence(){
_subdomainEntities);
Msg::Info("Reducing the Cell Complex...");
Msg::StatusBar(1, false, "Reducing...");
- double t1 = Cpu();
-
- printf("Cell Complex: \n %d volumes, %d faces, %d edges and %d vertices. \n",
- cellComplex->getSize(3), cellComplex->getSize(2),
- cellComplex->getSize(1), cellComplex->getSize(0));
+ double t1 = Cpu();
cellComplex->reduceComplex(false);
-
- printf(" %d volumes, %d faces, %d edges and %d vertices. \n",
- cellComplex->getSize(3), cellComplex->getSize(2),
- cellComplex->getSize(1), cellComplex->getSize(0));
-
- cellComplex->combine(3);
- cellComplex->reduction(2);
- cellComplex->combine(2);
- cellComplex->reduction(1);
- cellComplex->combine(1);
-
- printf(" %d volumes, %d faces, %d edges and %d vertices. \n",
- cellComplex->getSize(3), cellComplex->getSize(2),
- cellComplex->getSize(1), cellComplex->getSize(0));
-
- cellComplex->checkCoherence();
-
double t2 = Cpu();
+
Msg::Info("Cell Complex reduction complete (%g s).", t2 - t1);
Msg::Info("%d volumes, %d faces, %d edges and %d vertices.",
cellComplex->getSize(3), cellComplex->getSize(2),