diff --git a/Common/Views.cpp b/Common/Views.cpp index 0396689bb8fa7d13ac647c10ee42dac405c37648..74ea6f0dbe280aa1c1e0dd8f48b46b62ebda4cc3 100644 --- a/Common/Views.cpp +++ b/Common/Views.cpp @@ -1,4 +1,4 @@ -// $Id: Views.cpp,v 1.188 2006-04-18 09:57:32 remacle Exp $ +// $Id: Views.cpp,v 1.189 2006-08-08 01:10:04 geuzaine Exp $ // // Copyright (C) 1997-2006 C. Geuzaine, J.-F. Remacle // @@ -347,6 +347,21 @@ void Stat_Text(Post_View * v, List_T *D, List_T *C, int nb) if(nbtime > v->NbTimeStep) v->NbTimeStep = nbtime; } + + if(nb == 5){ + for(int i = 0; i < List_Nbr(D); i += nb){ + double x, y, z; + List_Read(D, i, &x); + List_Read(D, i+1, &y); + List_Read(D, i+2, &z); + if(x < v->BBox[0]) v->BBox[0] = x; + if(x > v->BBox[1]) v->BBox[1] = x; + if(y < v->BBox[2]) v->BBox[2] = y; + if(y > v->BBox[3]) v->BBox[3] = y; + if(z < v->BBox[4]) v->BBox[4] = z; + if(z > v->BBox[5]) v->BBox[5] = z; + } + } } void EndView(Post_View * v, int add_in_gui, char *file_name, char *name) diff --git a/Geo/GModelIO.cpp b/Geo/GModelIO.cpp index 092cc04114d26ecbe44d5b2de2244ddb9a7f447b..4419a2f4a6a549a7ddd48824cb98c835892d8d00 100644 --- a/Geo/GModelIO.cpp +++ b/Geo/GModelIO.cpp @@ -487,13 +487,13 @@ int GModel::writePOS(const std::string &name, double scalingFactor) "\"Gamma\", \"Eta\", \"Rho\"};\n", (1<<16)|(4<<8)); for(riter it = firstRegion(); it != lastRegion(); ++it) { for(unsigned int i = 0; i < (*it)->tetrahedra.size(); i++) - (*it)->tetrahedra[i]->writePOS(fp, scalingFactor); + (*it)->tetrahedra[i]->writePOS(fp, scalingFactor, (*it)->tag()); for(unsigned int i = 0; i < (*it)->hexahedra.size(); i++) - (*it)->hexahedra[i]->writePOS(fp, scalingFactor); + (*it)->hexahedra[i]->writePOS(fp, scalingFactor, (*it)->tag()); for(unsigned int i = 0; i < (*it)->prisms.size(); i++) - (*it)->prisms[i]->writePOS(fp, scalingFactor); + (*it)->prisms[i]->writePOS(fp, scalingFactor, (*it)->tag()); for(unsigned int i = 0; i < (*it)->pyramids.size(); i++) - (*it)->pyramids[i]->writePOS(fp, scalingFactor); + (*it)->pyramids[i]->writePOS(fp, scalingFactor, (*it)->tag()); } fprintf(fp, "};\n"); } @@ -504,9 +504,9 @@ int GModel::writePOS(const std::string &name, double scalingFactor) "\"Gamma\", \"Eta\", \"Rho\"};\n", (1<<16)|(4<<8)); for(fiter it = firstFace(); it != lastFace(); ++it) { for(unsigned int i = 0; i < (*it)->triangles.size(); i++) - (*it)->triangles[i]->writePOS(fp, scalingFactor); + (*it)->triangles[i]->writePOS(fp, scalingFactor, (*it)->tag()); for(unsigned int i = 0; i < (*it)->quadrangles.size(); i++) - (*it)->quadrangles[i]->writePOS(fp, scalingFactor); + (*it)->quadrangles[i]->writePOS(fp, scalingFactor, (*it)->tag()); } fprintf(fp, "};\n"); } @@ -517,7 +517,7 @@ int GModel::writePOS(const std::string &name, double scalingFactor) "\"Gamma\", \"Eta\", \"Rho\"};\n", (1<<16)|(4<<8)); for(eiter it = firstEdge(); it != lastEdge(); ++it) { for(unsigned int i = 0; i < (*it)->lines.size(); i++) - (*it)->lines[i]->writePOS(fp, scalingFactor); + (*it)->lines[i]->writePOS(fp, scalingFactor, (*it)->tag()); } fprintf(fp, "};\n"); } diff --git a/Geo/MElement.cpp b/Geo/MElement.cpp index 497f0e76d41951c25f5329616c928cd553c44701..dbc3231a529f8ab1c2069b26b73e1eed7891d48b 100644 --- a/Geo/MElement.cpp +++ b/Geo/MElement.cpp @@ -45,6 +45,33 @@ double MElement::rhoShapeMeasure() return 0.; } +double MTetrahedron::gammaShapeMeasure() +{ + double p0[3] = { _v[0]->x(), _v[0]->y(), _v[0]->z() }; + double p1[3] = { _v[1]->x(), _v[1]->y(), _v[1]->z() }; + double p2[3] = { _v[2]->x(), _v[2]->y(), _v[2]->z() }; + double p3[3] = { _v[3]->x(), _v[3]->y(), _v[3]->z() }; + double s1 = fabs(triangle_area(p0, p1, p2)); + double s2 = fabs(triangle_area(p0, p2, p3)); + double s3 = fabs(triangle_area(p0, p1, p3)); + double s4 = fabs(triangle_area(p1, p2, p3)); + double rhoin = 3. * fabs(getVolume()) / (s1 + s2 + s3 + s4); + return 12. * rhoin / (sqrt(6.) * maxEdge()); +} + +double MTetrahedron::etaShapeMeasure() +{ + double lij2 = 0.; + for(int i = 0; i <= 3; i++) { + for(int j = i + 1; j <= 3; j++) { + double lij = dist(_v[i], _v[j]); + lij2 += lij * lij; + } + } + double v = fabs(getVolume()); + return 12. * pow(0.9 * v * v, 1./3.) / lij2; +} + void MElement::cog(double &x, double &y, double &z) { x = y = z = 0.; @@ -98,7 +125,7 @@ void MElement::writeMSH(FILE *fp, double version, int num, int elementary, fprintf(fp, "\n"); } -void MElement::writePOS(FILE *fp, double scalingFactor) +void MElement::writePOS(FILE *fp, double scalingFactor, int elementary) { int n = getNumVertices(); double gamma = gammaShapeMeasure(); @@ -113,7 +140,7 @@ void MElement::writePOS(FILE *fp, double scalingFactor) } fprintf(fp, "){"); for(int i = 0; i < n; i++) - fprintf(fp, "%d,", getVertex(i)->onWhat()->tag()); + fprintf(fp, "%d,", elementary); for(int i = 0; i < n; i++) fprintf(fp, "%d,", getNum()); for(int i = 0; i < n; i++) diff --git a/Geo/MElement.h b/Geo/MElement.h index 4007a09b6d2f6a3bc37237214a3ab2cde38ae4a5..6adf9dafcab354959f2afe74b8f9a673c75cba15 100644 --- a/Geo/MElement.h +++ b/Geo/MElement.h @@ -94,13 +94,15 @@ class MElement // compute and change the orientation of 3D elements to get // positive volume + virtual double getVolume(){ return 0.; } virtual int getVolumeSign(){ return 1; } virtual void setVolumePositive(){} // IO routines virtual void writeMSH(FILE *fp, double version=1.0, int num=0, int elementary=1, int physical=1); - virtual void writePOS(FILE *fp, double scalingFactor=1.0); + virtual void writePOS(FILE *fp, double scalingFactor=1.0, + int elementary=1); virtual void writeSTL(FILE *fp, double scalingFactor=1.0); virtual void writeVRML(FILE *fp); virtual void writeUNV(FILE *fp, int type, int elementary); @@ -316,7 +318,7 @@ class MTetrahedron : public MElement { } int getTypeForMSH(){ return TET1; } char *getStringForPOS(){ return "SS"; } - virtual int getVolumeSign() + virtual double getVolume() { double mat[3][3]; mat[0][0] = _v[1]->x() - _v[0]->x(); @@ -328,8 +330,9 @@ class MTetrahedron : public MElement { mat[2][0] = _v[1]->z() - _v[0]->z(); mat[2][1] = _v[2]->z() - _v[0]->z(); mat[2][2] = _v[3]->z() - _v[0]->z(); - return sign(det3x3(mat)); + return det3x3(mat) / 6.; } + virtual int getVolumeSign(){ return sign(getVolume()); } void setVolumePositive() { if(getVolumeSign() < 0){ @@ -337,6 +340,8 @@ class MTetrahedron : public MElement { tmp = _v[0]; _v[0] = _v[1]; _v[1] = tmp; } } + virtual double gammaShapeMeasure(); + virtual double etaShapeMeasure(); }; // TODO: for MTetrahedron2 diff --git a/Numeric/Numeric.cpp b/Numeric/Numeric.cpp index c119ce751fa97996bc508d0a78d33880061c54a6..74d9b8f5092c03b17f51008a37e4702f9a77281e 100644 --- a/Numeric/Numeric.cpp +++ b/Numeric/Numeric.cpp @@ -1,4 +1,4 @@ -// $Id: Numeric.cpp,v 1.25 2006-01-17 17:09:05 geuzaine Exp $ +// $Id: Numeric.cpp,v 1.26 2006-08-08 01:10:05 geuzaine Exp $ // // Copyright (C) 1997-2006 C. Geuzaine, J.-F. Remacle // @@ -283,6 +283,22 @@ double angle_02pi(double A3) return A3; } +double triangle_area(double p0[3], double p1[3], double p2[3]) +{ + double a[3], b[3], c[3]; + + a[0] = p2[0] - p1[0]; + a[1] = p2[1] - p1[1]; + a[2] = p2[2] - p1[2]; + + b[0] = p0[0] - p1[0]; + b[1] = p0[1] - p1[1]; + b[2] = p0[2] - p1[2]; + + prodve(a, b, c); + return (0.5 * sqrt(c[0] * c[0] + c[1] * c[1] + c[2] * c[2])); +} + double InterpolateIso(double *X, double *Y, double *Z, double *Val, double V, int I1, int I2, double *XI, double *YI, double *ZI) diff --git a/Numeric/Numeric.h b/Numeric/Numeric.h index da87bf8344832ae411d859999e5129d0a8b9f50c..23a56da1368ddf0c9ae37503d6d30ce4f285949c 100644 --- a/Numeric/Numeric.h +++ b/Numeric/Numeric.h @@ -76,6 +76,7 @@ double trace3x3(double mat[3][3]); double trace2 (double mat[3][3]); double inv3x3(double mat[3][3], double inv[3][3]); double angle_02pi(double A3); +double triangle_area(double p0[3], double p1[3], double p2[3]); void eigenvalue(double mat[3][3], double re[3]); void FindCubicRoots(const double coeff[4], double re[3], double im[3]); void eigsort(double d[3]); diff --git a/Parser/OpenFile.cpp b/Parser/OpenFile.cpp index 52f0a33088f16e2684ccd4660c07492610a7e373..6b0784a9b1d3f2d234b422cd3d6cd0f9a0f9c16b 100644 --- a/Parser/OpenFile.cpp +++ b/Parser/OpenFile.cpp @@ -1,4 +1,4 @@ -// $Id: OpenFile.cpp,v 1.103 2006-08-07 21:04:05 geuzaine Exp $ +// $Id: OpenFile.cpp,v 1.104 2006-08-08 01:10:05 geuzaine Exp $ // // Copyright (C) 1997-2006 C. Geuzaine, J.-F. Remacle // @@ -149,7 +149,7 @@ void SetBoundingBox(void) if(bb.empty() && GMODEL) bb = GMODEL->recomputeBounds(); - if(List_Nbr(CTX.post.list)) { + if(bb.empty() && List_Nbr(CTX.post.list)) { for(int i = 0; i < List_Nbr(CTX.post.list); i++){ Post_View *v = *(Post_View **)List_Pointer(CTX.post.list, i); bb += SPoint3(v->BBox[0], v->BBox[2], v->BBox[4]);