diff --git a/Geo/GModel.cpp b/Geo/GModel.cpp
index 2c1cf48bb64154d7f21ab1901fe5b4aaa6fc9c52..66ebc4422f69136bf99f14d6a1ebd61932e63c02 100644
--- a/Geo/GModel.cpp
+++ b/Geo/GModel.cpp
@@ -969,14 +969,7 @@ void GModel::storeChain(int dim,
{
// create new discrete entities that have no associated MVertices
_storeElementsInEntities(entityMap);
- // give them physical tags
_storePhysicalTagsInEntities(dim, physicalMap);
- // For some weird reason, mesh element color by elementary/physical
- // entity causes crashes when drawing mesh elements of the new discrete
- // entities
- // (somehow related to MVertex.onWhat() as MVertices are unaware of the
- // new discrete entities)
- CTX::instance()->mesh.colorCarousel = 0; // color by element type
}
template<class T>
@@ -1113,6 +1106,23 @@ void GModel::_storeVerticesInEntities(std::vector<MVertex*> &vertices)
}
}
+void GModel::_pruneMeshVertexAssociations()
+{
+ std::vector<GEntity*> entities;
+ std::vector<MVertex*> vertices;
+ getEntities(entities);
+ for(unsigned int i = 0; i < entities.size(); i++) {
+ for(unsigned int j = 0; j < entities[i]->mesh_vertices.size(); j++) {
+ MVertex* v = entities[i]->mesh_vertices[j];
+ v->setEntity(0);
+ vertices.push_back(v);
+ }
+ entities[i]->mesh_vertices.clear();
+ }
+ _associateEntityWithMeshVertices();
+ _storeVerticesInEntities(vertices);
+}
+
void GModel::checkMeshCoherence(double tolerance)
{
int numEle = getNumMeshElements();
@@ -2752,7 +2762,7 @@ void GModel::computeHomology()
std::multimap<dpair, std::string>::iterator> itp =
_homologyRequests.equal_range(*it);
Homology* homology = new Homology(this, itp.first->first.first,
- itp.first->first.second);
+ itp.first->first.second, false);
CellComplex *cellcomplex = homology->createCellComplex();
if(cellcomplex->getSize(0)){
for(std::multimap<dpair, std::string>::iterator itt = itp.first;
@@ -2767,6 +2777,7 @@ void GModel::computeHomology()
else
Msg::Error("Unknown type of homology computation: %s", type.c_str());
}
+ _pruneMeshVertexAssociations();
}
delete cellcomplex;
delete homology;
diff --git a/Geo/GModel.h b/Geo/GModel.h
index 092994f17d54f968fefe6037f33ad5e6890a39ec..422e3c7173afb89d4d50dbcceee7ebdc02f8b850 100644
--- a/Geo/GModel.h
+++ b/Geo/GModel.h
@@ -102,6 +102,12 @@ class GModel
void _storeVerticesInEntities(std::map<int, MVertex*> &vertices);
void _storeVerticesInEntities(std::vector<MVertex*> &vertices);
+ // remove all mesh vertex associations to geometrical entities and
+ // remove vertices from geometrical entities, then
+ // _associateEntityWithMeshVertices and _storeVerticesInEntities
+ // are called to rebuild the associations
+ void _pruneMeshVertexAssociations();
+
// store the physical tags in the geometrical entities
void _storePhysicalTagsInEntities(int dim,
std::map<int, std::map<int, std::string> > &map);
diff --git a/Geo/Homology.cpp b/Geo/Homology.cpp
index 1877a1555d07ec1a35f0d623230c76a44525cd8e..ee3075550f291da1d1b229d9707799e87b35328c 100644
--- a/Geo/Homology.cpp
+++ b/Geo/Homology.cpp
@@ -15,10 +15,10 @@
#if defined(HAVE_KBIPACK)
Homology::Homology(GModel* model, std::vector<int> physicalDomain,
- std::vector<int> physicalSubdomain,
+ std::vector<int> physicalSubdomain, bool save0N,
bool combine, bool omit, bool smoothen) :
_model(model), _domain(physicalDomain), _subdomain(physicalSubdomain),
- _combine(combine), _omit(omit), _smoothen(smoothen)
+ _save0N(save0N), _combine(combine), _omit(omit), _smoothen(smoothen)
{
_fileName = "";
@@ -137,6 +137,8 @@ void Homology::findGenerators(CellComplex* cellComplex)
t2 = Cpu();
Msg::StatusBar(2, true, "Done computing homology spaces (%g s)", t2 - t1);
+ int dim = cellComplex->getDim();
+
int HRank[4];
for(int j = 0; j < 4; j++){
HRank[j] = 0;
@@ -165,7 +167,7 @@ void Homology::findGenerators(CellComplex* cellComplex)
// FIXME: Cell* and CellComplex* pointers should not outlive
// the objects, don't store Chains containing them for now
//_basisChains[chain->createPGroup()] = chain;
- chain->createPGroup();
+ if(!_save0N && (j != 0 && j != dim)) chain->createPGroup();
delete chain;
}
}
@@ -246,7 +248,7 @@ void Homology::findDualGenerators(CellComplex* cellComplex)
// FIXME: Cell* and CellComplex* pointers should not outlive
// the objects, don't store Chains containing them for now
//_basisChains[chain->createPGroup()] = chain;
- chain->createPGroup();
+ if(!_save0N && (j != 0 && j != dim)) chain->createPGroup();
delete chain;
}
}
@@ -266,114 +268,6 @@ void Homology::findDualGenerators(CellComplex* cellComplex)
HRank[0], HRank[1], HRank[2], HRank[3]);
}
-/*void Homology::findHomSequence(){
- CellComplex* cellComplex = createCellComplex(_domainEntities,
- _subdomainEntities);
- Msg::StatusBar(2, true, "Reducing cell complex...");
-
- double t1 = Cpu();
- cellComplex->reduceComplex();
- double t2 = Cpu();
-
- Msg::StatusBar(2, true, "Done reducing cell complex (%g s)", t2 - t1);
- Msg::Info("%d volumes, %d faces, %d edges and %d vertices",
- cellComplex->getSize(3), cellComplex->getSize(2),
- cellComplex->getSize(1), cellComplex->getSize(0));
-
- Msg::StatusBar(2, true, "Computing homology spaces...");
- t1 = Cpu();
-
- ChainComplex* subcomplex = new ChainComplex(cellComplex, 2);
- ChainComplex* complex = new ChainComplex(cellComplex, 1);
- ChainComplex* relcomplex = new ChainComplex(cellComplex, 0);
-
- subcomplex->computeHomology();
- complex->computeHomology();
- relcomplex->computeHomology();
-
- t2 = Cpu();
- Msg::StatusBar(2, true, "Done compuring homology spaces (%g s)", t2 - t1);
-
- Msg::StatusBar(2, true, "Computing homology sequence...");
- HomologySequence* seq = new HomologySequence(subcomplex,
- complex, relcomplex);
- t1 = Cpu();
- Msg::StatusBar(2, true, "Done computing homology sequence (%g s)", t1 - t2);
-
- for(int task = 0; task < 3; task++){
- ChainComplex* chains;
-
- std::string domainString = "";
- std::vector<int> empty;
- if(task == 0) {
- chains = subcomplex;
- domainString = getDomainString(_subdomain, empty);
- }
- else if(task == 1) {
- chains = complex;
- domainString = getDomainString(_domain, empty);
- }
- else{
- chains = relcomplex;
- domainString = getDomainString(_domain, _subdomain);
- }
-
- int HRank[4];
- for(int j = 0; j < 4; j++){
- HRank[j] = 0;
- std::string dimension = "";
- convert(j, dimension);
- for(int i = 1; i <= chains->getBasisSize(j, 3); i++){
-
- std::string generator = "";
- convert(i, generator);
-
- std::string name = "H" + dimension + domainString + generator;
- std::map<Cell*, int, Less_Cell> protoChain;
- chains->getBasisChain(protoChain, i, j, 3, true);
- Chain* chain = new Chain(protoChain, cellComplex, _model, name,
- chains->getTorsion(j,i));
- if(chain->getSize() == 0) {
- delete chain;
- continue;
- }
- HRank[j] = HRank[j] + 1;
- if(chain->getTorsion() != 1){
- Msg::Warning("H%d %d has torsion coefficient %d!",
- j, i, chain->getTorsion());
- }
- _generators.push_back(chain->createPGroup());
- delete chain;
- }
-
- }
- if(task == 0){
- Msg::Info("Ranks of relative homology spaces:");
- }
- if(task == 1){
- Msg::Info("Ranks of absolute homology spaces:");
- }
- if(task == 2){
- Msg::Info("Ranks of absolute homology spaces of relative subdomain:");
- }
- Msg::Info("H0 = %d", HRank[0]);
- Msg::Info("H1 = %d", HRank[1]);
- Msg::Info("H2 = %d", HRank[2]);
- Msg::Info("H3 = %d", HRank[3]);
-
- Msg::StatusBar(2, false, "H0: %d, H1: %d, H2: %d, H3: %d",
- HRank[0], HRank[1], HRank[2], HRank[3]);
- }
-
- delete subcomplex;
- delete complex;
- delete relcomplex;
-
- if(_fileName != "") writeGeneratorsMSH();
- delete cellComplex;
- delete seq;
-}
-*/
std::string Homology::getDomainString(const std::vector<int>& domain,
const std::vector<int>& subdomain)
{
diff --git a/Geo/Homology.h b/Geo/Homology.h
index 9b664a79e017ce89ba9ebd34945e8fab1af97974..de22d7d2b03b189e6e42abdddb6eaf3522727ee1 100644
--- a/Geo/Homology.h
+++ b/Geo/Homology.h
@@ -50,6 +50,9 @@ class Homology
// use chain smoothning
bool _smoothen;
+ // save chains of 0 and highest dimension
+ bool _save0N;
+
int _maxdomain;
int _maxnum;
@@ -61,7 +64,7 @@ class Homology
public:
Homology(GModel* model, std::vector<int> physicalDomain,
- std::vector<int> physicalSubdomain,
+ std::vector<int> physicalSubdomain, bool save0N=false,
bool combine=true, bool omit=true, bool smoothen=true);
~Homology();