diff --git a/FunctionSpace/FunctionSpace.cpp b/FunctionSpace/FunctionSpace.cpp
index fd3e2501bfde178f0fb1d6606f56932399bff76c..874cafc5b94980648b2a3914e0d718926c906745 100644
--- a/FunctionSpace/FunctionSpace.cpp
+++ b/FunctionSpace/FunctionSpace.cpp
@@ -19,14 +19,8 @@ FunctionSpace::~FunctionSpace(void){
delete basis;
// Dof //
- if(dof){
- set<const Dof*>::iterator dStop = dof->end();
- set<const Dof*>::iterator dIt = dof->begin();
-
- for(; dIt != dStop; dIt++)
- delete *dIt;
+ if(dof)
delete dof;
- }
// Group //
if(group){
@@ -96,7 +90,7 @@ void FunctionSpace::buildDof(void){
const vector<const MElement*>& element = goe->getAll();
// Init Struct //
- dof = new set<const Dof*, DofComparator>;
+ dof = new set<Dof>;
group = new vector<GroupOfDof*>(nElement);
eToGod = new map<const MElement*,
const GroupOfDof*,
@@ -109,13 +103,15 @@ void FunctionSpace::buildDof(void){
size_t nDof = myDof.size();
// Add Dof
- vector<const Dof*> trueDof(nDof);
-
for(size_t j = 0; j < nDof; j++)
- insertDof(myDof[j], trueDof, j);
+ dof->insert(myDof[j]);
+ //vector<Dof> trueDof(nDof);
+
+ //for(size_t j = 0; j < nDof; j++)
+ //insertDof(myDof[j], trueDof, j);
// Create new GroupOfDof
- GroupOfDof* god = new GroupOfDof(*(element[i]), trueDof);
+ GroupOfDof* god = new GroupOfDof(*(element[i]), myDof);
(*group)[i] = god;
// Map GOD
@@ -123,7 +119,7 @@ void FunctionSpace::buildDof(void){
(element[i], god));
}
}
-
+/*
void FunctionSpace::insertDof(Dof& d,
vector<const Dof*>& trueDof,
size_t index){
@@ -146,7 +142,7 @@ void FunctionSpace::insertDof(Dof& d,
trueDof[index] = *(p.first);
}
}
-
+*/
vector<Dof> FunctionSpace::getKeys(const MElement& elem) const{
// Const_Cast //
MElement& element = const_cast<MElement&>(elem);
@@ -165,6 +161,7 @@ vector<Dof> FunctionSpace::getKeys(const MElement& elem) const{
// New Element
MElementFactory factory;
+ //MElement* permElement = factory.create(elem.getLowOrderTypeForMSH(), vertex);
MElement* permElement = factory.create(elem.getTypeForMSH(), vertex);
// Edge & Face from Permuted Element //
diff --git a/FunctionSpace/FunctionSpace.h b/FunctionSpace/FunctionSpace.h
index 98ee56d6955f220776befbaef63536755dfdd4fb..bdaed37a29aed24bb4f0c9e3be1b53d88423e832 100644
--- a/FunctionSpace/FunctionSpace.h
+++ b/FunctionSpace/FunctionSpace.h
@@ -59,8 +59,8 @@ class FunctionSpace{
bool scalar;
// Dofs //
- std::set<const Dof*, DofComparator>* dof;
- std::vector<GroupOfDof*>* group;
+ std::set<Dof>* dof;
+ std::vector<GroupOfDof*>* group;
std::map<
const MElement*,
const GroupOfDof*, ElementComparator>* eToGod;
@@ -80,7 +80,7 @@ class FunctionSpace{
std::vector<Dof> getKeys(const MEdge& edge) const;
std::vector<Dof> getKeys(const MFace& face) const;
- const std::vector<const Dof*> getAllDofs(void) const;
+ const std::vector<Dof> getAllDofs(void) const;
const std::vector<GroupOfDof*>& getAllGroups(void) const;
const GroupOfDof& getGoDFromElement(const MElement& element) const;
@@ -209,8 +209,8 @@ inline size_t FunctionSpace::groupNumber(void) const{
return group->size();
}
-inline const std::vector<const Dof*> FunctionSpace::getAllDofs(void) const{
- return std::vector<const Dof*>(dof->begin(), dof->end());
+inline const std::vector<Dof> FunctionSpace::getAllDofs(void) const{
+ return std::vector<Dof>(dof->begin(), dof->end());
}
inline const std::vector<GroupOfDof*>& FunctionSpace::getAllGroups(void) const{
diff --git a/Geo/MElement.h b/Geo/MElement.h
index 328c3adf612e6d6ad558cfcaff327f680c72f209..1d086a18a1f9618c203c6171ca82448acb510bb2 100644
--- a/Geo/MElement.h
+++ b/Geo/MElement.h
@@ -347,6 +347,7 @@ class MElement
// info for specific IO formats (returning 0 means that the element
// is not implemented in that format)
virtual int getTypeForMSH() const { return 0; }
+ virtual int getLowOrderTypeForMSH() const { return 0; }
virtual int getTypeForUNV() const { return 0; }
virtual int getTypeForVTK() const { return 0; }
virtual const char *getStringForPOS() const { return 0; }
diff --git a/Geo/MLine.h b/Geo/MLine.h
index 1cb90314eb26197eb84d8ad77358abf465098a8b..28e341d59938f35f9137a702eed25e8359a4f8b9 100644
--- a/Geo/MLine.h
+++ b/Geo/MLine.h
@@ -65,6 +65,7 @@ class MLine : public MElement {
virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n){}
virtual int getType() const { return TYPE_LIN; }
virtual int getTypeForMSH() const { return MSH_LIN_2; }
+ virtual int getLowOrderTypeForMSH() const { return MSH_LIN_2; }
virtual int getTypeForUNV() const { return 21; } // linear beam
virtual int getTypeForVTK() const { return 3; }
virtual const char *getStringForPOS() const { return "SL"; }
diff --git a/Geo/MQuadrangle.h b/Geo/MQuadrangle.h
index 14e3bf691e6e20196ff77ee97c4e513d0f4a0479..28cdb5d515a00db2538a58f02b39877c2ed3f277 100644
--- a/Geo/MQuadrangle.h
+++ b/Geo/MQuadrangle.h
@@ -111,6 +111,7 @@ class MQuadrangle : public MElement {
}
virtual int getType() const { return TYPE_QUA; }
virtual int getTypeForMSH() const { return MSH_QUA_4; }
+ virtual int getLowOrderTypeForMSH() const { return MSH_QUA_4; }
virtual int getTypeForUNV() const { return 94; } // thin shell linear quadrilateral
virtual int getTypeForVTK() const { return 9; }
virtual const char *getStringForPOS() const { return "SQ"; }
diff --git a/Geo/MTetrahedron.h b/Geo/MTetrahedron.h
index 920a88bc2fcb1851e908e3a3417a3a3b9fc70a6f..f09ef8a8ee1c49ce76bbb663dde78c219ab0e96b 100644
--- a/Geo/MTetrahedron.h
+++ b/Geo/MTetrahedron.h
@@ -100,6 +100,7 @@ class MTetrahedron : public MElement {
}
virtual int getType() const { return TYPE_TET; }
virtual int getTypeForMSH() const { return MSH_TET_4; }
+ virtual int getLowOrderTypeForMSH() const { return MSH_TET_4; }
virtual int getTypeForUNV() const { return 111; } // solid linear tetrahedron
virtual int getTypeForVTK() const { return 10; }
virtual const char *getStringForPOS() const { return "SS"; }
diff --git a/Geo/MTriangle.h b/Geo/MTriangle.h
index eaab6370c4cd3b555b3e12b267759711c6cf3c35..345fe108368fe7b634beb7c00b3d06aab7e0ca83 100644
--- a/Geo/MTriangle.h
+++ b/Geo/MTriangle.h
@@ -114,6 +114,7 @@ class MTriangle : public MElement {
}
virtual int getType() const { return TYPE_TRI; }
virtual int getTypeForMSH() const { return MSH_TRI_3; }
+ virtual int getLowOrderTypeForMSH() const { return MSH_TRI_3; }
virtual int getTypeForUNV() const { return 91; } // thin shell linear triangle
virtual int getTypeForVTK() const { return 5; }
virtual const char *getStringForPOS() const { return "ST"; }
diff --git a/Numeric/BasisFactory.cpp b/Numeric/BasisFactory.cpp
index 60f8ea8aa10ea9847393232f2d583a63621fdff2..6f548388ee5bb6df3e41663f70a86fa56521342a 100644
--- a/Numeric/BasisFactory.cpp
+++ b/Numeric/BasisFactory.cpp
@@ -46,9 +46,17 @@ const nodalBasis* BasisFactory::getNodalBasis(int tag)
}
// FIXME: check if already exists to deallocate if necessary
- fs.insert(std::make_pair(tag, F));
+ std::pair<std::map<int, nodalBasis*>::const_iterator, bool> inserted;
- return F;
+ #pragma omp critical
+ {
+ inserted = fs.insert(std::make_pair(tag, F));
+
+ if (!inserted.second)
+ delete F;
+ }
+
+ return inserted.first->second;
}
const JacobianBasis* BasisFactory::getJacobianBasis(int tag)
diff --git a/Numeric/fullMatrix.h b/Numeric/fullMatrix.h
index 420e6e4df5e61da96123e3e298b55e5c89d7318c..8e1f8edb0dd3534347e8f41fd82631565ada6f0e 100644
--- a/Numeric/fullMatrix.h
+++ b/Numeric/fullMatrix.h
@@ -745,6 +745,25 @@ class fullMatrix
printf("};\n");
}
+
+ void printMatlab(const std::string name = "A") const
+ {
+ int ni = size1();
+ int nj = size2();
+
+ printf("%s = [", name.c_str());
+
+ for(int I = 0; I < ni; I++){
+ for(int J = 0; J < nj; J++){
+ printf("%+e,", (*this)(I, J));
+ }
+
+ printf(";");
+ }
+
+ printf("]\n");
+ }
+
// specific functions for dgshell
void mult_naiveBlock(const fullMatrix<scalar> &b, const int ncol, const int fcol,
const int alpha, const int beta, fullVector<scalar> &c,