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

* reduce (slightly) memory used by 3d algo by making sure we dealloc unused 

data asap
* finer estimate of mem usage
parent a4851849
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Mesh/meshGRegion.cpp
+ 43
− 44
View file @ 355ce8b8
// $Id: meshGRegion.cpp,v 1.43 2008-02-22 21:09:01 geuzaine Exp $ // $Id: meshGRegion.cpp,v 1.44 2008-03-06 14:19:01 geuzaine Exp $
// //
// Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
// //
...@@ -96,27 +96,25 @@ void buildTetgenStructure(GRegion *gr, tetgenio &in, std::vector<MVertex*> &numb ...@@ -96,27 +96,25 @@ void buildTetgenStructure(GRegion *gr, tetgenio &in, std::vector<MVertex*> &numb
for(unsigned int i = 0; i < gf->triangles.size(); i++){ for(unsigned int i = 0; i < gf->triangles.size(); i++){
MTriangle *t = gf->triangles[i]; MTriangle *t = gf->triangles[i];
tetgenio::facet *f = &in.facetlist[I]; tetgenio::facet *f = &in.facetlist[I];
tetgenio::init(f); tetgenio::init(f);
f->numberofholes = 0; f->numberofholes = 0;
f->numberofpolygons = 1; f->numberofpolygons = 1;
f->polygonlist = new tetgenio::polygon[f->numberofpolygons]; f->polygonlist = new tetgenio::polygon[f->numberofpolygons];
tetgenio::polygon *p = &f->polygonlist[0]; tetgenio::polygon *p = &f->polygonlist[0];
tetgenio::init(p); tetgenio::init(p);
p->numberofvertices = 3; p->numberofvertices = 3;
p->vertexlist = new int[p->numberofvertices]; p->vertexlist = new int[p->numberofvertices];
p->vertexlist[0] = t->getVertex(0)->getNum(); p->vertexlist[0] = t->getVertex(0)->getNum();
p->vertexlist[1] = t->getVertex(1)->getNum(); p->vertexlist[1] = t->getVertex(1)->getNum();
p->vertexlist[2] = t->getVertex(2)->getNum(); p->vertexlist[2] = t->getVertex(2)->getNum();
in.facetmarkerlist[I] = gf->tag(); in.facetmarkerlist[I] = gf->tag();
++I; ++I;
} }
++it; ++it;
} }
} }
void TransferTetgenMesh(GRegion *gr, void TransferTetgenMesh(GRegion *gr, tetgenio &in, tetgenio &out,
tetgenio &in,
tetgenio &out,
std::vector<MVertex*> &numberedV) std::vector<MVertex*> &numberedV)
{ {
for(int i = numberedV.size(); i < out.numberofpoints; i++){ for(int i = numberedV.size(); i < out.numberofpoints; i++){
...@@ -141,13 +139,13 @@ void TransferTetgenMesh(GRegion *gr, ...@@ -141,13 +139,13 @@ void TransferTetgenMesh(GRegion *gr,
gf->triangles.clear(); gf->triangles.clear();
gf->deleteVertexArrays(); gf->deleteVertexArrays();
++it; ++it;
} }
// TODO: re-create 1D mesh // TODO: re-create 1D mesh
for(int i = 0; i < out.numberofedges; i++){ for(int i = 0; i < out.numberofedges; i++){
MVertex *v[2]; MVertex *v[2];
v[0] = numberedV[out.edgelist[i * 2 + 0] -1]; v[0] = numberedV[out.edgelist[i * 2 + 0] - 1];
v[1] = numberedV[out.edgelist[i * 2 + 1] -1]; v[1] = numberedV[out.edgelist[i * 2 + 1] - 1];
} }
// re-create the triangular meshes FIXME: this can lead to hanging // re-create the triangular meshes FIXME: this can lead to hanging
...@@ -212,38 +210,40 @@ void MeshDelaunayVolume(std::vector<GRegion*> &regions) ...@@ -212,38 +210,40 @@ void MeshDelaunayVolume(std::vector<GRegion*> &regions)
allFaces.push_back(*it); allFaces.push_back(*it);
gr->set(allFaces); gr->set(allFaces);
// mesh with tetgen, possibly changing the mesh on boundaries // mesh with tetgen, possibly changing the mesh on boundaries (leave
tetgenio in, out; // this in block, so in/out are destroyed before we refine the mesh)
std::vector<MVertex*> numberedV; {
char opts[128]; tetgenio in, out;
buildTetgenStructure(gr, in, numberedV); std::vector<MVertex*> numberedV;
sprintf(opts, "pe%c", (CTX.verbosity < 3) ? 'Q': (CTX.verbosity > 6)? 'V': '\0'); char opts[128];
try{ buildTetgenStructure(gr, in, numberedV);
tetrahedralize(opts, &in, &out); sprintf(opts, "pe%c", (CTX.verbosity < 3) ? 'Q': (CTX.verbosity > 6)? 'V': '\0');
}
catch (int error){
Msg (WARNING, "Self intersecting Surface Mesh, computing intersections "
"(this could take a while)");
sprintf(opts, "dV");
try{ try{
tetrahedralize(opts, &in, &out); tetrahedralize(opts, &in, &out);
Msg(INFO,"%d faces self-intersect",out.numberoftrifaces); }
for (int i=0;i<out.numberoftrifaces;i++) catch (int error){
{ Msg (WARNING, "Self intersecting Surface Mesh, computing intersections "
"(this could take a while)");
sprintf(opts, "dV");
try{
tetrahedralize(opts, &in, &out);
Msg(INFO,"%d faces self-intersect",out.numberoftrifaces);
for (int i = 0; i < out.numberoftrifaces; i++){
Msg(INFO,"face (%d %d %d) on model face %d", Msg(INFO,"face (%d %d %d) on model face %d",
numberedV[out.trifacelist[i * 3 + 0] - 1]->getNum(), numberedV[out.trifacelist[i * 3 + 0] - 1]->getNum(),
numberedV[out.trifacelist[i * 3 + 1] - 1]->getNum(), numberedV[out.trifacelist[i * 3 + 1] - 1]->getNum(),
numberedV[out.trifacelist[i * 3 + 2] - 1]->getNum(), numberedV[out.trifacelist[i * 3 + 2] - 1]->getNum(),
out.trifacemarkerlist[i]); out.trifacemarkerlist[i]);
} }
}
catch (int error2){
Msg(GERROR, "Surface Mesh is wrong, cannot do the 3D mesh");
}
gr->set(faces);
return;
} }
catch (int error2){ TransferTetgenMesh(gr, in, out, numberedV);
Msg (GERROR, "Surface Mesh is wrong, cannot do the 3D mesh");
}
gr->set(faces);
return;
} }
TransferTetgenMesh(gr, in, out, numberedV);
// sort triangles in all model faces in order to be able to search in vectors // sort triangles in all model faces in order to be able to search in vectors
std::list<GFace*>::iterator itf = allFaces.begin(); std::list<GFace*>::iterator itf = allFaces.begin();
...@@ -257,7 +257,7 @@ void MeshDelaunayVolume(std::vector<GRegion*> &regions) ...@@ -257,7 +257,7 @@ void MeshDelaunayVolume(std::vector<GRegion*> &regions)
gr->set(faces); gr->set(faces);
// now do insertion of points // now do insertion of points
insertVerticesInRegion(gr); insertVerticesInRegion(gr);
Msg(INFO, "Gmsh 3D Delaunay has generated %d tets", gr->tetrahedra.size()); Msg(INFO, "Gmsh 3D Delaunay has generated %d tets", gr->tetrahedra.size());
#endif #endif
} }
...@@ -306,13 +306,13 @@ Ng_Mesh *buildNetgenStructure(GRegion *gr, bool importVolumeMesh, ...@@ -306,13 +306,13 @@ Ng_Mesh *buildNetgenStructure(GRegion *gr, bool importVolumeMesh,
std::list<GFace*> faces = gr->faces(); std::list<GFace*> faces = gr->faces();
std::list<GFace*>::iterator it = faces.begin(); std::list<GFace*>::iterator it = faces.begin();
while(it != faces.end()){ while(it != faces.end()){
GFace *gf = (*it); GFace *gf = (*it);
for(unsigned int i = 0; i< gf->triangles.size(); i++){ for(unsigned int i = 0; i< gf->triangles.size(); i++){
MTriangle *t = gf->triangles[i]; MTriangle *t = gf->triangles[i];
int tmp[3]; int tmp[3];
tmp[0] = t->getVertex(0)->getNum(); tmp[0] = t->getVertex(0)->getNum();
tmp[1] = t->getVertex(1)->getNum(); tmp[1] = t->getVertex(1)->getNum();
tmp[2] = t->getVertex(2)->getNum(); tmp[2] = t->getVertex(2)->getNum();
Ng_AddSurfaceElement(ngmesh, NG_TRIG, tmp); Ng_AddSurfaceElement(ngmesh, NG_TRIG, tmp);
} }
++it; ++it;
...@@ -326,8 +326,8 @@ Ng_Mesh *buildNetgenStructure(GRegion *gr, bool importVolumeMesh, ...@@ -326,8 +326,8 @@ Ng_Mesh *buildNetgenStructure(GRegion *gr, bool importVolumeMesh,
int tmp[4]; int tmp[4];
tmp[0] = t->getVertex(0)->getNum(); tmp[0] = t->getVertex(0)->getNum();
tmp[1] = t->getVertex(1)->getNum(); tmp[1] = t->getVertex(1)->getNum();
tmp[2] = t->getVertex(2)->getNum(); tmp[2] = t->getVertex(2)->getNum();
tmp[3] = t->getVertex(3)->getNum(); tmp[3] = t->getVertex(3)->getNum();
Ng_AddVolumeElement(ngmesh, NG_TET, tmp); Ng_AddVolumeElement(ngmesh, NG_TET, tmp);
} }
} }
...@@ -416,10 +416,10 @@ int intersect_line_triangle(double X[3], double Y[3], double Z[3] , ...@@ -416,10 +416,10 @@ int intersect_line_triangle(double X[3], double Y[3], double Z[3] ,
b[1] = P[1] - Y[0]; b[1] = P[1] - Y[0];
b[2] = P[2] - Z[0]; b[2] = P[2] - Z[0];
if(!sys3x3_with_tol(mat, b, res, &det)) if(!sys3x3_with_tol(mat, b, res, &det))
return 0; return 0;
// Msg(INFO,"going there %g %g %g",res[0],res[1],res[2]); // Msg(INFO, "going there %g %g %g", res[0], res[1], res[2]);
if(res[0] >= eps_prec && res[0] <= 1.0 - eps_prec && if(res[0] >= eps_prec && res[0] <= 1.0 - eps_prec &&
res[1] >= eps_prec && res[1] <= 1.0 - eps_prec && res[1] >= eps_prec && res[1] <= 1.0 - eps_prec &&
...@@ -454,7 +454,7 @@ void meshNormalsPointOutOfTheRegion(GRegion *gr) ...@@ -454,7 +454,7 @@ void meshNormalsPointOutOfTheRegion(GRegion *gr)
setRand(rrr); setRand(rrr);
while(it != faces.end()){ while(it != faces.end()){
GFace *gf = (*it); GFace *gf = (*it);
int nb_intersect = 0; int nb_intersect = 0;
for(unsigned int i = 0; i < gf->triangles.size(); i++){ for(unsigned int i = 0; i < gf->triangles.size(); i++){
MTriangle *t = gf->triangles[i]; MTriangle *t = gf->triangles[i];
...@@ -464,7 +464,7 @@ void meshNormalsPointOutOfTheRegion(GRegion *gr) ...@@ -464,7 +464,7 @@ void meshNormalsPointOutOfTheRegion(GRegion *gr)
double P[3] = {(X[0]+X[1]+X[2])/3., (Y[0]+Y[1]+Y[2])/3., (Z[0]+Z[1]+Z[2])/3.}; double P[3] = {(X[0]+X[1]+X[2])/3., (Y[0]+Y[1]+Y[2])/3., (Z[0]+Z[1]+Z[2])/3.};
double v1[3] = {X[0]-X[1], Y[0]-Y[1], Z[0]-Z[1]}; double v1[3] = {X[0]-X[1], Y[0]-Y[1], Z[0]-Z[1]};
double v2[3] = {X[2]-X[1], Y[2]-Y[1], Z[2]-Z[1]}; double v2[3] = {X[2]-X[1], Y[2]-Y[1], Z[2]-Z[1]};
double N[3] ; double N[3];
prodve(v1, v2, N); prodve(v1, v2, N);
norme(v1); norme(v1);
norme(v2); norme(v2);
...@@ -611,7 +611,6 @@ void optimizeMeshGRegionGmsh::operator() (GRegion *gr) ...@@ -611,7 +611,6 @@ void optimizeMeshGRegionGmsh::operator() (GRegion *gr)
// import mesh into netgen, including volume tets // import mesh into netgen, including volume tets
gmshOptimizeMesh (gr, QMTET_2); gmshOptimizeMesh (gr, QMTET_2);
} }
bool buildFaceSearchStructure(GModel *model, fs_cont &search) bool buildFaceSearchStructure(GModel *model, fs_cont &search)
......
Mesh/meshGRegionDelaunayInsertion.cpp
+ 31
− 35
View file @ 355ce8b8
// $Id: meshGRegionDelaunayInsertion.cpp,v 1.38 2008-03-04 08:51:14 geuzaine Exp $ // $Id: meshGRegionDelaunayInsertion.cpp,v 1.39 2008-03-06 14:19:01 geuzaine Exp $
// //
// Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
// //
...@@ -76,7 +76,7 @@ struct faceXtet{ ...@@ -76,7 +76,7 @@ struct faceXtet{
}; };
template <class ITER> template <class ITER>
void connectTets ( ITER beg, ITER end) void connectTets(ITER beg, ITER end)
{ {
std::set<faceXtet> conn; std::set<faceXtet> conn;
while (beg != end){ while (beg != end){
...@@ -131,11 +131,9 @@ void recurFindCavity(std::list<faceXtet> & shell, ...@@ -131,11 +131,9 @@ void recurFindCavity(std::list<faceXtet> & shell,
} }
} }
bool insertVertex(MVertex *v, bool insertVertex(MVertex *v, MTet4 *t, MTet4Factory &myFactory,
MTet4 *t, std::set<MTet4*, compareTet4Ptr> &allTets,
MTet4Factory &myFactory, std::vector<double> &vSizes)
std::set<MTet4*,compareTet4Ptr> &allTets,
std::vector<double> & vSizes)
{ {
std::list<faceXtet> shell; std::list<faceXtet> shell;
std::list<MTet4*> cavity; std::list<MTet4*> cavity;
...@@ -175,7 +173,7 @@ bool insertVertex(MVertex *v, ...@@ -175,7 +173,7 @@ bool insertVertex(MVertex *v,
++ittet; ++ittet;
} }
// fprintf(ff2,"};\n"); // fprintf(ff2,"};\n");
// fclose (ff2); // fclose(ff2);
// Msg(INFO,"cavity of size %d volume %g",cavity.size(),oldVolume); // Msg(INFO,"cavity of size %d volume %g",cavity.size(),oldVolume);
// create new tetrahedron using faces that are // create new tetrahedron using faces that are
// on the border of the cavity // on the border of the cavity
...@@ -186,7 +184,7 @@ bool insertVertex(MVertex *v, ...@@ -186,7 +184,7 @@ bool insertVertex(MVertex *v,
// char name[245]; // char name[245];
// sprintf(name,"test%d.pos",III); // sprintf(name,"test%d.pos",III);
// FILE *ff = fopen (name,"w"); // FILE *ff = fopen(name,"w");
// fprintf(ff,"View\"test\"{\n"); // fprintf(ff,"View\"test\"{\n");
MTet4** newTets = new MTet4*[shell.size()];; MTet4** newTets = new MTet4*[shell.size()];;
...@@ -195,7 +193,7 @@ bool insertVertex(MVertex *v, ...@@ -195,7 +193,7 @@ bool insertVertex(MVertex *v,
std::list<faceXtet>::iterator it = shell.begin(); std::list<faceXtet>::iterator it = shell.begin();
while (it != shell.end()){ while (it != shell.end()){
MTetrahedron *tr = new MTetrahedron (it->v[0], it->v[1], it->v[2], v); MTetrahedron *tr = new MTetrahedron(it->v[0], it->v[1], it->v[2], v);
// Msg(INFO,"shell %d %d %d",it->v[0]->getNum(),it->v[1]->getNum(),it->v[2]->getNum()); // Msg(INFO,"shell %d %d %d",it->v[0]->getNum(),it->v[1]->getNum(),it->v[2]->getNum());
// fprintf(ff,"ST(%g,%g,%g,%g,%g,%g,%g,%g,%g) {0,0,0};\n", // fprintf(ff,"ST(%g,%g,%g,%g,%g,%g,%g,%g,%g) {0,0,0};\n",
// it->v[0]->x(), // it->v[0]->x(),
...@@ -229,8 +227,8 @@ bool insertVertex(MVertex *v, ...@@ -229,8 +227,8 @@ bool insertVertex(MVertex *v,
// Msg(INFO,"new cavity of vol %g (%d boundaries)",newVolume,shell.size()); // Msg(INFO,"new cavity of vol %g (%d boundaries)",newVolume,shell.size());
// OK, the cavity is star shaped // OK, the cavity is star shaped
if (fabs(oldVolume - newVolume) < 1.e-10 * oldVolume){ if (fabs(oldVolume - newVolume) < 1.e-10 * oldVolume){
connectTets ( new_cavity.begin(),new_cavity.end()); connectTets(new_cavity.begin(), new_cavity.end());
allTets.insert(newTets,newTets+shell.size()); allTets.insert(newTets, newTets + shell.size());
// ittet = cavity.begin(); // ittet = cavity.begin();
// ittete = cavity.end(); // ittete = cavity.end();
...@@ -254,7 +252,7 @@ bool insertVertex(MVertex *v, ...@@ -254,7 +252,7 @@ bool insertVertex(MVertex *v,
} }
} }
static void setLcs(MTetrahedron *t, std::map<MVertex*,double> &vSizes) static void setLcs(MTetrahedron *t, std::map<MVertex*, double> &vSizes)
{ {
for (int i = 0; i < 4; i++){ for (int i = 0; i < 4; i++){
for (int j = i + 1; j < 4; j++){ for (int j = i + 1; j < 4; j++){
...@@ -319,12 +317,9 @@ GRegion *getRegionFromBoundingFaces(GModel *model, ...@@ -319,12 +317,9 @@ GRegion *getRegionFromBoundingFaces(GModel *model,
return 0; return 0;
} }
void recur_classify(MTet4 *t , void recur_classify(MTet4 *t, std::list<MTet4*> &theRegion,
std::list<MTet4*> &theRegion, std::set<GFace*> &faces_bound, GRegion *bidon,
std::set<GFace *> &faces_bound, GModel *model, const fs_cont &search)
GRegion *bidon ,
GModel *model,
const fs_cont &search)
{ {
if (!t) Msg (GERROR,"a tet is not connected by a boundary face"); if (!t) Msg (GERROR,"a tet is not connected by a boundary face");
if (t->onWhat()) return; // should never return here... if (t->onWhat()) return; // should never return here...
...@@ -368,7 +363,7 @@ void adaptMeshGRegion::operator () (GRegion *gr) ...@@ -368,7 +363,7 @@ void adaptMeshGRegion::operator () (GRegion *gr)
typedef std::list<MTet4 *> CONTAINER ; typedef std::list<MTet4 *> CONTAINER ;
CONTAINER allTets; CONTAINER allTets;
for(unsigned int i=0;i<gr->tetrahedra.size();i++){ for(unsigned int i = 0; i < gr->tetrahedra.size(); i++){
allTets.push_back(new MTet4(gr->tetrahedra[i], qm)); allTets.push_back(new MTet4(gr->tetrahedra[i], qm));
} }
gr->tetrahedra.clear(); gr->tetrahedra.clear();
...@@ -377,8 +372,8 @@ void adaptMeshGRegion::operator () (GRegion *gr) ...@@ -377,8 +372,8 @@ void adaptMeshGRegion::operator () (GRegion *gr)
double t1 = Cpu(); double t1 = Cpu();
std::vector<MTet4*> illegals; std::vector<MTet4*> illegals;
const int nbRanges=10; const int nbRanges = 10;
int quality_ranges [nbRanges]; int quality_ranges[nbRanges];
{ {
double totalVolumeb = 0.0; double totalVolumeb = 0.0;
double worst = 1.0; double worst = 1.0;
...@@ -545,7 +540,7 @@ void adaptMeshGRegion::operator () (GRegion *gr) ...@@ -545,7 +540,7 @@ void adaptMeshGRegion::operator () (GRegion *gr)
} }
//template <class CONTAINER, class DATA> //template <class CONTAINER, class DATA>
void gmshOptimizeMesh (GRegion *gr, const gmshQualityMeasure4Tet &qm) void gmshOptimizeMesh(GRegion *gr, const gmshQualityMeasure4Tet &qm)
{ {
typedef std::list<MTet4 *> CONTAINER ; typedef std::list<MTet4 *> CONTAINER ;
CONTAINER allTets; CONTAINER allTets;
...@@ -726,19 +721,20 @@ void insertVerticesInRegion (GRegion *gr) ...@@ -726,19 +721,20 @@ void insertVerticesInRegion (GRegion *gr)
//printf("sizeof MTet4 = %d sizeof MTetrahedron %d sizeof(MVertex) %d\n", //printf("sizeof MTet4 = %d sizeof MTetrahedron %d sizeof(MVertex) %d\n",
// sizeof(MTet4), sizeof(MTetrahedron), sizeof(MVertex)); // sizeof(MTet4), sizeof(MTetrahedron), sizeof(MVertex));
std::map<MVertex*,double> vSizesMap;
std::vector<double> vSizes; std::vector<double> vSizes;
MTet4Factory myFactory(1600000); MTet4Factory myFactory(1600000);
std::set<MTet4*,compareTet4Ptr> &allTets = myFactory.getAllTets(); std::set<MTet4*, compareTet4Ptr> &allTets = myFactory.getAllTets();
for(unsigned int i = 0; i < gr->tetrahedra.size(); i++)
setLcs(gr->tetrahedra[i], vSizesMap);
int NUM = 0; int NUM = 0;
for(std::map<MVertex*, double>::iterator it = vSizesMap.begin();
it != vSizesMap.end(); ++it){ { // leave this in a block so the map gets deallocated directly
it->first->setNum(NUM++); std::map<MVertex*, double> vSizesMap;
vSizes.push_back(it->second); for(unsigned int i = 0; i < gr->tetrahedra.size(); i++)
setLcs(gr->tetrahedra[i], vSizesMap);
for(std::map<MVertex*, double>::iterator it = vSizesMap.begin();
it != vSizesMap.end(); ++it){
it->first->setNum(NUM++);
vSizes.push_back(it->second);
}
} }
for(unsigned int i = 0; i < gr->tetrahedra.size(); i++) for(unsigned int i = 0; i < gr->tetrahedra.size(); i++)
...@@ -751,7 +747,7 @@ void insertVerticesInRegion (GRegion *gr) ...@@ -751,7 +747,7 @@ void insertVerticesInRegion (GRegion *gr)
// Msg (INFO,"reclassifying %d tets", allTets.size()); // Msg (INFO,"reclassifying %d tets", allTets.size());
fs_cont search; fs_cont search;
buildFaceSearchStructure(gr->model(), search ); buildFaceSearchStructure(gr->model(), search);
for(MTet4Factory::iterator it = allTets.begin(); it != allTets.end(); ++it){ for(MTet4Factory::iterator it = allTets.begin(); it != allTets.end(); ++it){
if(!(*it)->onWhat()){ if(!(*it)->onWhat()){
...@@ -812,7 +808,7 @@ void insertVerticesInRegion (GRegion *gr) ...@@ -812,7 +808,7 @@ void insertVerticesInRegion (GRegion *gr)
double uvw[3]; double uvw[3];
bool inside = worst->tet()->invertmapping(center, uvw); bool inside = worst->tet()->invertmapping(center, uvw);
if(inside){ if(inside){
MVertex *v = new MVertex(center[0],center[1],center[2], worst->onWhat()); MVertex *v = new MVertex(center[0], center[1], center[2], worst->onWhat());
v->setNum(NUM++); v->setNum(NUM++);
double lc1 = double lc1 =
(1 - uvw[0] - uvw[1] - uvw[2]) * vSizes[worst->tet()->getVertex(0)->getNum()] + (1 - uvw[0] - uvw[1] - uvw[2]) * vSizes[worst->tet()->getVertex(0)->getNum()] +
......
Mesh/meshGRegionDelaunayInsertion.h
+ 26
− 14
View file @ 355ce8b8
...@@ -34,16 +34,28 @@ class GModel; ...@@ -34,16 +34,28 @@ class GModel;
class MTet4Factory; class MTet4Factory;
// sizeof(MTet4) = 36 Bytes and sizeof(MTetrahedron) = 28 Bytes, so // Memory usage for 1 million tets:
// normally it should take 36+28 = 64 MByte per million tets. We also //
// store a rb tree containing all pointers sorted with respect to tet // * sizeof(MTet4) = 36 Bytes and sizeof(MTetrahedron) = 28 Bytes
// radius. Each bucket of the tree contain 4 pointers, i.e. 16 Bytes // -> 64 MB
// plus the data -> Extra cost of 20 Bytes/Tet, i.e., 84 MB per // * rb tree containing all pointers sorted with respect to tet
// million tets. A MVertex has a cost of 44 Bytes and there are about // radius: each bucket of the tree contains 4 pointers (16 Bytes)
// 200000 of them per million tet, i.e., a new cost of 9MB per million // plus the data -> 20 MB
// tets. // * sizeof(MVertex) = 44 Bytes and there are about 200000 verts per
// million tet -> 9MB
// Grand total should be 92 MB per million tet (I observe 160M MB!) // * vector of char lengths per vertex -> 1.6Mb
// * vectors in GEntities to store the element and vertex pointers
// -> 5Mb
//
// Grand total should thus be about 100 MB.
//
// The observed mem usage with "demos/cube.geo -clscale 0.61" is
// 157MB. Where do the extra 57 MB come from?
//
// * surface mesh + all other overheads (model, etc.) is 19Mb
// * tetgen initial mesh is about 20Mb, but it is deleted before mesh
// refinement.
// * ?
class MTet4 class MTet4
{ {
...@@ -92,8 +104,8 @@ class MTet4 ...@@ -92,8 +104,8 @@ class MTet4
} }
inline GRegion *onWhat () const { return gr; } inline GRegion *onWhat () const { return gr; }
inline void setOnWhat (GRegion *g) { gr = g; } inline void setOnWhat (GRegion *g) { gr = g; }
bool isDeleted () const { return deleted; } inline bool isDeleted () const { return deleted; }
void forceRadius (double r){ circum_radius = r; } inline void forceRadius (double r){ circum_radius = r; }
inline double getRadius () const { return circum_radius; } inline double getRadius () const { return circum_radius; }
inline double getQuality () const { return circum_radius; } inline double getQuality () const { return circum_radius; }
inline void setQuality (const double &q){ circum_radius = q; } inline void setQuality (const double &q){ circum_radius = q; }
...@@ -111,7 +123,7 @@ class MTet4 ...@@ -111,7 +123,7 @@ class MTet4
{ {
return inCircumSphere(v->x(), v->y(), v->z()); return inCircumSphere(v->x(), v->y(), v->z());
} }
double getVolume() const { return base->getVolume(); } inline double getVolume() const { return base->getVolume(); }
inline void setDeleted(bool d) inline void setDeleted(bool d)
{ {
deleted = d; deleted = d;
...@@ -142,7 +154,7 @@ class compareTet4Ptr ...@@ -142,7 +154,7 @@ class compareTet4Ptr
{ {
if (a->getRadius() > b->getRadius()) return true; if (a->getRadius() > b->getRadius()) return true;
if (a->getRadius() < b->getRadius()) return false; if (a->getRadius() < b->getRadius()) return false;
return a<b; return a < b;
} }
}; };
......
benchmarks/3d/periodic.geo
+ 1
− 1
View file @ 355ce8b8
...@@ -46,7 +46,7 @@ Physical Line(18) = {2,3}; ...@@ -46,7 +46,7 @@ Physical Line(18) = {2,3};
Physical Surface(19) = {15}; Physical Surface(19) = {15};
If(use_prisms) If(use_prisms)
Extrude Surface {15, {0.,0.,2.*R}}{Layers{nb_layers,83,1}; Recombine; }; Extrude Surface {15, {0.,0.,2.*R}}{ Layers{nb_layers}; Recombine; };
EndIf EndIf
If(!use_prisms) If(!use_prisms)
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment