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

Cleaned some interdepencies of classes.

parent 3ead45a7
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Geo/Cell.cpp
+ 1
− 9
View file @ 2aeaa427
......@@ -103,7 +103,7 @@ int Cell::getFacetOri(std::vector<MVertex*> &v)
void Cell::printCell()
{
printf("%d-cell %d: \n" , getDim(), getNum());
printf("%d-cell: \n" , getDim());
printf("Vertices: ");
for(int i = 0; i < this->getNumVertices(); i++){
printf("%d ", this->getSortedVertex(i));
......@@ -205,14 +205,6 @@ bool Cell::hasCoboundary(Cell* cell, bool org)
}
}
void Cell::makeDualCell()
{
std::map<Cell*, int, Less_Cell > temp = _boundary;
_boundary = _coboundary;
_coboundary = temp;
_dim = 3-_dim;
}
void Cell::printBoundary(bool org)
{
for(biter it = firstBoundary(org); it != lastBoundary(org); it++){
......
Geo/Cell.h
+ 3
− 10
View file @ 2aeaa427
......@@ -87,10 +87,6 @@ class Cell
virtual int getDim() const { return _dim; };
virtual int getIndex() const { return _index; };
virtual void setIndex(int index) { _index = index; };
virtual int getNum() const { return _image->getNum(); }
virtual int getType() const { return _image->getType(); }
virtual int getTypeForMSH() const { return _image->getTypeForMSH(); }
virtual int getPartition() const { return _image->getPartition(); }
virtual void setImmune(bool immune) { _immune = immune; };
virtual bool getImmune() const { return _immune; };
virtual bool inSubdomain() const { return _inSubdomain; }
......@@ -146,10 +142,7 @@ class Cell
virtual void clearBoundary() { _boundary.clear(); }
virtual void clearCoboundary() { _coboundary.clear(); }
// algebraic dual of the cell
virtual void makeDualCell();
// print cell info
// print cell debug info
virtual void printCell();
virtual void printBoundary(bool org=false);
virtual void printCoboundary(bool org=false);
......@@ -160,14 +153,14 @@ class Cell
virtual void getFacetVertices(const int num, std::vector<MVertex*> &v) const;
// get boundary cell orientation
virtual int getFacetOri(std::vector<MVertex*> &v);
virtual int getFacetOri(Cell* cell) {
std::vector<MVertex*> v;
for(int i = 0; i < cell->getNumVertices(); i++) {
v.push_back(cell->getVertex(i));
}
return getFacetOri(v);
}
virtual int getFacetOri(std::vector<MVertex*> &v);
}
// tools for combined cells
virtual bool isCombined() { return _combined; }
......
Geo/CellComplex.cpp
+ 61
− 91
View file @ 2aeaa427
......@@ -9,37 +9,21 @@
#if defined(HAVE_KBIPACK)
CellComplex::CellComplex( std::vector<GEntity*> domain,
std::vector<GEntity*> subdomain)
CellComplex::CellComplex( std::vector<MElement*>& domainElements,
std::vector<MElement*>& subdomainElements)
{
_domain = domain;
_subdomain = subdomain;
_dim = 0;
_simplicial = true;
_multidim = false;
int dim = 0;
for(unsigned int i = 0; i < domain.size(); i++){
GEntity* entity = domain.at(i);
if(i == 0) dim = entity->dim();
if(dim != entity->dim()){
_multidim = true;
//printf("Warning: Domain is not a manifold.");
break;
}
}
// insert cells into cell complex
// subdomain need to be inserted first!
if(!insert_cells(_subdomain, true)){ panic_exit(); return; }
if(!insert_cells(_domain, false)) { panic_exit(); return; }
// subdomain needs to be inserted first!
if(!insert_cells(subdomainElements, true)){ panic_exit(); return; }
if(!insert_cells(domainElements, false)) { panic_exit(); return; }
for(int i = 0; i < 4; i++){
_ocells[i] = _cells[i];
if(getSize(i) > _dim) _dim = i;
}
}
}
void CellComplex::panic_exit(){
......@@ -53,68 +37,50 @@ void CellComplex::panic_exit(){
}
for(unsigned int i = 0; i < _newcells.size(); i++) delete _newcells.at(i);
_newcells.clear();
_domain.clear();
_subdomain.clear();
}
bool CellComplex::insert_cells(std::vector<GEntity*>& domain, bool subdomain)
bool CellComplex::insert_cells(std::vector<MElement*>& elements,
bool subdomain)
{
std::vector<MVertex*> vertices;
std::pair<citer, bool> insertInfo;
// add highest dimensional cells
for(unsigned int j=0; j < domain.size(); j++) {
for(unsigned int i=0; i < domain.at(j)->getNumMeshElements(); i++){
vertices.clear();
MElement* element = domain.at(j)->getMeshElement(i);
for(int k=0; k < element->getNumVertices(); k++){
MVertex* vertex = element->getVertex(k);
vertices.push_back(vertex);
}
int dim = element->getDim();
int type = element->getTypeForMSH();
for(unsigned int i=0; i < elements.size(); i++){
vertices.clear();
MElement* element = elements.at(i);
for(int k=0; k < element->getNumVertices(); k++){
MVertex* vertex = element->getVertex(k);
vertices.push_back(vertex);
}
Cell* cell;
// simplex types
if(type == MSH_LIN_2 || type == MSH_TRI_3 || type == MSH_TET_4
|| type == MSH_LIN_3 || type == MSH_TRI_6 || type == MSH_TET_10
|| type == MSH_PNT || type == MSH_TRI_9 || type == MSH_TRI_10
|| type == MSH_TRI_12 || type == MSH_TRI_15 || type == MSH_TRI_15I
|| type == MSH_TRI_21 || type == MSH_LIN_4 || type == MSH_LIN_5
|| type == MSH_LIN_6 || type == MSH_TET_20 || type == MSH_TET_35
|| type == MSH_TET_56
|| type == MSH_TET_34 || type == MSH_TET_52 ){
cell = new Cell(element);
}
else if(type == MSH_QUA_4 || type == MSH_QUA_8 || type == MSH_QUA_9){
cell = new Cell(element);
_simplicial = false;
}
else if(type == MSH_HEX_8 || type == MSH_HEX_27 || type == MSH_HEX_20){
cell = new Cell(element);
_simplicial = false;
}
else if(type == MSH_PRI_6 || type == MSH_PRI_18 || type == MSH_PRI_15){
cell = new Cell(element);
_simplicial = false;
}/* FIXME: no getFaceInfo methods for these MElements
else if(type == MSH_PYR_5 || type == MSH_PYR_14 || type == MSH_PYR_13){
cell = new Cell(element, subdomain, boundary);
_simplicial = false;
}*/
else {
//printf("Error: mesh element %d not implemented yet! \n", type);
return false;
}
cell->setImmune(false);
cell->setInSubdomain(subdomain);
insertInfo = _cells[dim].insert(cell);
if(!insertInfo.second) delete cell;
int dim = element->getDim();
int type = element->getTypeForMSH();
Cell* cell = new Cell(element);
// simplex types
if( !(type == MSH_PNT
|| type == MSH_LIN_2 || type == MSH_TRI_3 || type == MSH_TET_4
|| type == MSH_LIN_3 || type == MSH_TRI_6 || type == MSH_TET_10
|| type == MSH_PNT || type == MSH_TRI_9 || type == MSH_TRI_10
|| type == MSH_TRI_12 || type == MSH_TRI_15 || type == MSH_TRI_15I
|| type == MSH_TRI_21 || type == MSH_LIN_4 || type == MSH_LIN_5
|| type == MSH_LIN_6 || type == MSH_TET_20 || type == MSH_TET_35
|| type == MSH_TET_56 || type == MSH_TET_34 || type == MSH_TET_52) ){
_simplicial = false;
}
/* FIXME: no getFaceInfo methods for these MElements */
if(type == MSH_PYR_5 || type == MSH_PYR_14 || type == MSH_PYR_13){
//printf("Error: mesh element %d not implemented yet! \n", type);
return false;
}
cell->setImmune(false);
cell->setInSubdomain(subdomain);
insertInfo = _cells[dim].insert(cell);
if(!insertInfo.second) delete cell;
}
// add lower dimensional cells recursively
MElementFactory factory;
for (int dim = 3; dim > 0; dim--){
......@@ -123,7 +89,7 @@ bool CellComplex::insert_cells(std::vector<GEntity*>& domain, bool subdomain)
std::vector<MVertex*> vertices;
for(int i = 0; i < cell->getNumFacets(); i++){
cell->getFacetVertices(i, vertices);
int type = cell->getTypeForMSH();
int type = cell->getImageMElement()->getTypeForMSH();
int newtype = 0;
//FIXME: add missing boundary cell type relations
//FIXME: high order meshes don't work
......@@ -131,30 +97,34 @@ bool CellComplex::insert_cells(std::vector<GEntity*>& domain, bool subdomain)
if(type == MSH_TET_4) newtype = MSH_TRI_3;
/*else if(type == MSH_TET_10) newtype = MSH_TRI_6;
else if(type == MSH_TET_20) newtype = MSH_TRI_9;*/
else if(type == MSH_HEX_8) newtype = MSH_QUA_4;
/*else if(type == MSH_HEX_8) newtype = MSH_QUA_4;*/
/*else if(type == MSH_HEX_20) newtype = MSH_QUA_8;
else if(type == MSH_HEX_27) newtype = MSH_QUA_9;*/
else if(type == MSH_PRI_6 && vertices.size() == 3) newtype = MSH_TRI_3;
else if(type == MSH_PRI_6 && vertices.size() == 4) newtype = MSH_QUA_4;
else if(type == MSH_PRI_6
&& vertices.size() == 3) newtype = MSH_TRI_3;
else if(type == MSH_PRI_6
&& vertices.size() == 4) newtype = MSH_QUA_4;
}
else if(dim == 2){
if(type == MSH_TRI_3 || type == MSH_QUA_4) newtype = MSH_LIN_2;
/*else if(type == MSH_TRI_6 || type == MSH_QUA_8) newtype = MSH_LIN_3;
else if(type == MSH_TRI_9) newtype = MSH_LIN_4;
else if(type == MSH_QUA_9) newtype = MSH_LIN_3;*/
}
else if(dim == 1){
if(type == MSH_LIN_2) newtype = MSH_PNT;
/*else if(type == MSH_LIN_3 || type == MSH_LIN_4 ||
type == MSH_LIN_5 || type == MSH_LIN_6) newtype = MSH_PNT;*/
}
if(newtype == 0){
if(type == MSH_TRI_3 || type == MSH_QUA_4) newtype = MSH_LIN_2;
/*else if(type == MSH_TRI_6
|| type == MSH_QUA_8) newtype = MSH_LIN_3;
else if(type == MSH_TRI_9) newtype = MSH_LIN_4;
else if(type == MSH_QUA_9) newtype = MSH_LIN_3;*/
}
else if(dim == 1){
if(type == MSH_LIN_2) newtype = MSH_PNT;
/*else if(type == MSH_LIN_3 || type == MSH_LIN_4 ||
type == MSH_LIN_5 || type == MSH_LIN_6) newtype = MSH_PNT;*/
}
if(newtype == 0){
//printf("Error: mesh element %d not implemented yet! \n", type);
return false;
}
MElement* element = factory.create(newtype, vertices, 0,
cell->getPartition());
cell->getImageMElement()->
getPartition());
Cell* newCell = new Cell(element);
newCell->setImmune(cell->getImmune());
newCell->setInSubdomain(subdomain);
......
Geo/CellComplex.h
+ 4
− 20
View file @ 2aeaa427
......@@ -20,7 +20,7 @@
#include "Cell.h"
#include "MElement.h"
#include "ChainComplex.h"
#include "GModel.h"
//#include "GModel.h"
class Cell;
......@@ -28,16 +28,6 @@ class Cell;
class CellComplex
{
private:
// the domain in the model which this cell complex covers
std::vector<GEntity*> _domain;
// a subdomain of the given domain
// used in relative homology computation, may be empty
std::vector<GEntity*> _subdomain;
// entities on the boundary of the homology computation domain
std::vector<GEntity*> _boundary;
// sorted containers of unique cells in this cell complex
// one for each dimension
std::set<Cell*, Less_Cell> _cells[4];
......@@ -56,9 +46,6 @@ class CellComplex
// is the cell complex simplicial
bool _simplicial;
// does the domain contain entities of different dimensions
bool _multidim;
// enqueue cells in queue if they are not there already
void enqueueCells(std::map<Cell*, int, Less_Cell>& cells,
std::queue<Cell*>& Q, std::set<Cell*, Less_Cell>& Qset);
......@@ -66,7 +53,7 @@ class CellComplex
void removeCellQset(Cell* cell, std::set<Cell*, Less_Cell>& Qset);
// for constructor
bool insert_cells(std::vector<GEntity*>& domain, bool subdomain);
bool insert_cells(std::vector<MElement*>& elements, bool subdomain);
void panic_exit();
// insert/remove a cell from this cell complex
......@@ -78,13 +65,10 @@ class CellComplex
public:
CellComplex( std::vector<GEntity*> domain, std::vector<GEntity*> subdomain);
//CellComplex(CellComplex* cellComplex)
CellComplex( std::vector<MElement*>& domainElements,
std::vector<MElement*>& subdomainElements);
~CellComplex();
std::vector<GEntity*> getDomain() const { return _domain; }
std::vector<GEntity*> getSubdomain() const { return _subdomain; }
// get the number of certain dimensional cells
int getSize(int dim, bool org=false){
if(!org) return _cells[dim].size();
......
Geo/Homology.cpp
+ 23
− 7
View file @ 2aeaa427
......@@ -69,8 +69,24 @@ CellComplex* Homology::createCellComplex(std::vector<GEntity*>& domainEntities,
if(domainEntities.empty()) Msg::Error("Domain is empty.");
if(subdomainEntities.empty()) Msg::Info("Subdomain is empty.");
CellComplex* cellComplex = new CellComplex(domainEntities,
subdomainEntities);
std::vector<MElement*> domainElements;
std::vector<MElement*> subdomainElements;
for(unsigned int j=0; j < domainEntities.size(); j++) {
for(unsigned int i=0; i < domainEntities.at(j)->getNumMeshElements(); i++){
MElement* element = domainEntities.at(j)->getMeshElement(i);
domainElements.push_back(element);
}
}
for(unsigned int j=0; j < subdomainEntities.size(); j++) {
for(unsigned int i=0; i < subdomainEntities.at(j)->getNumMeshElements();
i++){
MElement* element = subdomainEntities.at(j)->getMeshElement(i);
subdomainElements.push_back(element);
}
}
CellComplex* cellComplex = new CellComplex(domainElements,
subdomainElements);
if(cellComplex->getSize(0) == 0){
Msg::Error("Cell Complex is empty!");
......@@ -94,7 +110,8 @@ Homology::~Homology()
void Homology::findGenerators()
{
CellComplex* cellComplex = createCellComplex(_domainEntities, _subdomainEntities);
CellComplex* cellComplex = createCellComplex(_domainEntities,
_subdomainEntities);
std::string domainString = getDomainString(_domain, _subdomain);
Msg::Info("Reducing the Cell Complex...");
......@@ -320,9 +337,8 @@ void Homology::findDualGenerators()
}
void Homology::findHomSequence(){
CellComplex* cellComplex = new CellComplex(_domainEntities,
_subdomainEntities);
CellComplex* cellComplex = createCellComplex(_domainEntities,
_subdomainEntities);
Msg::Info("Reducing the Cell Complex...");
Msg::StatusBar(1, false, "Reducing...");
double t1 = Cpu();
......@@ -718,7 +734,7 @@ int Chain::writeChainMSH(const std::string &name)
for(citer cit = _cells.begin(); cit != _cells.end(); cit++){
Cell* cell = (*cit).first;
int coeff = (*cit).second;
fprintf(fp, "%d %d \n", cell->getNum(), coeff );
fprintf(fp, "%d %d \n", cell->getImageMElement()->getNum(), coeff );
}
fprintf(fp, "$EndElementData\n");
......
Geo/Homology.h
+ 1
− 0
View file @ 2aeaa427
......@@ -12,6 +12,7 @@
#include "CellComplex.h"
#include "ChainComplex.h"
#include "OS.h"
#include "GModel.h"
#include "PView.h"
#include "Options.h"
......
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