Newer
Older
// Gmsh - Copyright (C) 1997-2017 C. Geuzaine, J.-F. Remacle
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@onelab.info>.
#include "GmshConfig.h"
Christophe Geuzaine
committed
#include "GmshMessage.h"
#include "Numeric.h"
#include "Context.h"
#include "MLine.h"
#include "MTriangle.h"
#include "MQuadrangle.h"
#include "MTetrahedron.h"
#include "MHexahedron.h"
#include "MPrism.h"
#include "MPyramid.h"
#include "BoundaryLayers.h"
#include "GFaceCompound.h"
#include "Field.h"
#include "yamakawa.h"
Paul-Emile Bernard
committed
#include "pointInsertion.h"
#if defined(HAVE_OPTHOM)
#include "OptHomRun.h"
#include "OptHomElastic.h"
#endif
Christophe Geuzaine
committed
#include "PView.h"
#include "PViewData.h"
Christophe Geuzaine
committed
#endif
class TEST_IF_MESH_IS_COMPATIBLE_WITH_EMBEDDED_ENTITIES {
public:
void operator () (GRegion *gr) {
std::list<GEdge*> e = gr->embeddedEdges();
std::list<GFace*> f = gr->embeddedFaces();
if (e.empty() && f.empty())return;
std::map<MEdge,GEdge*,Less_Edge> edges;
std::map<MFace,GFace*,Less_Face> faces;
std::list<GEdge*>::iterator it = e.begin();
std::list<GFace*>::iterator itf = f.begin();
for ( ; it != e.end() ; ++it){
for (unsigned int i=0;i<(*it)->lines.size(); ++i){
if (distance ((*it)->lines[i]->getVertex(0),(*it)->lines[i]->getVertex(1)) > 1.e-12)
edges.insert(std::make_pair(MEdge((*it)->lines[i]->getVertex(0),(*it)->lines[i]->getVertex(1)),*it));
}
}
for ( ; itf != f.end() ; ++itf){
for (unsigned int i=0;i<(*itf)->triangles.size(); ++i){
faces.insert(std::make_pair(MFace((*itf)->triangles[i]->getVertex(0),(*itf)->triangles[i]->getVertex(1),(*itf)->triangles[i]->getVertex(2)),*itf));
}
}
Msg::Info ("Searching for %d embedded mesh edges and %d embedded mesh faces in region %d", edges.size(), faces.size(), gr->tag());
for (unsigned int k=0;k<gr->getNumMeshElements();k++){
for (int j=0;j<gr->getMeshElement(k)->getNumEdges();j++){
edges.erase (gr->getMeshElement(k)->getEdge(j));
}
for (int j=0;j<gr->getMeshElement(k)->getNumFaces();j++){
faces.erase (gr->getMeshElement(k)->getFace(j));
}
}
if (edges.empty() && faces.empty()) {
Msg::Info ("All embedded edges and faces are present in the final mesh");
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
}
if (edges.size()) {
char name[256];
sprintf(name,"missingEdgesOnRegion%d.pos",gr->tag());
Msg::Error("Region %d : %d mesh edges that should be embedded are missing in the final mesh",gr->tag(), (int)edges.size());
Msg::Error("Saving the missing edges in file %s",name);
FILE *f = fopen(name,"w");
fprintf(f,"View \" \" {\n");
for (std::map<MEdge,GEdge*,Less_Edge>::iterator it = edges.begin() ; it != edges.end(); ++it){
MVertex *v1 = it->first.getVertex(0);
MVertex *v2 = it->first.getVertex(1);
fprintf(f,"SL(%g,%g,%g,%g,%g,%g){%d,%d};\n",v1->x(),v1->y(),v1->z(),v2->x(),v2->y(),v2->z(), it->second->tag(),it->second->tag());
}
fprintf(f,"};\n");
fclose(f);
}
if (faces.size()) {
char name[256];
sprintf(name,"missingFacesOnRegion%d.pos",gr->tag());
Msg::Error("Region %d : %d mesh faces that should be embedded are missing in the final mesh",gr->tag(), (int)faces.size());
Msg::Error("Saving the missing faces in file %s",name);
FILE *f = fopen(name,"w");
fprintf(f,"View \" \" {\n");
for (std::map<MFace,GFace*,Less_Face>::iterator it = faces.begin() ; it != faces.end(); ++it){
MVertex *v1 = it->first.getVertex(0);
MVertex *v2 = it->first.getVertex(1);
MVertex *v3 = it->first.getVertex(2);
fprintf(f,"ST(%g,%g,%g,%g,%g,%g,%g,%g,%g){%d,%d,%d};\n",v1->x(),v1->y(),v1->z(),v2->x(),v2->y(),v2->z(),
v3->x(),v3->y(),v3->z(),it->second->tag(),it->second->tag(),it->second->tag());
}
fprintf(f,"};\n");
fclose(f);
}
}
};
Christophe Geuzaine
committed
template<class T>
static void GetQualityMeasure(std::vector<T*> &ele,
double &gamma, double &gammaMin, double &gammaMax,
double &minSICN, double &minSICNMin, double &minSICNMax,
double quality[3][100])
for(unsigned int i = 0; i < ele.size(); i++){
double g = ele[i]->gammaShapeMeasure();
gamma += g;
gammaMin = std::min(gammaMin, g);
double s = ele[i]->minSICNShapeMeasure();
minSICN += s;
minSICNMin = std::min(minSICNMin, s);
minSICNMax = std::max(minSICNMax, s);
rho += r;
rhoMin = std::min(rhoMin, r);
rhoMax = std::max(rhoMax, r);
for(int j = 0; j < 100; j++){
if(s > (2*j-100) / 100. && s <= (2*j-98) / 100.) quality[0][j]++;
if(g > j / 100. && g <= (j + 1) / 100.) quality[1][j]++;
if(r > j / 100. && r <= (j + 1) / 100.) quality[2][j]++;
}
}
void GetStatistics(double stat[50], double quality[3][100])
for(int i = 0; i < 50; i++) stat[i] = 0.;
stat[0] = m->getNumVertices();
stat[1] = m->getNumEdges();
stat[2] = m->getNumFaces();
stat[3] = m->getNumRegions();
std::map<int, std::vector<GEntity*> > physicals[4];
stat[45] = physicals[0].size() + physicals[1].size() +
physicals[2].size() + physicals[3].size();
for(GModel::eiter it = m->firstEdge(); it != m->lastEdge(); ++it)
for(GModel::fiter it = m->firstFace(); it != m->lastFace(); ++it){
// TODO: make this an option! if((*it)->getVisibility()){
stat[5] += (*it)->mesh_vertices.size();
stat[7] += (*it)->triangles.size();
stat[8] += (*it)->quadrangles.size();
for(GModel::riter it = m->firstRegion(); it != m->lastRegion(); ++it){
stat[6] += (*it)->mesh_vertices.size();
stat[9] += (*it)->tetrahedra.size();
stat[10] += (*it)->hexahedra.size();
stat[11] += (*it)->prisms.size();
stat[12] += (*it)->pyramids.size();
stat[14] = CTX::instance()->meshTimer[0];
stat[15] = CTX::instance()->meshTimer[1];
stat[16] = CTX::instance()->meshTimer[2];
if(quality){
for(int i = 0; i < 3; i++)
for(int j = 0; j < 100; j++)
double minSICN = 0., minSICNMin = 1., minSICNMax = -1.;
double gamma = 0., gammaMin = 1., gammaMax = 0.;
double rho = 0., rhoMin = 1., rhoMax = 0.;
double N = stat[9] + stat[10] + stat[11] + stat[12] + stat[13];
for(GModel::riter it = m->firstRegion(); it != m->lastRegion(); ++it){
GetQualityMeasure((*it)->tetrahedra, gamma, gammaMin, gammaMax,
minSICN, minSICNMin, minSICNMax, rho, rhoMin, rhoMax, quality);
Loading
Loading full blame...