Skip to content
Snippets Groups Projects
Commit bbc7ddfd authored by Matti Pellika's avatar Matti Pellika
Browse files 

Bugfixes.
parent 91fd43c1
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Geo/CellComplex.cpp
+ 58
− 10
View file @ bbc7ddfd
......@@ -272,6 +272,8 @@ void CellComplex::removeCell(Cell* cell){
bdCell->removeCoboundaryCell(cell);
}
//cell->clearBoundary();
//cell->clearCoboundary();
}
......@@ -283,7 +285,8 @@ void CellComplex::enqueueCells(std::list<Cell*>& cells, std::queue<Cell*>& Q, st
for(std::list<Cell*>::iterator cit = cells.begin(); cit != cells.end(); cit++){
Cell* cell = *cit;
citer it = Qset.find(cell);
if(it == Qset.end()){
//citer it2 = _cells[cell->getDim()].find(cell);
if(it == Qset.end()){// && it2 != _cells[cell->getDim()].end()){
Qset.insert(cell);
Q.push(cell);
}
......@@ -519,6 +522,23 @@ void CellComplex::replaceCells(Cell* c1, Cell* c2, Cell* newCell, bool orMatch,
cbdCell->removeBoundaryCell(c1);
cbdCell->addBoundaryCell(ori, newCell, true);
}
/*
for(std::list< std::pair<int, Cell*> >::iterator it = coboundary1.begin(); it != coboundary1.end(); it++){
Cell* cbdCell = (*it).second;
int ori = (*it).first;
cbdCell->removeBoundaryCell(c1);
if(!co){
bool old = false;
for(std::list< std::pair<int, Cell* > >::iterator it2 = coboundary2.begin(); it2 != coboundary2.end(); it2++){
Cell* cell2 = (*it2).second;
if(*cell2 == *cbdCell) old = true;
}
if(!old) cbdCell->addBoundaryCell(ori, newCell, true);
}
else cbdCell->addBoundaryCell(ori, newCell, true);
}
*/
for(std::list< std::pair<int, Cell*> >::iterator it = coboundary2.begin(); it != coboundary2.end(); it++){
Cell* cbdCell = (*it).second;
int ori = (*it).first;
......@@ -535,12 +555,31 @@ void CellComplex::replaceCells(Cell* c1, Cell* c2, Cell* newCell, bool orMatch,
else cbdCell->addBoundaryCell(ori, newCell, true);
}
for(std::list< std::pair<int, Cell* > >::iterator it = boundary1.begin(); it != boundary1.end(); it++){
Cell* bdCell = (*it).second;
int ori = (*it).first;
bdCell->removeCoboundaryCell(c1);
bdCell->addCoboundaryCell(ori, newCell, true);
}
/*
for(std::list< std::pair<int, Cell* > >::iterator it = boundary1.begin(); it != boundary1.end(); it++){
Cell* bdCell = (*it).second;
int ori = (*it).first;
bdCell->removeCoboundaryCell(c1);
if(co){
bool old = false;
for(std::list< std::pair<int, Cell* > >::iterator it2 = boundary2.begin(); it2 != boundary2.end(); it2++){
Cell* cell2 = (*it2).second;
if(*cell2 == *bdCell) old = true;
}
if(!old) bdCell->addCoboundaryCell(ori, newCell, true);
}
else bdCell->addCoboundaryCell(ori, newCell, true);
}
*/
for(std::list< std::pair<int, Cell* > >::iterator it = boundary2.begin(); it != boundary2.end(); it++){
Cell* bdCell = (*it).second;
int ori = (*it).first;
......@@ -560,6 +599,8 @@ void CellComplex::replaceCells(Cell* c1, Cell* c2, Cell* newCell, bool orMatch,
_cells[dim].erase(c1);
_cells[dim].erase(c2);
//removeCell(c1);
//removeCell(c2);
//_trash.insert(c1);
//_trash.insert(c2);
_cells[dim].insert(newCell);
......@@ -601,19 +642,21 @@ int CellComplex::cocombine(int dim){
int or2 = bd_c.back().first;
Cell* c1 = bd_c.front().second;
Cell* c2 = bd_c.back().second;
if( c1->getNumVertices() < getSize(dim) && c2->getNumVertices() < getSize(dim) ){
removeCell(s);
cbd_c = c1->getCoboundary();
enqueueCells(cbd_c, Q, Qset);
cbd_c = c2->getCoboundary();
enqueueCells(cbd_c, Q, Qset);
removeCell(s);
CombinedCell* newCell = new CombinedCell(c1, c2, (or1 != or2), true );
replaceCells(c1, c2, newCell, (or1 != or2), true);
cit = firstCell(dim);
count++;
}
}
removeCellQset(s, Qset);
......@@ -653,23 +696,28 @@ int CellComplex::combine(int dim){
cbd_c = s->getOrientedCoboundary();
if(cbd_c.size() == 2 && !(*(cbd_c.front().second) == *(cbd_c.back().second))
&& inSameDomain(s, cbd_c.front().second) && inSameDomain(s, cbd_c.back().second) ){
&& inSameDomain(s, cbd_c.front().second) && inSameDomain(s, cbd_c.back().second)){
int or1 = cbd_c.front().first;
int or2 = cbd_c.back().first;
Cell* c1 = cbd_c.front().second;
Cell* c2 = cbd_c.back().second;
if( c1->getNumVertices() < getSize(dim) && c2->getNumVertices() < getSize(dim)){
removeCell(s);
bd_c = c1->getBoundary();
enqueueCells(bd_c, Q, Qset);
bd_c = c2->getBoundary();
enqueueCells(bd_c, Q, Qset);
removeCell(s);
CombinedCell* newCell = new CombinedCell(c1, c2, (or1 != or2) );
replaceCells(c1, c2, newCell, (or1 != or2));
//removeCell(s);
cit = firstCell(dim);
//cit++;
count++;
}
}
removeCellQset(s, Qset);
......@@ -807,10 +855,10 @@ int CellComplex::writeComplexMSH(const std::string &name){
void CellComplex::printComplex(int dim){
for (citer cit = firstCell(dim); cit != lastCell(dim); cit++){
Cell* cell = *cit;
for(int i = 0; i < cell->getNumVertices(); i++){
printf("%d ", cell->getVertex(i)->getNum());
}
printf("\n");
cell->printCell();
cell->printBoundary();
cell->printCoboundary();
printf("--- \n");
}
}
......
Geo/CellComplex.h
+ 72
− 6
View file @ bbc7ddfd
......@@ -141,7 +141,15 @@ class Cell
}
return count;
}
virtual bool hasBoundary(Cell* cell){
for(std::list< std::pair<int, Cell*> >::iterator it = _boundary.begin(); it != _boundary.end(); it++){
Cell* cell2 = (*it).second;
if(*cell2 == *cell) return true;
}
return false;
}
virtual void clearBoundary() { _boundary.clear(); }
virtual void clearCoboundary() { _coboundary.clear(); }
......@@ -160,6 +168,14 @@ class Cell
cell2->printCell();
}
}
virtual void printCoboundary() {
for(std::list< std::pair<int, Cell*> >::iterator it = _coboundary.begin(); it != _coboundary.end(); it++){
printf("Coboundary cell orientation: %d ", (*it).first);
Cell* cell2 = (*it).second;
cell2->printCell();
}
}
virtual bool inSubdomain() const { return _inSubdomain; }
......@@ -435,6 +451,7 @@ class CombinedCell : public Cell{
std::vector<MVertex*> _v;
std::vector<int> _vs;
std::list< std::pair<int, Cell*> > _cells;
int _num;
public:
......@@ -442,10 +459,15 @@ class CombinedCell : public Cell{
_index = c1->getIndex();
_tag = c1->getTag();
_dim = c1->getDim();
_num = c1->getNum();
_inSubdomain = c1->inSubdomain();
_onDomainBoundary = c1->onDomainBoundary();
_combined = true;
_v.reserve(c1->getNumVertices() + c2->getNumVertices());
_vs.reserve(c1->getNumVertices() + c2->getNumVertices());
_v = c1->getVertexVector();
for(int i = 0; i < c2->getNumVertices(); i++){
if(!this->hasVertex(c2->getVertex(i)->getNum())) _v.push_back(c2->getVertex(i));
......@@ -461,14 +483,31 @@ class CombinedCell : public Cell{
_cells.push_back(*it);
}
_boundary = c1->getOrientedBoundary();
if(!co) _boundary = c1->getOrientedBoundary();
std::list< std::pair<int, Cell*> > c1Boundary = c1->getOrientedBoundary();
std::list< std::pair<int, Cell*> > c2Boundary = c2->getOrientedBoundary();
/*
for(std::list< std::pair<int, Cell*> >::iterator it = c1Boundary.begin(); it != c1Boundary.end(); it++){
Cell* cell = (*it).second;
if(co){
bool old = false;
for(std::list< std::pair<int, Cell* > >::iterator it2 = c2Boundary.begin(); it2 != c2Boundary.end(); it2++){
Cell* cell2 = (*it2).second;
if(*cell2 == *cell) old = true;
}
if(!old) _boundary.push_back(*it);
}
else _boundary.push_back(*it);
}
*/
for(std::list< std::pair<int, Cell*> >::iterator it = c2Boundary.begin(); it != c2Boundary.end(); it++){
if(!orMatch) (*it).first = -1*(*it).first;
Cell* cell = (*it).second;
if(co){
bool old = false;
for(std::list< std::pair<int, Cell* > >::iterator it2 = _boundary.begin(); it2 != _boundary.end(); it2++){
for(std::list< std::pair<int, Cell* > >::iterator it2 = c1Boundary.begin(); it2 != c1Boundary.end(); it2++){
Cell* cell2 = (*it2).second;
if(*cell2 == *cell) old = true;
}
......@@ -477,14 +516,32 @@ class CombinedCell : public Cell{
else _boundary.push_back(*it);
}
_coboundary = c1->getOrientedCoboundary();
if(co) _coboundary = c1->getOrientedCoboundary();
std::list< std::pair<int, Cell*> > c1Coboundary = c1->getOrientedCoboundary();
std::list< std::pair<int, Cell*> > c2Coboundary = c2->getOrientedCoboundary();
/*
for(std::list< std::pair<int, Cell* > >::iterator it = c1Coboundary.begin(); it != c1Coboundary.end(); it++){
Cell* cell = (*it).second;
if(!co){
bool old = false;
for(std::list< std::pair<int, Cell* > >::iterator it2 = c2Coboundary.begin(); it2 != c2Coboundary.end(); it2++){
Cell* cell2 = (*it2).second;
if(*cell2 == *cell) old = true;
}
if(!old) _coboundary.push_back(*it);
}
else _coboundary.push_back(*it);
}
*/
for(std::list< std::pair<int, Cell* > >::iterator it = c2Coboundary.begin(); it != c2Coboundary.end(); it++){
if(!orMatch) (*it).first = -1*(*it).first;
Cell* cell = (*it).second;
if(!co){
bool old = false;
for(std::list< std::pair<int, Cell* > >::iterator it2 = _coboundary.begin(); it2 != _coboundary.end(); it2++){
for(std::list< std::pair<int, Cell* > >::iterator it2 = c1Coboundary.begin(); it2 != c1Coboundary.end(); it2++){
Cell* cell2 = (*it2).second;
if(*cell2 == *cell) old = true;
}
......@@ -497,6 +554,7 @@ class CombinedCell : public Cell{
~CombinedCell(){}
int getNum() { return _num; }
int getNumVertices() const { return _v.size(); }
MVertex* getVertex(int vertex) const { return _v.at(vertex); }
int getSortedVertex(int vertex) const { return _vs.at(vertex); }
......@@ -506,6 +564,9 @@ class CombinedCell : public Cell{
// true if this cell has given vertex
bool hasVertex(int vertex) const {
//std::vector<int>::iterator it = std::find(_vs.begin(), _vs.end(), vertex);
//if (it == _vs.end()) return false;
//else return true;
for(int i = 0; i < _v.size(); i++){
if(_v.at(i)->getNum() == vertex) return true;
}
......@@ -522,7 +583,6 @@ class CombinedCell : public Cell{
std::list< std::pair<int, Cell*> > getCells() { return _cells; }
int getNumCells() {return _cells.size();}
};
......@@ -626,11 +686,17 @@ class CellComplex
//}
for(int i = 0; i < 4; i++){
for(citer cit = _cells[i].begin(); cit != _cells[i].end(); cit++){
Cell* cell = *cit;
if(cell->combined()) delete cell;
}
_cells[i].clear();
for(citer cit = _cells2[i].begin(); cit != _cells2[i].end(); cit++){
Cell* cell = *cit;
delete cell;
}
_cells2[i].clear();
}
}
......
Geo/ChainComplex.cpp
+ 1
− 1
View file @ bbc7ddfd
......@@ -87,7 +87,7 @@ ChainComplex::ChainComplex(CellComplex* cellComplex){
mpz_set_si(elem, old_elem + (*it).first);
if( (old_elem + (*it).first) > 1 || (old_elem + (*it).first) < -1 ){
printf("Warning: Invalid incidence index: %d! HMatrix: %d. \n", (old_elem + (*it).first), dim);
mpz_set_si(elem, 0);
//mpz_set_si(elem, 0);
}
gmp_matrix_set_elem(elem, bdCell->getIndex(), cell->getIndex(), _HMatrix[dim]);
}
......
Geo/Homology.cpp
+ 4
− 0
View file @ bbc7ddfd
......@@ -72,9 +72,12 @@ Homology::Homology(GModel* model, std::vector<int> physicalDomain, std::vector<i
void Homology::findGenerators(std::string fileName){
Msg::Info("Reducing Cell Complex...");
double t1 = Cpu();
int omitted = _cellComplex->reduceComplex(true);
//_cellComplex->checkCoherence();
if(getCombine()){
_cellComplex->combine(3);
......@@ -90,6 +93,7 @@ void Homology::findGenerators(std::string fileName){
Msg::Info("%d volumes, %d faces, %d edges and %d vertices.",
_cellComplex->getSize(3), _cellComplex->getSize(2), _cellComplex->getSize(1), _cellComplex->getSize(0));
//for(int i = 0; i < 4; i++) _cellComplex->printComplex(i);
_cellComplex->writeComplexMSH(fileName);
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment