Skip to content
Snippets Groups Projects
Commit 0196d1c7 authored by Christophe Geuzaine's avatar Christophe Geuzaine
Browse files

coverity fixes

parent d837e8e4
Branches
Tags
No related merge requests found
...@@ -483,6 +483,10 @@ static bool skipVertex(GVertex *gv) ...@@ -483,6 +483,10 @@ static bool skipVertex(GVertex *gv)
int GModel::writeGEO(const std::string &name, bool printLabels, bool onlyPhysicals) int GModel::writeGEO(const std::string &name, bool printLabels, bool onlyPhysicals)
{ {
FILE *fp = Fopen(name.c_str(), "w"); FILE *fp = Fopen(name.c_str(), "w");
if(!fp){
Msg::Error("Could not open file '%s'", name.c_str());
return 0;
}
std::map<double, std::string> meshSizeParameters; std::map<double, std::string> meshSizeParameters;
int cpt = 0; int cpt = 0;
...@@ -533,6 +537,6 @@ int GModel::writeGEO(const std::string &name, bool printLabels, bool onlyPhysica ...@@ -533,6 +537,6 @@ int GModel::writeGEO(const std::string &name, bool printLabels, bool onlyPhysica
if(getFields()->getBackgroundField() > 0) if(getFields()->getBackgroundField() > 0)
fprintf(fp, "Background Field = %i;\n", getFields()->getBackgroundField()); fprintf(fp, "Background Field = %i;\n", getFields()->getBackgroundField());
if(fp) fclose(fp); fclose(fp);
return 1; return 1;
} }
...@@ -60,13 +60,21 @@ void normalizeAngle(double &angle) ...@@ -60,13 +60,21 @@ void normalizeAngle(double &angle)
void backgroundMesh2D::create_face_mesh() void backgroundMesh2D::create_face_mesh()
{ {
quadsToTriangles(dynamic_cast<GFace*>(gf), 100000); GFace *face = dynamic_cast<GFace*>(gf);
if(!face){
Msg::Error("Entity is not a face in background mesh");
return;
}
quadsToTriangles(face, 100000);
// storing the initial mesh from GFace // storing the initial mesh from GFace
tempTR.clear(); tempTR.clear();
GFace *face = dynamic_cast<GFace*>(gf);
for(unsigned int i = 0; i < face->triangles.size(); i++) for(unsigned int i = 0; i < face->triangles.size(); i++)
tempTR.push_back(new MTriangle(face->triangles[i]->getVertex(0), face->triangles[i]->getVertex(1), face->triangles[i]->getVertex(2))); tempTR.push_back(new MTriangle(face->triangles[i]->getVertex(0),
face->triangles[i]->getVertex(1),
face->triangles[i]->getVertex(2)));
// avoid computing curvatures on the fly : only on the // avoid computing curvatures on the fly : only on the
// BGM computes once curvatures at each node // BGM computes once curvatures at each node
...@@ -121,19 +129,23 @@ void backgroundMesh2D::reset(bool erase_2D3D) ...@@ -121,19 +129,23 @@ void backgroundMesh2D::reset(bool erase_2D3D)
} }
} }
void backgroundMesh2D::unset()
void backgroundMesh2D::unset(){ {
for (unsigned int i = 0; i < vertices.size(); i++) delete vertices[i]; for (unsigned int i = 0; i < vertices.size(); i++) delete vertices[i];
for (unsigned int i = 0; i < getNumMeshElements(); i++) delete elements[i]; for (unsigned int i = 0; i < getNumMeshElements(); i++) delete elements[i];
if (octree)delete octree; if (octree)delete octree;
octree=NULL; octree=NULL;
} }
void backgroundMesh2D::create_mesh_copy()
void backgroundMesh2D::create_mesh_copy(){ {
// TODO: useful to extend it to other elements ??? // TODO: useful to extend it to other elements ???
//std::set<SPoint2> myBCNodes; //std::set<SPoint2> myBCNodes;
GFace *face = dynamic_cast<GFace*>(gf); GFace *face = dynamic_cast<GFace*>(gf);
if(!face){
Msg::Error("Entity is not a face in background mesh");
return;
}
for (unsigned int i = 0; i < face->triangles.size(); i++){ for (unsigned int i = 0; i < face->triangles.size(); i++){
MTriangle *e = face->triangles[i]; MTriangle *e = face->triangles[i];
MVertex *news[3]; MVertex *news[3];
...@@ -158,8 +170,14 @@ void backgroundMesh2D::create_mesh_copy(){ ...@@ -158,8 +170,14 @@ void backgroundMesh2D::create_mesh_copy(){
} }
GPoint backgroundMesh2D::get_GPoint_from_MVertex(const MVertex *v)const{ GPoint backgroundMesh2D::get_GPoint_from_MVertex(const MVertex *v)const
return dynamic_cast<GFace*>(gf)->point(SPoint2(v->x(),v->y())); {
GFace *face = dynamic_cast<GFace*>(gf);
if(!face){
Msg::Error("Entity is not a face in background mesh");
return GPoint();
}
return face->point(SPoint2(v->x(),v->y()));
} }
...@@ -171,6 +189,9 @@ backgroundMesh2D::backgroundMesh2D(GFace *_gf, bool erase_2D3D):BGMBase(2,_gf),s ...@@ -171,6 +189,9 @@ backgroundMesh2D::backgroundMesh2D(GFace *_gf, bool erase_2D3D):BGMBase(2,_gf),s
// now, the new mesh has been copied in local in backgroundMesh2D, deleting the mesh // now, the new mesh has been copied in local in backgroundMesh2D, deleting the mesh
// from GFace, back to the previous one ! // from GFace, back to the previous one !
GFace *face = dynamic_cast<GFace*>(gf); GFace *face = dynamic_cast<GFace*>(gf);
if(!face)
Msg::Error("Entity is not a face in background mesh");
else
face->triangles = tempTR; face->triangles = tempTR;
} }
} }
...@@ -181,7 +202,10 @@ backgroundMesh2D::~backgroundMesh2D() ...@@ -181,7 +202,10 @@ backgroundMesh2D::~backgroundMesh2D()
} }
void backgroundMesh2D::propagateValues(DoubleStorageType &dirichlet, simpleFunction<double> &eval_diffusivity, bool in_parametric_plane){ void backgroundMesh2D::propagateValues(DoubleStorageType &dirichlet,
simpleFunction<double> &eval_diffusivity,
bool in_parametric_plane)
{
#if defined(HAVE_SOLVER) #if defined(HAVE_SOLVER)
linearSystem<double> *_lsys = 0; linearSystem<double> *_lsys = 0;
#if defined(HAVE_PETSC) && !defined(HAVE_TAUCS) #if defined(HAVE_PETSC) && !defined(HAVE_TAUCS)
...@@ -207,12 +231,16 @@ void backgroundMesh2D::propagateValues(DoubleStorageType &dirichlet, simpleFunct ...@@ -207,12 +231,16 @@ void backgroundMesh2D::propagateValues(DoubleStorageType &dirichlet, simpleFunct
// Number vertices // Number vertices
std::set<MVertex*> vs; std::set<MVertex*> vs;
GFace *face = dynamic_cast<GFace*>(gf); GFace *face = dynamic_cast<GFace*>(gf);
if(!face){
Msg::Error("Entity is not a face in background mesh");
delete _lsys;
return;
}
for (unsigned int k = 0; k < face->triangles.size(); k++) for (unsigned int k = 0; k < face->triangles.size(); k++)
for (int j=0;j<3;j++)vs.insert(face->triangles[k]->getVertex(j)); for (int j=0;j<3;j++)vs.insert(face->triangles[k]->getVertex(j));
for (unsigned int k = 0; k < face->quadrangles.size(); k++) for (unsigned int k = 0; k < face->quadrangles.size(); k++)
for (int j=0;j<4;j++)vs.insert(face->quadrangles[k]->getVertex(j)); for (int j=0;j<4;j++)vs.insert(face->quadrangles[k]->getVertex(j));
std::map<MVertex*,SPoint3> theMap; std::map<MVertex*,SPoint3> theMap;
if ( in_parametric_plane) { if ( in_parametric_plane) {
for (std::set<MVertex*>::iterator it = vs.begin(); it != vs.end(); ++it){ for (std::set<MVertex*>::iterator it = vs.begin(); it != vs.end(); ++it){
...@@ -254,10 +282,15 @@ void backgroundMesh2D::propagateValues(DoubleStorageType &dirichlet, simpleFunct ...@@ -254,10 +282,15 @@ void backgroundMesh2D::propagateValues(DoubleStorageType &dirichlet, simpleFunct
#endif #endif
} }
void backgroundMesh2D::computeSizeField()
void backgroundMesh2D::computeSizeField(){ {
GFace *face = dynamic_cast<GFace*>(gf); GFace *face = dynamic_cast<GFace*>(gf);
list<GEdge*> e;// = face->edges(); if(!face){
Msg::Error("Entity is not a face in background mesh");
return;
}
list<GEdge*> e;
replaceMeshCompound(face, e); replaceMeshCompound(face, e);
list<GEdge*>::const_iterator it = e.begin(); list<GEdge*>::const_iterator it = e.begin();
DoubleStorageType sizes; DoubleStorageType sizes;
...@@ -303,7 +336,8 @@ inline double myAngle (const SVector3 &a, const SVector3 &b, const SVector3 &d){ ...@@ -303,7 +336,8 @@ inline double myAngle (const SVector3 &a, const SVector3 &b, const SVector3 &d){
} }
void backgroundMesh2D::updateSizes(){ void backgroundMesh2D::updateSizes()
{
DoubleStorageType::iterator itv = sizeField.begin(); DoubleStorageType::iterator itv = sizeField.begin();
for ( ; itv != sizeField.end(); ++itv){ for ( ; itv != sizeField.end(); ++itv){
SPoint2 p; SPoint2 p;
...@@ -319,6 +353,10 @@ void backgroundMesh2D::updateSizes(){ ...@@ -319,6 +353,10 @@ void backgroundMesh2D::updateSizes(){
} }
else{ else{
GFace *face = dynamic_cast<GFace*>(gf); GFace *face = dynamic_cast<GFace*>(gf);
if(!face){
Msg::Error("Entity is not a face in background mesh");
return;
}
reparamMeshVertexOnFace(v, face, p); reparamMeshVertexOnFace(v, face, p);
lc = sizeFactor * BGM_MeshSize(face, p.x(), p.y(), v->x(), v->y(), v->z()); lc = sizeFactor * BGM_MeshSize(face, p.x(), p.y(), v->x(), v->y(), v->z());
} }
...@@ -363,6 +401,9 @@ frameFieldBackgroundMesh2D::frameFieldBackgroundMesh2D(GFace *_gf):backgroundMes ...@@ -363,6 +401,9 @@ frameFieldBackgroundMesh2D::frameFieldBackgroundMesh2D(GFace *_gf):backgroundMes
// now, the new mesh has been copied in local in backgroundMesh2D, deleting the mesh // now, the new mesh has been copied in local in backgroundMesh2D, deleting the mesh
// from GFace, back to the previous one ! // from GFace, back to the previous one !
GFace *face = dynamic_cast<GFace*>(gf); GFace *face = dynamic_cast<GFace*>(gf);
if(!face)
Msg::Error("Entity is not a face in background mesh");
else
face->triangles = tempTR; face->triangles = tempTR;
} }
...@@ -437,6 +478,11 @@ void frameFieldBackgroundMesh2D::computeCrossField(simpleFunction<double> &eval_ ...@@ -437,6 +478,11 @@ void frameFieldBackgroundMesh2D::computeCrossField(simpleFunction<double> &eval_
list<GEdge*> e; list<GEdge*> e;
GFace *face = dynamic_cast<GFace*>(gf); GFace *face = dynamic_cast<GFace*>(gf);
if(!face){
Msg::Error("Entity is not a face in background mesh");
return;
}
replaceMeshCompound(face, e); replaceMeshCompound(face, e);
list<GEdge*>::const_iterator it = e.begin(); list<GEdge*>::const_iterator it = e.begin();
...@@ -523,6 +569,10 @@ void frameFieldBackgroundMesh2D::eval_crossfield(MVertex *vert, STensor3 &cf) ...@@ -523,6 +569,10 @@ void frameFieldBackgroundMesh2D::eval_crossfield(MVertex *vert, STensor3 &cf)
{ {
SPoint2 parampoint; SPoint2 parampoint;
GFace *face = dynamic_cast<GFace*>(gf); GFace *face = dynamic_cast<GFace*>(gf);
if(!face){
Msg::Error("Entity is not a face in background mesh");
return;
}
reparamMeshVertexOnFace(vert, face, parampoint); reparamMeshVertexOnFace(vert, face, parampoint);
return eval_crossfield(parampoint[0], parampoint[1], cf); return eval_crossfield(parampoint[0], parampoint[1], cf);
} }
...@@ -596,10 +646,16 @@ void frameFieldBackgroundMesh2D::exportCrossField(const std::string &filename) ...@@ -596,10 +646,16 @@ void frameFieldBackgroundMesh2D::exportCrossField(const std::string &filename)
} }
// returns the cross field as a pair of othogonal vectors (NOT in parametric coordinates, but real 3D coordinates) // returns the cross field as a pair of othogonal vectors (NOT in parametric coordinates, but real 3D coordinates)
Pair<SVector3, SVector3> frameFieldBackgroundMesh2D::compute_crossfield_directions(double u,double v, double angle_current) Pair<SVector3, SVector3> frameFieldBackgroundMesh2D::compute_crossfield_directions(double u, double v,
double angle_current)
{ {
// get the unit normal at that point // get the unit normal at that point
GFace *face = dynamic_cast<GFace*>(gf); GFace *face = dynamic_cast<GFace*>(gf);
if(!face){
Msg::Error("Entity is not a face in background mesh");
return Pair<SVector3,SVector3>(SVector3(), SVector3());
}
Pair<SVector3, SVector3> der = face->firstDer(SPoint2(u,v)); Pair<SVector3, SVector3> der = face->firstDer(SPoint2(u,v));
SVector3 s1 = der.first(); SVector3 s1 = der.first();
SVector3 s2 = der.second(); SVector3 s2 = der.second();
...@@ -623,7 +679,6 @@ Pair<SVector3, SVector3> frameFieldBackgroundMesh2D::compute_crossfield_directio ...@@ -623,7 +679,6 @@ Pair<SVector3, SVector3> frameFieldBackgroundMesh2D::compute_crossfield_directio
bool frameFieldBackgroundMesh2D::compute_RK_infos(double u,double v, double x, double y, double z, RK_form &infos) bool frameFieldBackgroundMesh2D::compute_RK_infos(double u,double v, double x, double y, double z, RK_form &infos)
{ {
// check if point is in domain // check if point is in domain
if (!inDomain(u,v)) return false; if (!inDomain(u,v)) return false;
...@@ -635,6 +690,11 @@ bool frameFieldBackgroundMesh2D::compute_RK_infos(double u,double v, double x, d ...@@ -635,6 +690,11 @@ bool frameFieldBackgroundMesh2D::compute_RK_infos(double u,double v, double x, d
// get the unit normal at that point // get the unit normal at that point
GFace *face = dynamic_cast<GFace*>(gf); GFace *face = dynamic_cast<GFace*>(gf);
if(!face){
Msg::Error("Entity is not a face in background mesh");
return false;
}
Pair<SVector3, SVector3> der = face->firstDer(SPoint2(u,v)); Pair<SVector3, SVector3> der = face->firstDer(SPoint2(u,v));
SVector3 s1 = der.first(); SVector3 s1 = der.first();
SVector3 s2 = der.second(); SVector3 s2 = der.second();
......
...@@ -27,12 +27,12 @@ using namespace std; ...@@ -27,12 +27,12 @@ using namespace std;
*/ */
class backgroundMesh2D : public BGMBase { class backgroundMesh2D : public BGMBase {
protected: protected:
virtual void propagateValues(DoubleStorageType &dirichlet, simpleFunction<double> &eval_diffusivity, bool in_parametric_plane = false); virtual void propagateValues(DoubleStorageType &dirichlet,
simpleFunction<double> &eval_diffusivity,
bool in_parametric_plane = false);
virtual void computeSizeField(); virtual void computeSizeField();
virtual inline unsigned int getNumMeshElements()const{return elements.size();} virtual inline unsigned int getNumMeshElements()const{return elements.size();}
virtual const MElement* getElement(unsigned int i)const; virtual const MElement* getElement(unsigned int i)const;
virtual GPoint get_GPoint_from_MVertex(const MVertex *) const; virtual GPoint get_GPoint_from_MVertex(const MVertex *) const;
// only 2D: // only 2D:
...@@ -43,20 +43,15 @@ class backgroundMesh2D : public BGMBase { ...@@ -43,20 +43,15 @@ class backgroundMesh2D : public BGMBase {
void create_face_mesh(); void create_face_mesh();
double sizeFactor; double sizeFactor;
std::vector<MTriangle*> tempTR; std::vector<MTriangle*> tempTR;
vector<MElement*> elements; vector<MElement*> elements;
vector<MVertex*> vertices; vector<MVertex*> vertices;
map<MVertex*,MVertex*> _3Dto2D; map<MVertex*,MVertex*> _3Dto2D;
map<MVertex*,MVertex*> _2Dto3D; map<MVertex*,MVertex*> _2Dto3D;
public: public:
backgroundMesh2D(GFace *, bool erase_2D3D=true); backgroundMesh2D(GFace *, bool erase_2D3D=true);
virtual ~backgroundMesh2D(); virtual ~backgroundMesh2D();
virtual MElementOctree* getOctree(); virtual MElementOctree* getOctree();
// TODO: only 2D // TODO: only 2D
...@@ -66,7 +61,6 @@ class backgroundMesh2D : public BGMBase { ...@@ -66,7 +61,6 @@ class backgroundMesh2D : public BGMBase {
// not used !!!! TODO !!! // not used !!!! TODO !!!
void setSizeFactor (double s) {sizeFactor = s;} void setSizeFactor (double s) {sizeFactor = s;}
virtual vector<MVertex*>::iterator beginvertices(){ return vertices.begin(); } virtual vector<MVertex*>::iterator beginvertices(){ return vertices.begin(); }
virtual vector<MVertex*>::iterator endvertices(){ return vertices.end(); } virtual vector<MVertex*>::iterator endvertices(){ return vertices.end(); }
virtual vector<MElement*>::iterator beginelements(){ return elements.begin(); } virtual vector<MElement*>::iterator beginelements(){ return elements.begin(); }
...@@ -75,26 +69,13 @@ class backgroundMesh2D : public BGMBase { ...@@ -75,26 +69,13 @@ class backgroundMesh2D : public BGMBase {
virtual vector<MVertex*>::const_iterator endvertices() const { return vertices.end(); } virtual vector<MVertex*>::const_iterator endvertices() const { return vertices.end(); }
virtual vector<MElement*>::const_iterator beginelements() const { return elements.begin(); } virtual vector<MElement*>::const_iterator beginelements() const { return elements.begin(); }
virtual vector<MElement*>::const_iterator endelements() const { return elements.end(); } virtual vector<MElement*>::const_iterator endelements() const { return elements.end(); }
}; };
class RK_form{ // informations for RK at one point class RK_form{ // informations for RK at one point
public: public:
RK_form(){ RK_form() : angle(0.), localsize(0.)
{
} }
// RK_form(const RK_form &other){
// t1 = other.t1;
// t2 = other.t2;
// h.first = other.h.first;
// h.second = other.h.second;
// paramh.first = other.paramh.first;
// paramh.second = other.paramh.second;
// paramt1 = other.paramt1;
// paramt2 = other.paramt2;
// angle = other.angle;
// }
SMetric3 metricField; SMetric3 metricField;
SVector3 t1, t2;// 3D cross field directions SVector3 t1, t2;// 3D cross field directions
SVector3 normal;// 3D cross field directions SVector3 normal;// 3D cross field directions
...@@ -104,14 +85,11 @@ class RK_form{// informations for RK at one point ...@@ -104,14 +85,11 @@ class RK_form{// informations for RK at one point
double angle, localsize; double angle, localsize;
}; };
class frameFieldBackgroundMesh2D : public backgroundMesh2D { class frameFieldBackgroundMesh2D : public backgroundMesh2D {
private: private:
// specification for cross field // specification for cross field
DoubleStorageType angles; // an attached angle DoubleStorageType angles; // an attached angle
DoubleStorageType smoothness; DoubleStorageType smoothness;
void computeCrossField(simpleFunction<double> &eval_diffusivity); void computeCrossField(simpleFunction<double> &eval_diffusivity);
void computeSmoothness(); void computeSmoothness();
...@@ -131,10 +109,14 @@ class frameFieldBackgroundMesh2D : public backgroundMesh2D{ ...@@ -131,10 +109,14 @@ class frameFieldBackgroundMesh2D : public backgroundMesh2D{
void eval_crossfield(MVertex *vert, STensor3 &cf); void eval_crossfield(MVertex *vert, STensor3 &cf);
void exportCrossField(const string &filename); void exportCrossField(const string &filename);
Pair<SVector3, SVector3> compute_crossfield_directions(double u,double v, double angle_current); Pair<SVector3, SVector3> compute_crossfield_directions(double u,double v,
double angle_current);
bool compute_RK_infos(double u,double v, double x, double y, double z,RK_form &infos); bool compute_RK_infos(double u,double v, double x, double y, double z,RK_form &infos);
inline void exportSmoothness(const string &filename) const{export_scalar(filename,smoothness);}; inline void exportSmoothness(const string &filename) const
{
export_scalar(filename,smoothness);
}
}; };
#endif #endif
...@@ -61,6 +61,10 @@ void backgroundMesh3D::computeSizeField() ...@@ -61,6 +61,10 @@ void backgroundMesh3D::computeSizeField()
// fills dirichlet BC // fills dirichlet BC
GRegion *region = dynamic_cast<GRegion*>(gf); GRegion *region = dynamic_cast<GRegion*>(gf);
if(!region){
Msg::Error("Entity is not a region in background mesh");
return;
}
list<GFace*> faces = region->faces(); list<GFace*> faces = region->faces();
MVertex *v; MVertex *v;
MElement *e; MElement *e;
...@@ -68,6 +72,10 @@ void backgroundMesh3D::computeSizeField() ...@@ -68,6 +72,10 @@ void backgroundMesh3D::computeSizeField()
for (list<GFace*>::iterator it=faces.begin();it!=faces.end();it++){// for all GFace for (list<GFace*>::iterator it=faces.begin();it!=faces.end();it++){// for all GFace
GFace *face = *it; GFace *face = *it;
frameFieldBackgroundMesh2D *bgm2d = dynamic_cast<frameFieldBackgroundMesh2D*>(BGMManager::get(face)); frameFieldBackgroundMesh2D *bgm2d = dynamic_cast<frameFieldBackgroundMesh2D*>(BGMManager::get(face));
if(!bgm2d){
Msg::Error("Background mesh is not a 2D frame field");
return;
}
for (unsigned int i=0;i<face->getNumMeshElements();i++){// for all elements for (unsigned int i=0;i<face->getNumMeshElements();i++){// for all elements
e = face->getMeshElement(i); e = face->getMeshElement(i);
if (e->getDim()!=2) continue;// of dim=2 if (e->getDim()!=2) continue;// of dim=2
...@@ -94,6 +102,10 @@ void backgroundMesh3D::propagateValues(DoubleStorageType &dirichlet, ...@@ -94,6 +102,10 @@ void backgroundMesh3D::propagateValues(DoubleStorageType &dirichlet,
{ {
// same as Size_field::solve() // same as Size_field::solve()
GRegion *gr = dynamic_cast<GRegion*>(gf); GRegion *gr = dynamic_cast<GRegion*>(gf);
if(!gr){
Msg::Error("Entity is not a region in background mesh");
return;
}
#if defined(HAVE_PETSC) #if defined(HAVE_PETSC)
linearSystem<double>* system = 0; linearSystem<double>* system = 0;
...@@ -195,6 +207,10 @@ MElementOctree* backgroundMesh3D::getOctree() ...@@ -195,6 +207,10 @@ MElementOctree* backgroundMesh3D::getOctree()
{ {
if(!octree){ if(!octree){
GRegion *gr = dynamic_cast<GRegion*>(gf); GRegion *gr = dynamic_cast<GRegion*>(gf);
if(!gr){
Msg::Error("Entity is not a region in background mesh");
return 0;
}
Msg::Debug("Rebuilding BackgroundMesh element octree"); Msg::Debug("Rebuilding BackgroundMesh element octree");
std::vector<MElement*> copy;// !!! std::vector<MElement*> copy;// !!!
for (std::vector<MTetrahedron*>::iterator it = gr->tetrahedra.begin(); for (std::vector<MTetrahedron*>::iterator it = gr->tetrahedra.begin();
...@@ -578,6 +594,10 @@ void frameFieldBackgroundMesh3D::initiate_crossfield() ...@@ -578,6 +594,10 @@ void frameFieldBackgroundMesh3D::initiate_crossfield()
// first, for boundaries : // first, for boundaries :
GRegion *gr = dynamic_cast<GRegion*>(gf); GRegion *gr = dynamic_cast<GRegion*>(gf);
if(!gr){
Msg::Error("Entity is not a region in background mesh");
return;
}
list<GFace*> faces = gr->faces(); list<GFace*> faces = gr->faces();
// here, not using the gm2D since we are interested by the new 2D vertices, // here, not using the gm2D since we are interested by the new 2D vertices,
// not the old (now erased) ones... alternative would be to reset the 2DBGM... // not the old (now erased) ones... alternative would be to reset the 2DBGM...
...@@ -585,6 +605,10 @@ void frameFieldBackgroundMesh3D::initiate_crossfield() ...@@ -585,6 +605,10 @@ void frameFieldBackgroundMesh3D::initiate_crossfield()
GFace *face = *it; GFace *face = *it;
frameFieldBackgroundMesh2D *bgm2d = frameFieldBackgroundMesh2D *bgm2d =
dynamic_cast<frameFieldBackgroundMesh2D*>(BGMManager::get(face)); dynamic_cast<frameFieldBackgroundMesh2D*>(BGMManager::get(face));
if(!bgm2d){
Msg::Error("Background mesh is not a 2D frame field");
return;
}
// initializing the vertices on the faces // initializing the vertices on the faces
for (unsigned int i=0;i<face->getNumMeshElements();i++){// for all elements for (unsigned int i=0;i<face->getNumMeshElements();i++){// for all elements
MElement *e = face->getMeshElement(i); MElement *e = face->getMeshElement(i);
...@@ -721,6 +745,10 @@ double frameFieldBackgroundMesh3D::get_vectorial_smoothness(const int idir, ...@@ -721,6 +745,10 @@ double frameFieldBackgroundMesh3D::get_vectorial_smoothness(const int idir,
void frameFieldBackgroundMesh3D::build_vertex_to_element_table() void frameFieldBackgroundMesh3D::build_vertex_to_element_table()
{ {
GRegion *gr = dynamic_cast<GRegion*>(gf); GRegion *gr = dynamic_cast<GRegion*>(gf);
if(!gr){
Msg::Error("Entity is not a region in background mesh");
return;
}
MElement *e; MElement *e;
MVertex *v; MVertex *v;
...@@ -746,12 +774,20 @@ void frameFieldBackgroundMesh3D::build_vertex_to_element_table() ...@@ -746,12 +774,20 @@ void frameFieldBackgroundMesh3D::build_vertex_to_element_table()
const MElement* backgroundMesh3D::getElement(unsigned int i)const const MElement* backgroundMesh3D::getElement(unsigned int i)const
{ {
GRegion *gr = dynamic_cast<GRegion*>(gf); GRegion *gr = dynamic_cast<GRegion*>(gf);
if(!gr){
Msg::Error("Entity is not a region in background mesh");
return 0;
}
return gr->getMeshElement(i); return gr->getMeshElement(i);
} }
unsigned int backgroundMesh3D::getNumMeshElements()const unsigned int backgroundMesh3D::getNumMeshElements()const
{ {
GRegion *gr = dynamic_cast<GRegion*>(gf); GRegion *gr = dynamic_cast<GRegion*>(gf);
if(!gr){
Msg::Error("Entity is not a region in background mesh");
return 0;
}
return gr->getNumMeshElements(); return gr->getNumMeshElements();
} }
......
...@@ -106,23 +106,29 @@ private: ...@@ -106,23 +106,29 @@ private:
SVector3 *vec; SVector3 *vec;
bool ok; bool ok;
public: public:
Wrapper3D():ok(true){}; Wrapper3D()
Wrapper3D(MVertex* _i,MVertex* _p):individual(_i), parent(_p),ok(true){}; : individual(0), parent(0), size(0), cf(0), vec(0), ok(true)
~Wrapper3D(){}; {
void set_ok(bool b){ok=b;}; }
void set_individual(MVertex *vertex){individual=vertex;}; Wrapper3D(MVertex* _i,MVertex* _p)
void set_parent(MVertex *vertex){parent=vertex;}; : individual(_i), parent(_p), size(0), cf(0), vec(0), ok(true)
void set_size(double *h){size=h;}; {
void set_crossfield(STensor3 *_cf){cf=_cf;}; }
void set_direction(SVector3 *_v){vec=_v;}; ~Wrapper3D(){}
bool get_ok(){return ok;}; void set_ok(bool b){ ok = b; }
void set_bgm(frameFieldBackgroundMesh3D *bgm){bgmesh = bgm;}; void set_individual(MVertex *vertex){ individual = vertex; }
frameFieldBackgroundMesh3D * bgm(){return bgmesh;}; void set_parent(MVertex *vertex){ parent = vertex; }
MVertex* get_individual(){return individual;}; void set_size(double *h){ size = h; }
MVertex* get_parent(){return parent;}; void set_crossfield(STensor3 *_cf){ cf = _cf; }
STensor3* get_crossfield(){return cf;}; void set_direction(SVector3 *_v){ vec = _v; }
SVector3* get_direction(){return vec;}; bool get_ok(){ return ok; }
double* get_size(){return size;}; void set_bgm(frameFieldBackgroundMesh3D *bgm){ bgmesh = bgm; }
frameFieldBackgroundMesh3D * bgm(){ return bgmesh; }
MVertex* get_individual(){ return individual; }
MVertex* get_parent(){ return parent; }
STensor3* get_crossfield(){ return cf; }
SVector3* get_direction(){ return vec; }
double* get_size(){ return size; }
}; };
extern double infinity_distance_3D(const MVertex *v1,const MVertex *v2,STensor3 &cf); extern double infinity_distance_3D(const MVertex *v1,const MVertex *v2,STensor3 &cf);
...@@ -134,8 +140,8 @@ extern void fill_min_max(double x,double y,double z,double h,double *min,double ...@@ -134,8 +140,8 @@ extern void fill_min_max(double x,double y,double z,double h,double *min,double
// listOfPoints AND in RTree: larger memory footprint but less CPU time... // listOfPoints AND in RTree: larger memory footprint but less CPU time...
class smoothness_vertex_pair{ class smoothness_vertex_pair{
public: public:
smoothness_vertex_pair(){}; smoothness_vertex_pair() : rank(0.), size(0.), v(0), dir(0), layer(0) {}
~smoothness_vertex_pair(){}; ~smoothness_vertex_pair(){}
STensor3 cf; STensor3 cf;
SVector3 direction; SVector3 direction;
double rank, size; double rank, size;
...@@ -180,30 +186,30 @@ public: ...@@ -180,30 +186,30 @@ public:
class listOfPointsScalarSmoothness : public listOfPoints{ class listOfPointsScalarSmoothness : public listOfPoints{
public: public:
listOfPointsScalarSmoothness(){ }; listOfPointsScalarSmoothness(){ }
virtual ~listOfPointsScalarSmoothness() virtual ~listOfPointsScalarSmoothness()
{ {
while (!empty()) while (!empty())
erase_first(); erase_first();
}; }
virtual void insert(smoothness_vertex_pair *svp){ points.insert(svp); }; virtual void insert(smoothness_vertex_pair *svp){ points.insert(svp); }
virtual unsigned int size(){ return points.size(); }; virtual unsigned int size(){ return points.size(); }
virtual MVertex* get_first_vertex(){ return (*points.begin())->v; }; virtual MVertex* get_first_vertex(){ return (*points.begin())->v; }
virtual STensor3 get_first_crossfield(){ return (*points.begin())->cf; }; virtual STensor3 get_first_crossfield(){ return (*points.begin())->cf; }
virtual double get_first_size(){ return (*points.begin())->size; }; virtual double get_first_size(){ return (*points.begin())->size; }
virtual int get_first_layer(){ return (*points.begin())->layer; }; virtual int get_first_layer(){ return (*points.begin())->layer; }
virtual SVector3 get_first_direction() virtual SVector3 get_first_direction()
{ {
Msg::Error("listOfPointsScalarSmoothness::get_first_direction NOT applicable"); Msg::Error("listOfPointsScalarSmoothness::get_first_direction NOT applicable");
return SVector3(0.); return SVector3(0.);
}; }
virtual void erase_first() virtual void erase_first()
{ {
smoothness_vertex_pair *ptr = *(points.begin()); smoothness_vertex_pair *ptr = *(points.begin());
points.erase(points.begin()); points.erase(points.begin());
delete ptr; delete ptr;
}; }
virtual bool empty(){ return points.empty(); }; virtual bool empty(){ return points.empty(); }
protected: protected:
std::set<smoothness_vertex_pair*, compareSmoothnessVertexPairs> points; std::set<smoothness_vertex_pair*, compareSmoothnessVertexPairs> points;
...@@ -211,12 +217,12 @@ protected: ...@@ -211,12 +217,12 @@ protected:
class listOfPointsVectorialSmoothness : public listOfPointsScalarSmoothness{ class listOfPointsVectorialSmoothness : public listOfPointsScalarSmoothness{
public: public:
listOfPointsVectorialSmoothness(){}; listOfPointsVectorialSmoothness(){}
virtual ~listOfPointsVectorialSmoothness(){ virtual ~listOfPointsVectorialSmoothness(){
while (!empty()) while (!empty())
erase_first(); erase_first();
}; }
virtual SVector3 get_first_direction(){ return (*points.begin())->direction; }; virtual SVector3 get_first_direction(){ return (*points.begin())->direction; }
protected: protected:
std::set<smoothness_vertex_pair*, compareSmoothnessVertexPairs> points; std::set<smoothness_vertex_pair*, compareSmoothnessVertexPairs> points;
}; };
...@@ -228,24 +234,24 @@ public: ...@@ -228,24 +234,24 @@ public:
while (!empty()) while (!empty())
erase_first(); erase_first();
}; };
virtual void insert(smoothness_vertex_pair *svp){ points.push(svp); }; virtual void insert(smoothness_vertex_pair *svp){ points.push(svp); }
virtual unsigned int size(){ return points.size(); }; virtual unsigned int size(){ return points.size(); }
virtual MVertex* get_first_vertex(){ return (points.front())->v; }; virtual MVertex* get_first_vertex(){ return (points.front())->v; }
virtual STensor3 get_first_crossfield(){ return (points.front())->cf; }; virtual STensor3 get_first_crossfield(){ return (points.front())->cf; }
virtual double get_first_size(){ return (points.front())->size; }; virtual double get_first_size(){ return (points.front())->size; }
virtual int get_first_layer(){ return (points.front())->layer; }; virtual int get_first_layer(){ return (points.front())->layer; }
virtual SVector3 get_first_direction() virtual SVector3 get_first_direction()
{ {
Msg::Error("listOfPointsFifo::get_first_direction NOT applicable"); Msg::Error("listOfPointsFifo::get_first_direction NOT applicable");
return SVector3(0.); return SVector3(0.);
}; }
virtual void erase_first() virtual void erase_first()
{ {
smoothness_vertex_pair *ptr = points.front(); smoothness_vertex_pair *ptr = points.front();
points.pop(); points.pop();
delete ptr; delete ptr;
}; }
virtual bool empty(){ return points.empty(); }; virtual bool empty(){ return points.empty(); }
protected: protected:
std::queue<smoothness_vertex_pair*> points; std::queue<smoothness_vertex_pair*> points;
......
...@@ -110,6 +110,10 @@ void GMSH_DistancePlugin::printView(std::vector<GEntity*> _entities, ...@@ -110,6 +110,10 @@ void GMSH_DistancePlugin::printView(std::vector<GEntity*> _entities,
Msg::Info("Writing %s", _fileName.c_str()); Msg::Info("Writing %s", _fileName.c_str());
FILE *fName = Fopen(_fileName.c_str(),"w"); FILE *fName = Fopen(_fileName.c_str(),"w");
if(!fName){
Msg::Error("Could not open file '%s'", _fileName.c_str());
return;
}
fprintf(fName, "View \"distance \"{\n"); fprintf(fName, "View \"distance \"{\n");
for (unsigned int ii=0; ii<_entities.size(); ii++) { for (unsigned int ii=0; ii<_entities.size(); ii++) {
...@@ -496,6 +500,15 @@ PView *GMSH_DistancePlugin::execute(PView *v) ...@@ -496,6 +500,15 @@ PView *GMSH_DistancePlugin::execute(PView *v)
Msg::Info("Writing orthogonal.pos"); Msg::Info("Writing orthogonal.pos");
FILE *f5 = Fopen("orthogonal.pos","w"); FILE *f5 = Fopen("orthogonal.pos","w");
if(!f5){
Msg::Error("Could not open file 'orthogonal.pos'");
#if defined(HAVE_SOLVER)
delete lsys;
delete dofView;
#endif
return view;
}
fprintf(f5,"View \"orthogonal\"{\n"); fprintf(f5,"View \"orthogonal\"{\n");
for (std::vector<MElement* >::iterator it = allElems.begin(); for (std::vector<MElement* >::iterator it = allElems.begin();
it != allElems.end(); it++){ it != allElems.end(); it++){
......
...@@ -20,7 +20,7 @@ GETDP=svn ...@@ -20,7 +20,7 @@ GETDP=svn
#GMSH=2.8.5 #GMSH=2.8.5
#GETDP=2.4.4 #GETDP=2.4.4
MODELS='machines relay inductor indheat magnetometer antennas acoustic_scattering time_reversal shielding waveguides transfo_simple ddm_wave_simple bloch_periodic_waveguides magnets thermal_conduction' MODELS='machines relay inductor indheat magnetometer antennas acoustic_scattering time_reversal shielding waveguides transfo_simple ddm_wave_simple bloch_periodic_waveguides magnets thermal_conduction magnetostriction'
# get onelab models # get onelab models
mkdir /tmp/models mkdir /tmp/models
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment