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

Added a reduction heuristic

parent 8725626e
No related branches found
No related tags found
No related merge requests found
......@@ -18,6 +18,8 @@ CellComplex::CellComplex(GModel* model,
_model(model), _dim(0), _simplicial(true), _saveorig(saveOriginalComplex),
_relative(false)
{
_smallestCell.second = -1.;
_biggestCell.second = -1.;
_deleteCount = 0;
_createCount = 0;
_insertCells(subdomainElements, 1);
......@@ -62,8 +64,17 @@ CellComplex::CellComplex(GModel* model,
bool CellComplex::_insertCells(std::vector<MElement*>& elements,
int domain)
{
std::pair<Cell*, double> smallestElement[4];
std::pair<Cell*, double> biggestElement[4];
for(int i = 0; i < 4; i++) {
smallestElement[i].second = -1.;
biggestElement[i].second = -1.;
}
_dim = 0;
for(unsigned int i=0; i < elements.size(); i++){
MElement* element = elements.at(i);
int dim = element->getDim();
int type = element->getType();
if(type == TYPE_PYR || type == TYPE_PRI ||
type == TYPE_POLYG || type == TYPE_POLYH) {
......@@ -77,6 +88,8 @@ bool CellComplex::_insertCells(std::vector<MElement*>& elements,
delete maybeCell.first;
continue;
}
if(_dim < dim) _dim = dim;
Cell* cell = maybeCell.first;
std::pair<citer, bool> insert =
_cells[cell->getDim()].insert(cell);
......@@ -86,7 +99,17 @@ bool CellComplex::_insertCells(std::vector<MElement*>& elements,
if(domain) cell->setDomain(domain);
}
else _createCount++;
if(domain == 0) {
double size = fabs(element->getVolume());
if(smallestElement[dim].second < 0. || smallestElement[dim].second > size)
smallestElement[dim] = std::make_pair(cell, size);
if(biggestElement[dim].second < 0. || biggestElement[dim].second < size)
biggestElement[dim] = std::make_pair(cell, size);
}
}
_smallestCell = smallestElement[_dim];
_biggestCell = biggestElement[_dim];
for (int dim = 3; dim > 0; dim--){
for(citer cit = firstCell(dim); cit != lastCell(dim); cit++){
......@@ -109,6 +132,10 @@ bool CellComplex::_insertCells(std::vector<MElement*>& elements,
if(domain == 0) {
int ori = cell->findBdCellOrientation(newCell, i);
cell->addBoundaryCell( ori, newCell, true);
if(_smallestCell.first == cell)
_smallestCell = std::make_pair(newCell, _smallestCell.second);
if(_biggestCell.first == cell)
_biggestCell = std::make_pair(newCell, _biggestCell.second);
}
}
}
......@@ -508,7 +535,7 @@ void CellComplex::removeCells(int dim)
_reduced = true;
}
int CellComplex::coreduceComplex(int combine, bool omit)
int CellComplex::coreduceComplex(int combine, bool omit, int heuristic)
{
Msg::Debug("Cell Complex coreduction:");
Msg::Debug(" %d volumes, %d faces, %d edges, and %d vertices",
......@@ -537,9 +564,23 @@ int CellComplex::coreduceComplex(int combine, bool omit)
std::vector<Cell*> newCells;
while (getSize(0) != 0){
citer cit = firstCell(0);
Cell* cell = *cit;
if(heuristic == -1 && _smallestCell.second != 0. &&
hasCell(_smallestCell.first)) {
Msg::Debug("Omitted a cell in the smallest mesh element with volume %g",
_smallestCell.second);
cell = _smallestCell.first;
}
else if(heuristic == 1 && _biggestCell.second != 0. &&
hasCell(_biggestCell.first)) {
Msg::Debug("Omitted a cell in the biggest mesh element with volume %g",
_biggestCell.second);
cell = _biggestCell.first;
}
newCells.push_back(_omitCell(cell, true));
}
for(unsigned int i = 0; i < newCells.size(); i++){
......@@ -599,7 +640,7 @@ int CellComplex::combine(int dim)
inSameDomain(s, c1) && inSameDomain(s, c2) &&
c1->getImmune() == c2->getImmune() ) {
removeCell(s, true, true);
removeCell(s, true, false);
c1->getBoundary(bd_c);
enqueueCells(bd_c, Q, Qset);
......@@ -615,10 +656,6 @@ int CellComplex::combine(int dim)
cit = firstCell(dim);
count++;
if(s->isCombined()) {
delete s;
_deleteCount++;
}
if(c1->isCombined()) {
delete c1;
_deleteCount++;
......@@ -656,6 +693,7 @@ int CellComplex::cocombine(int dim)
for(citer cit = firstCell(dim); cit != lastCell(dim); cit++){
Cell* cell = *cit;
cell->getCoboundary(cbd_c);
enqueueCells(cbd_c, Q, Qset);
......@@ -674,7 +712,8 @@ int CellComplex::cocombine(int dim)
if(!(*c1 == *c2) && abs(or1) == abs(or2)
&& inSameDomain(s, c1) && inSameDomain(s, c2)){
removeCell(s, true, true);
removeCell(s, true, false);
c1->getCoboundary(cbd_c);
enqueueCells(cbd_c, Q, Qset);
......@@ -691,10 +730,6 @@ int CellComplex::cocombine(int dim)
cit = firstCell(dim);
count++;
if(s->isCombined()) {
delete s;
_deleteCount++;
}
if(c1->isCombined()) {
delete c1;
_deleteCount++;
......
......@@ -25,6 +25,10 @@ class CellComplex
{
private:
std::pair<Cell*, double> _smallestCell;
std::pair<Cell*, double> _biggestCell;
GModel* _model;
// sorted containers of unique cells in this cell complex
......@@ -158,8 +162,11 @@ class CellComplex
// (combine = 1 -> with combining)
// (omit = true -> with highest dimensional cell omitting?)
// (homseq = true -> homology sequence splitting possible after reduction)
// (heuristic = 0 -> no heuristic, let mesh indexing determine
// 1 -> omit 0-cell in biggest element
// -1 -> omit 0-cell in smallest element)
int reduceComplex(int combine=1, bool omit=true, bool homseq=false);
int coreduceComplex(int combine=1, bool omit=true);
int coreduceComplex(int combine=1, bool omit=true, int heuristic=0);
bool isReduced() const { return _reduced; }
......
......@@ -19,10 +19,11 @@ Homology::Homology(GModel* model,
std::vector<int> physicalSubdomain,
std::vector<int> physicalImdomain,
bool saveOrig,
int combine, bool omit, bool smoothen) :
int combine, bool omit, bool smoothen, int heuristic) :
_model(model), _domain(physicalDomain), _subdomain(physicalSubdomain),
_imdomain(physicalImdomain), _saveOrig(saveOrig),
_combine(combine), _omit(omit), _smoothen(smoothen), _cellComplex(NULL)
_combine(combine), _omit(omit), _smoothen(smoothen),
_heuristic(heuristic), _cellComplex(NULL)
{
_fileName = "";
......@@ -49,6 +50,9 @@ Homology::Homology(GModel* model,
_cohomologyComputed[i] = false;
_betti[i] = -1;
}
if(abs(_heuristic) > 1) _heuristic = 0;
}
void Homology::_getEntities(const std::vector<int>& physicalGroups,
......@@ -257,7 +261,7 @@ void Homology::findCohomologyBasis(std::vector<int> dim)
double t1 = Cpu();
int omitted = _cellComplex->coreduceComplex(_combine, _omit);
int omitted = _cellComplex->coreduceComplex(_combine, _omit, _heuristic);
std::sort(dim.begin(), dim.end());
if(!dim.empty() && _combine > 1) {
......
......@@ -49,6 +49,8 @@ class Homology
bool _omit;
// use chain smoothning
bool _smoothen;
// corecution heuristic
int _heuristic;
// file name to store the results
std::string _fileName;
......@@ -93,7 +95,8 @@ class Homology
std::vector<int> physicalSubdomain,
std::vector<int> physicalIm,
bool saveOrig=true,
int combine=2, bool omit=true, bool smoothen=true);
int combine=2, bool omit=true, bool smoothen=true,
int heuristic=1);
~Homology();
GModel* getModel() const { return _model; }
......
......@@ -25,7 +25,8 @@ StringXNumber HomologyComputationOptions_Number[] = {
{GMSH_FULLRC, "CreatePostProcessingViews", NULL, 1.},
{GMSH_FULLRC, "ReductionOmit", NULL, 1.},
{GMSH_FULLRC, "ReductionCombine", NULL, 2.},
{GMSH_FULLRC, "PostProcessSimplify", NULL, 1.}
{GMSH_FULLRC, "PostProcessSimplify", NULL, 1.},
{GMSH_FULLRC, "ReductionHeuristic", NULL, 1.}
};
StringXString HomologyComputationOptions_String[] = {
......@@ -111,6 +112,7 @@ PView *GMSH_HomologyComputationPlugin::execute(PView *v)
bool omit = (bool)HomologyComputationOptions_Number[5].def;
int combine = (int)HomologyComputationOptions_Number[6].def;
bool smoothen = (bool)HomologyComputationOptions_Number[7].def;
int heuristic = (int)HomologyComputationOptions_Number[8].def;
std::vector<int> domain;
std::vector<int> subdomain;
......@@ -124,7 +126,7 @@ PView *GMSH_HomologyComputationPlugin::execute(PView *v)
GModel* m = GModel::current();
Homology* homology = new Homology(m, domain, subdomain, imdomain,
true, combine, omit, smoothen);
true, combine, omit, smoothen, heuristic);
homology->setFileName(fileName);
if(hom != 0) homology->findHomologyBasis(dimsave);
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment