diff --git a/Geo/MHexahedron.cpp b/Geo/MHexahedron.cpp index 9077b1e26a984b1ecf57c434693b0fb7a3ca2ff5..184dd63ac0793f24853a83678adb4ede8e90ae9d 100644 --- a/Geo/MHexahedron.cpp +++ b/Geo/MHexahedron.cpp @@ -14,6 +14,42 @@ #include "qualityMeasures.h" #endif +void MHexahedron::getEdgeRep(bool curved, int num, double *x, double *y, double *z, + SVector3 *n) +{ + // don't use MElement::_getEdgeRep: it's slow due to the creation of MFaces + MVertex *v0 = _v[edges_hexa(num, 0)]; + MVertex *v1 = _v[edges_hexa(num, 1)]; + x[0] = v0->x(); y[0] = v0->y(); z[0] = v0->z(); + x[1] = v1->x(); y[1] = v1->y(); z[1] = v1->z(); + if(CTX::instance()->mesh.lightLines > 1){ + static const int vv[12] = {2, 2, 3, 3, 2, 0, 1, 0, 6, 5, 4, 4}; + MVertex *v2 = _v[vv[num]]; + SVector3 t1(x[1] - x[0], y[1] - y[0], z[1] - z[0]); + SVector3 t2(v2->x() - x[0], v2->y() - y[0], v2->z() - z[0]); + SVector3 normal = crossprod(t1, t2); + normal.normalize(); + n[0] = n[1] = normal; + } + else{ + n[0] = n[1] = SVector3(0., 0., 1.); + } +} + +void MHexahedron::getFaceRep(bool curved, int num, double *x, double *y, double *z, + SVector3 *n) +{ + static const int f[12][3] = { + {0, 3, 2}, {0, 2, 1}, + {0, 1, 5}, {0, 5, 4}, + {0, 4, 7}, {0, 7, 3}, + {1, 2, 6}, {1, 6, 5}, + {2, 3, 7}, {2, 7, 6}, + {4, 5, 6}, {4, 6, 7} + }; + _getFaceRep(_v[f[num][0]], _v[f[num][1]], _v[f[num][2]], x, y, z, n); +} + int MHexahedron::getVolumeSign() { double mat[3][3]; diff --git a/Geo/MHexahedron.h b/Geo/MHexahedron.h index 30a25d63c633e61763f632786d09915bd6c618e8..d7609e79033bd170b249396f3a096b6d33302334 100644 --- a/Geo/MHexahedron.h +++ b/Geo/MHexahedron.h @@ -70,12 +70,7 @@ class MHexahedron : public MElement { return MEdge(_v[edges_hexa(num, 0)], _v[edges_hexa(num, 1)]); } virtual int getNumEdgesRep(bool curved){ return 12; } - virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n) - { - static const int f[12] = {0, 0, 1, 0, 1, 0, 3, 2, 1, 2, 3, 4}; - MEdge e(getEdge(num)); - _getEdgeRep(e.getVertex(0), e.getVertex(1), x, y, z, n, f[num]); - } + virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n); virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const { v.resize(2); @@ -93,19 +88,7 @@ class MHexahedron : public MElement { virtual double angleShapeMeasure(); virtual void getFaceInfo (const MFace & face, int &ithFace, int &sign, int &rot)const; virtual int getNumFacesRep(bool curved){ return 12; } - virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n) - { - static const int f[12][3] = { - {0, 3, 2}, {0, 2, 1}, - {0, 1, 5}, {0, 5, 4}, - {0, 4, 7}, {0, 7, 3}, - {1, 2, 6}, {1, 6, 5}, - {2, 3, 7}, {2, 7, 6}, - {4, 5, 6}, {4, 6, 7} - }; - _getFaceRep(getVertex(f[num][0]), getVertex(f[num][1]), getVertex(f[num][2]), - x, y, z, n); - } + virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n); virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const { v.resize(4); diff --git a/Geo/MQuadrangle.cpp b/Geo/MQuadrangle.cpp index 7a354659b7119fb5049c86378a3ec2ca822e83e7..89f345394e9c8ce544ff57451545a908c5a831e6 100644 --- a/Geo/MQuadrangle.cpp +++ b/Geo/MQuadrangle.cpp @@ -18,6 +18,37 @@ #define SQU(a) ((a)*(a)) +void MQuadrangle::getEdgeRep(bool curved, int num, double *x, double *y, double *z, + SVector3 *n) +{ + // don't use MElement::_getEdgeRep: it's slow due to the creation of MFace + MVertex *v0 = _v[edges_quad(num, 0)]; + MVertex *v1 = _v[edges_quad(num, 1)]; + x[0] = v0->x(); y[0] = v0->y(); z[0] = v0->z(); + x[1] = v1->x(); y[1] = v1->y(); z[1] = v1->z(); + if(CTX::instance()->mesh.lightLines > 1){ + static const int vv[4] = {2, 3, 0, 1}; + MVertex *v2 = _v[vv[num]]; + SVector3 t1(x[1] - x[0], y[1] - y[0], z[1] - z[0]); + SVector3 t2(v2->x() - x[0], v2->y() - y[0], v2->z() - z[0]); + SVector3 normal = crossprod(t1, t2); + normal.normalize(); + n[0] = n[1] = normal; + } + else{ + n[0] = n[1] = SVector3(0., 0., 1.); + } +} + +void MQuadrangle::getFaceRep(bool curved, int num, double *x, double *y, double *z, + SVector3 *n) +{ + static const int f[2][3] = { + {0, 1, 2}, {0, 2, 3} + }; + _getFaceRep(_v[f[num][0]], _v[f[num][1]], _v[f[num][2]], x, y, z, n); +} + int MQuadrangleN::getNumEdgesRep(bool curved) { return curved ? 4 * CTX::instance()->mesh.numSubEdges : 4; @@ -336,10 +367,10 @@ void MQuadrangle::reorient(int rot, bool swap) { std::memcpy(_v,tmp,4*sizeof(MVertex*)); } -void MQuadrangle8::reorient(int rot, bool swap) { - +void MQuadrangle8::reorient(int rot, bool swap) { + if (rot == 0 && !swap) return; - + MQuadrangle::reorient(rot,swap); MVertex* tmp[4]; if (swap) for (int i=0;i<4;i++) tmp[i] = _vs[(7-i+rot)%4]; @@ -347,10 +378,10 @@ void MQuadrangle8::reorient(int rot, bool swap) { std::memcpy(_vs,tmp,4*sizeof(MVertex*)); } -void MQuadrangle9::reorient(int rot, bool swap) { +void MQuadrangle9::reorient(int rot, bool swap) { if (rot == 0 && !swap) return; - + MQuadrangle::reorient(rot,swap); MVertex* tmp[4]; if (swap) for (int i=0;i<4;i++) tmp[i] = _vs[(7-i+rot)%4]; // edge swapped @@ -381,7 +412,7 @@ void MQuadrangleN::reorient(int rot, bool swap) { for (int i=0;i<nbEdgePts;i++) tmp.push_back(_vs[edgeIdx+i]); } } - + idx += 4*nbEdgePts; if (_vs.size() >= idx) { @@ -395,7 +426,7 @@ void MQuadrangleN::reorient(int rot, bool swap) { else for (int i=0;i<4;i++) tmp.push_back(_vs[idx + (4+i-rot)%4]); idx += 4; if (order > 4) Msg::Error("Reorientation of quad not supported above order 4"); - } + } } tmp.push_back(_vs[idx]); } diff --git a/Geo/MQuadrangle.h b/Geo/MQuadrangle.h index a4d744609f5a95fd816930a46c82f4291dbc4a15..092e44cbc525f33567947efca8d1df1aa562f165 100644 --- a/Geo/MQuadrangle.h +++ b/Geo/MQuadrangle.h @@ -83,11 +83,7 @@ class MQuadrangle : public MElement { Msg::Error("Could not get edge information for quadranglee %d", getNum()); } virtual int getNumEdgesRep(bool curved){ return 4; } - virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n) - { - MEdge e(getEdge(num)); - _getEdgeRep(e.getVertex(0), e.getVertex(1), x, y, z, n, 0); - } + virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n); virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const { v.resize(2); @@ -96,14 +92,7 @@ class MQuadrangle : public MElement { virtual int getNumFaces(){ return 1; } virtual MFace getFace(int num){ return MFace(_v[0], _v[1], _v[2], _v[3]); } virtual int getNumFacesRep(bool curved){ return 2; } - virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n) - { - static const int f[2][3] = { - {0, 1, 2}, {0, 2, 3} - }; - _getFaceRep(getVertex(f[num][0]), getVertex(f[num][1]), getVertex(f[num][2]), - x, y, z, n); - } + virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n); virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const { v.resize(4); diff --git a/Geo/MTetrahedron.cpp b/Geo/MTetrahedron.cpp index 4e76eb969336f83d975a84bc13467a807d12c5df..0bb4f764f861901c6f3562e56fff86b661d80257 100644 --- a/Geo/MTetrahedron.cpp +++ b/Geo/MTetrahedron.cpp @@ -42,18 +42,10 @@ void MTetrahedron::getEdgeRep(bool curved, int num, double *x, double *y, double void MTetrahedron::getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n) { - // don't use MElement::_getFaceRep: it's slow due to the creation of MFaces MVertex *v0 = _v[faces_tetra(num, 0)]; MVertex *v1 = _v[faces_tetra(num, 1)]; MVertex *v2 = _v[faces_tetra(num, 2)]; - x[0] = v0->x(); x[1] = v1->x(); x[2] = v2->x(); - y[0] = v0->y(); y[1] = v1->y(); y[2] = v2->y(); - z[0] = v0->z(); z[1] = v1->z(); z[2] = v2->z(); - SVector3 t1(x[1] - x[0], y[1] - y[0], z[1] - z[0]); - SVector3 t2(x[2] - x[0], y[2] - y[0], z[2] - z[0]); - SVector3 normal = crossprod(t1, t2); - normal.normalize(); - for(int i = 0; i < 3; i++) n[i] = normal; + _getFaceRep(v0, v1, v2, x, y, z, n); } SPoint3 MTetrahedron::circumcenter() diff --git a/Geo/MTriangle.cpp b/Geo/MTriangle.cpp index 836aebd428790931e9d156b4a063c87e8ed85668..203758f5766f29de7cf81fe8c4dbb66a1f7c31f0 100644 --- a/Geo/MTriangle.cpp +++ b/Geo/MTriangle.cpp @@ -42,15 +42,7 @@ void MTriangle::getEdgeRep(bool curved, int num, double *x, double *y, double *z void MTriangle::getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n) { - // don't use MElement::_getFaceRep: it's slow due to the creation of MFaces - x[0] = _v[0]->x(); x[1] = _v[1]->x(); x[2] = _v[2]->x(); - y[0] = _v[0]->y(); y[1] = _v[1]->y(); y[2] = _v[2]->y(); - z[0] = _v[0]->z(); z[1] = _v[1]->z(); z[2] = _v[2]->z(); - SVector3 t1(x[1] - x[0], y[1] - y[0], z[1] - z[0]); - SVector3 t2(x[2] - x[0], y[2] - y[0], z[2] - z[0]); - SVector3 normal = crossprod(t1, t2); - normal.normalize(); - for(int i = 0; i < 3; i++) n[i] = normal; + _getFaceRep(_v[0], _v[1], _v[2], x, y, z, n); } SPoint3 MTriangle::circumcenter()