Skip to content
Snippets Groups Projects
Commit 2344885c authored by Amaury Johnen's avatar Amaury Johnen
Browse files

Simplification

parent 2649f2ec
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Geo/MHexahedron.cpp
+ 11
− 49
View file @ 2344885c
...@@ -47,7 +47,7 @@ double MHexahedron::angleShapeMeasure() ...@@ -47,7 +47,7 @@ double MHexahedron::angleShapeMeasure()
std::vector<MVertex*> vv; std::vector<MVertex*> vv;
vv.push_back(getFace(i).getVertex(0)); vv.push_back(getFace(i).getVertex(0));
vv.push_back(getFace(i).getVertex(1)); vv.push_back(getFace(i).getVertex(1));
vv.push_back(getFace(i).getVertex(2)); vv.push_back(getFace(i).getVertex(2));
vv.push_back(getFace(i).getVertex(3)); vv.push_back(getFace(i).getVertex(3));
// MVertex *v0 = new MVertex(0, 0, 0); vv.push_back(v0); // MVertex *v0 = new MVertex(0, 0, 0); vv.push_back(v0);
// MVertex *v1 = new MVertex(1., 0, 0);vv.push_back(v1); // MVertex *v1 = new MVertex(1., 0, 0);vv.push_back(v1);
...@@ -63,7 +63,7 @@ double MHexahedron::angleShapeMeasure() ...@@ -63,7 +63,7 @@ double MHexahedron::angleShapeMeasure()
//printf("angle max =%g min =%g \n", angleMax*180/M_PI, angleMin*180/M_PI); //printf("angle max =%g min =%g \n", angleMax*180/M_PI, angleMin*180/M_PI);
} }
zeta = 1.-std::max((angleMax-0.5*M_PI)/(0.5*M_PI),(0.5*M_PI-angleMin)/(0.5*M_PI)); zeta = 1.-std::max((angleMax-0.5*M_PI)/(0.5*M_PI),(0.5*M_PI-angleMin)/(0.5*M_PI));
return zeta; return zeta;
#else #else
return 1.; return 1.;
#endif #endif
...@@ -189,30 +189,12 @@ int MHexahedronN::getNumEdgesRep() ...@@ -189,30 +189,12 @@ int MHexahedronN::getNumEdgesRep()
{ {
return 12 * CTX::instance()->mesh.numSubEdges; return 12 * CTX::instance()->mesh.numSubEdges;
} }
const nodalBasis* MHexahedron::getFunctionSpace(int o) const
{
int order = (o == -1) ? getPolynomialOrder() : o;
int nv = getNumVolumeVertices(); const nodalBasis* MHexahedron::getFunctionSpace(int order) const
{
if (order == -1) return BasisFactory::getNodalBasis(getTypeForMSH());
if ((nv == 0) && (o == -1)) { switch (order) {
switch (order) {
case 0: return BasisFactory::getNodalBasis(MSH_HEX_1);
case 1: return BasisFactory::getNodalBasis(MSH_HEX_8);
case 2: return BasisFactory::getNodalBasis(MSH_HEX_20);
case 3: return BasisFactory::getNodalBasis(MSH_HEX_32);
case 4: return BasisFactory::getNodalBasis(MSH_HEX_44);
case 5: return BasisFactory::getNodalBasis(MSH_HEX_56);
case 6: return BasisFactory::getNodalBasis(MSH_HEX_68);
case 7: return BasisFactory::getNodalBasis(MSH_HEX_80);
case 8: return BasisFactory::getNodalBasis(MSH_HEX_92);
case 9: return BasisFactory::getNodalBasis(MSH_HEX_104);
default: Msg::Error("Order %d hex function space not implemented", order); break;
}
}
else {
switch (order) {
case 0: return BasisFactory::getNodalBasis(MSH_HEX_1); case 0: return BasisFactory::getNodalBasis(MSH_HEX_1);
case 1: return BasisFactory::getNodalBasis(MSH_HEX_8); case 1: return BasisFactory::getNodalBasis(MSH_HEX_8);
case 2: return BasisFactory::getNodalBasis(MSH_HEX_27); case 2: return BasisFactory::getNodalBasis(MSH_HEX_27);
...@@ -224,34 +206,15 @@ const nodalBasis* MHexahedron::getFunctionSpace(int o) const ...@@ -224,34 +206,15 @@ const nodalBasis* MHexahedron::getFunctionSpace(int o) const
case 8: return BasisFactory::getNodalBasis(MSH_HEX_729); case 8: return BasisFactory::getNodalBasis(MSH_HEX_729);
case 9: return BasisFactory::getNodalBasis(MSH_HEX_1000); case 9: return BasisFactory::getNodalBasis(MSH_HEX_1000);
default: Msg::Error("Order %d hex function space not implemented", order); break; default: Msg::Error("Order %d hex function space not implemented", order); break;
}
} }
return 0; return NULL;
} }
const JacobianBasis* MHexahedron::getJacobianFuncSpace(int o) const const JacobianBasis* MHexahedron::getJacobianFuncSpace(int order) const
{ {
int order = (o == -1) ? getPolynomialOrder() : o; if (order == -1) return BasisFactory::getJacobianBasis(getTypeForMSH());
int nv = getNumVolumeVertices();
if ((nv == 0) && (o == -1)) { switch (order) {
switch (order) {
case 0: return BasisFactory::getJacobianBasis(MSH_HEX_1);
case 1: return BasisFactory::getJacobianBasis(MSH_HEX_8);
case 2: return BasisFactory::getJacobianBasis(MSH_HEX_20);
case 3: return BasisFactory::getJacobianBasis(MSH_HEX_32);
case 4: return BasisFactory::getJacobianBasis(MSH_HEX_44);
case 5: return BasisFactory::getJacobianBasis(MSH_HEX_56);
case 6: return BasisFactory::getJacobianBasis(MSH_HEX_68);
case 7: return BasisFactory::getJacobianBasis(MSH_HEX_80);
case 8: return BasisFactory::getJacobianBasis(MSH_HEX_92);
case 9: return BasisFactory::getJacobianBasis(MSH_HEX_104);
default: Msg::Error("Order %d hex incomplete Jacobian function space not implemented", order); break;
}
}
else {
switch (order) {
case 0: return BasisFactory::getJacobianBasis(MSH_HEX_1); case 0: return BasisFactory::getJacobianBasis(MSH_HEX_1);
case 1: return BasisFactory::getJacobianBasis(MSH_HEX_8); case 1: return BasisFactory::getJacobianBasis(MSH_HEX_8);
case 2: return BasisFactory::getJacobianBasis(MSH_HEX_27); case 2: return BasisFactory::getJacobianBasis(MSH_HEX_27);
...@@ -263,9 +226,8 @@ const JacobianBasis* MHexahedron::getJacobianFuncSpace(int o) const ...@@ -263,9 +226,8 @@ const JacobianBasis* MHexahedron::getJacobianFuncSpace(int o) const
case 8: return BasisFactory::getJacobianBasis(MSH_HEX_729); case 8: return BasisFactory::getJacobianBasis(MSH_HEX_729);
case 9: return BasisFactory::getJacobianBasis(MSH_HEX_1000); case 9: return BasisFactory::getJacobianBasis(MSH_HEX_1000);
default: Msg::Error("Order %d hex Jacobian function space not implemented", order); break; default: Msg::Error("Order %d hex Jacobian function space not implemented", order); break;
}
} }
return 0; return NULL;
} }
......
Geo/MPrism.cpp
+ 29
− 37
View file @ 2344885c
...@@ -9,7 +9,7 @@ ...@@ -9,7 +9,7 @@
#include "Context.h" #include "Context.h"
int MPrism::getVolumeSign() int MPrism::getVolumeSign()
{ {
double mat[3][3]; double mat[3][3];
mat[0][0] = _v[1]->x() - _v[0]->x(); mat[0][0] = _v[1]->x() - _v[0]->x();
mat[0][1] = _v[2]->x() - _v[0]->x(); mat[0][1] = _v[2]->x() - _v[0]->x();
...@@ -32,52 +32,44 @@ void MPrism::getIntegrationPoints(int pOrder, int *npts, IntPt **pts) ...@@ -32,52 +32,44 @@ void MPrism::getIntegrationPoints(int pOrder, int *npts, IntPt **pts)
*pts = getGQPriPts(pOrder); *pts = getGQPriPts(pOrder);
} }
const nodalBasis* MPrism::getFunctionSpace(int o) const const nodalBasis* MPrism::getFunctionSpace(int order) const
{ {
int order = (o == -1) ? getPolynomialOrder() : o; if (order == -1) return BasisFactory::getNodalBasis(getTypeForMSH());
int nv = getNumVolumeVertices();
if ((nv == 0) && (o == -1)) { switch (order) {
switch (order) {
case 0: return BasisFactory::getNodalBasis(MSH_PRI_1); case 0: return BasisFactory::getNodalBasis(MSH_PRI_1);
case 1: return BasisFactory::getNodalBasis(MSH_PRI_6); case 1: return BasisFactory::getNodalBasis(MSH_PRI_6);
case 2: return BasisFactory::getNodalBasis(MSH_PRI_18); case 2: return BasisFactory::getNodalBasis(MSH_PRI_18);
case 3: return BasisFactory::getNodalBasis(MSH_PRI_40);
case 4: return BasisFactory::getNodalBasis(MSH_PRI_75);
case 5: return BasisFactory::getNodalBasis(MSH_PRI_126);
case 6: return BasisFactory::getNodalBasis(MSH_PRI_196);
case 7: return BasisFactory::getNodalBasis(MSH_PRI_288);
case 8: return BasisFactory::getNodalBasis(MSH_PRI_405);
case 9: return BasisFactory::getNodalBasis(MSH_PRI_550);
default: Msg::Error("Order %d prism function space not implemented", order); default: Msg::Error("Order %d prism function space not implemented", order);
}
}
else {
switch (order) {
case 0: return BasisFactory::getNodalBasis(MSH_PRI_1);
case 1: return BasisFactory::getNodalBasis(MSH_PRI_6);
case 2: return BasisFactory::getNodalBasis(MSH_PRI_18);
default: Msg::Error("Order %d prism function space not implemented", order);
}
} }
return 0; return NULL;
} }
const JacobianBasis* MPrism::getJacobianFuncSpace(int o) const const JacobianBasis* MPrism::getJacobianFuncSpace(int order) const
{ {
int order = (o == -1) ? getPolynomialOrder() : o; if (order == -1) return BasisFactory::getJacobianBasis(getTypeForMSH());
int nv = getNumVolumeVertices(); switch (order) {
case 0: return BasisFactory::getJacobianBasis(MSH_PRI_1);
if ((nv == 0) && (o == -1)) {
switch (order) {
case 1: return BasisFactory::getJacobianBasis(MSH_PRI_6); case 1: return BasisFactory::getJacobianBasis(MSH_PRI_6);
case 2: return BasisFactory::getJacobianBasis(MSH_PRI_18); case 2: return BasisFactory::getJacobianBasis(MSH_PRI_18);
case 3: return BasisFactory::getJacobianBasis(MSH_PRI_40);
case 4: return BasisFactory::getJacobianBasis(MSH_PRI_75);
case 5: return BasisFactory::getJacobianBasis(MSH_PRI_126);
case 6: return BasisFactory::getJacobianBasis(MSH_PRI_196);
case 7: return BasisFactory::getJacobianBasis(MSH_PRI_288);
case 8: return BasisFactory::getJacobianBasis(MSH_PRI_405);
case 9: return BasisFactory::getJacobianBasis(MSH_PRI_550);
default: Msg::Error("Order %d prism function space not implemented", order); default: Msg::Error("Order %d prism function space not implemented", order);
}
}
else {
switch (order) {
case 1: return BasisFactory::getJacobianBasis(MSH_PRI_6);
case 2: return BasisFactory::getJacobianBasis(MSH_PRI_18);
default: Msg::Error("Order %d prism function space not implemented", order);
}
} }
return 0; return NULL;
} }
double MPrism::getInnerRadius() double MPrism::getInnerRadius()
...@@ -123,7 +115,7 @@ void MPrism::getFaceInfo(const MFace &face, int &ithFace, int &sign, int &rot) c ...@@ -123,7 +115,7 @@ void MPrism::getFaceInfo(const MFace &face, int &ithFace, int &sign, int &rot) c
} }
else { else {
MVertex *v3 = _v[faces_prism(ithFace, 3)]; MVertex *v3 = _v[faces_prism(ithFace, 3)];
if (v0 == face.getVertex(0) && v1 == face.getVertex(1) && if (v0 == face.getVertex(0) && v1 == face.getVertex(1) &&
v2 == face.getVertex(2) && v3 == face.getVertex(3)){ v2 == face.getVertex(2) && v3 == face.getVertex(3)){
sign = 1; rot = 0; return; sign = 1; rot = 0; return;
} }
...@@ -135,7 +127,7 @@ void MPrism::getFaceInfo(const MFace &face, int &ithFace, int &sign, int &rot) c ...@@ -135,7 +127,7 @@ void MPrism::getFaceInfo(const MFace &face, int &ithFace, int &sign, int &rot) c
v2 == face.getVertex(0) && v3 == face.getVertex(1)){ v2 == face.getVertex(0) && v3 == face.getVertex(1)){
sign = 1; rot = 2; return; sign = 1; rot = 2; return;
} }
if (v0 == face.getVertex(3) && v1 == face.getVertex(0) && if (v0 == face.getVertex(3) && v1 == face.getVertex(0) &&
v2 == face.getVertex(1) && v3 == face.getVertex(2)){ v2 == face.getVertex(1) && v3 == face.getVertex(2)){
sign = 1; rot = 3; return; sign = 1; rot = 3; return;
} }
...@@ -147,7 +139,7 @@ void MPrism::getFaceInfo(const MFace &face, int &ithFace, int &sign, int &rot) c ...@@ -147,7 +139,7 @@ void MPrism::getFaceInfo(const MFace &face, int &ithFace, int &sign, int &rot) c
v2 == face.getVertex(3) && v3 == face.getVertex(2)){ v2 == face.getVertex(3) && v3 == face.getVertex(2)){
sign = -1; rot = 1; return; sign = -1; rot = 1; return;
} }
if (v0 == face.getVertex(2) && v1 == face.getVertex(1) && if (v0 == face.getVertex(2) && v1 == face.getVertex(1) &&
v2 == face.getVertex(0) && v3 == face.getVertex(3)){ v2 == face.getVertex(0) && v3 == face.getVertex(3)){
sign = -1; rot = 2; return; sign = -1; rot = 2; return;
} }
...@@ -156,7 +148,7 @@ void MPrism::getFaceInfo(const MFace &face, int &ithFace, int &sign, int &rot) c ...@@ -156,7 +148,7 @@ void MPrism::getFaceInfo(const MFace &face, int &ithFace, int &sign, int &rot) c
sign = -1; rot = 3; return; sign = -1; rot = 3; return;
} }
} }
} }
Msg::Error("Could not get face information for prism %d", getNum()); Msg::Error("Could not get face information for prism %d", getNum());
} }
...@@ -451,4 +443,4 @@ int MPrism15::getNumFacesRep() { ...@@ -451,4 +443,4 @@ int MPrism15::getNumFacesRep() {
int MPrism18::getNumFacesRep() { int MPrism18::getNumFacesRep() {
return 4 * (CTX::instance()->mesh.numSubEdges * CTX::instance()->mesh.numSubEdges * 2); return 4 * (CTX::instance()->mesh.numSubEdges * CTX::instance()->mesh.numSubEdges * 2);
} }
\ No newline at end of file
Geo/MPyramid.cpp
+ 5
− 23
View file @ 2344885c
...@@ -33,29 +33,11 @@ void MPyramid::getIntegrationPoints(int pOrder, int *npts, IntPt **pts) ...@@ -33,29 +33,11 @@ void MPyramid::getIntegrationPoints(int pOrder, int *npts, IntPt **pts)
*pts = getGQPyrPts(pOrder); *pts = getGQPyrPts(pOrder);
} }
const nodalBasis* MPyramid::getFunctionSpace(int o) const const nodalBasis* MPyramid::getFunctionSpace(int order) const
{ {
int order = (o == -1) ? getPolynomialOrder() : o; if (order == -1) return BasisFactory::getNodalBasis(getTypeForMSH());
int nv = getNumVolumeVertices();
if ((nv == 0) && (o == -1)) { switch (order) {
switch (order) {
case 1: return BasisFactory::getNodalBasis(MSH_PYR_5);
case 2: return BasisFactory::getNodalBasis(MSH_PYR_14);
case 0: return BasisFactory::getNodalBasis(MSH_PYR_1);
case 3: return BasisFactory::getNodalBasis(MSH_PYR_21);
case 4: return BasisFactory::getNodalBasis(MSH_PYR_29);
case 5: return BasisFactory::getNodalBasis(MSH_PYR_37);
case 6: return BasisFactory::getNodalBasis(MSH_PYR_45);
case 7: return BasisFactory::getNodalBasis(MSH_PYR_53);
case 8: return BasisFactory::getNodalBasis(MSH_PYR_61);
case 9: return BasisFactory::getNodalBasis(MSH_PYR_69);
default: Msg::Error("Order %d pyramid function space not implemented", order);
}
}
else {
switch (order) {
case 0: return BasisFactory::getNodalBasis(MSH_PYR_1); case 0: return BasisFactory::getNodalBasis(MSH_PYR_1);
case 1: return BasisFactory::getNodalBasis(MSH_PYR_5); case 1: return BasisFactory::getNodalBasis(MSH_PYR_5);
case 2: return BasisFactory::getNodalBasis(MSH_PYR_14); case 2: return BasisFactory::getNodalBasis(MSH_PYR_14);
...@@ -67,13 +49,13 @@ const nodalBasis* MPyramid::getFunctionSpace(int o) const ...@@ -67,13 +49,13 @@ const nodalBasis* MPyramid::getFunctionSpace(int o) const
case 8: return BasisFactory::getNodalBasis(MSH_PYR_285); case 8: return BasisFactory::getNodalBasis(MSH_PYR_285);
case 9: return BasisFactory::getNodalBasis(MSH_PYR_385); case 9: return BasisFactory::getNodalBasis(MSH_PYR_385);
default: Msg::Error("Order %d pyramid function space not implemented", order); default: Msg::Error("Order %d pyramid function space not implemented", order);
}
} }
return 0; return NULL;
} }
const JacobianBasis* MPyramid::getJacobianFuncSpace(int o) const const JacobianBasis* MPyramid::getJacobianFuncSpace(int o) const
{ {
// FIXME add other order and see MPyramid::getFunctionSpace for 'design'
int order = (o == -1) ? getPolynomialOrder() : o; int order = (o == -1) ? getPolynomialOrder() : o;
switch (order) { switch (order) {
......
Geo/MQuadrangle.cpp
+ 12
− 46
View file @ 2344885c
...@@ -17,27 +17,10 @@ ...@@ -17,27 +17,10 @@
#define SQU(a) ((a)*(a)) #define SQU(a) ((a)*(a))
const nodalBasis* MQuadrangle::getFunctionSpace(int o) const const nodalBasis* MQuadrangle::getFunctionSpace(int order) const
{ {
int order = (o == -1) ? getPolynomialOrder() : o; if (order == -1) return BasisFactory::getNodalBasis(getTypeForMSH());
int nf = getNumFaceVertices();
if ((nf == 0) && (o == -1)) {
switch (order) {
case 0: return BasisFactory::getNodalBasis(MSH_QUA_1);
case 1: return BasisFactory::getNodalBasis(MSH_QUA_4);
case 2: return BasisFactory::getNodalBasis(MSH_QUA_8);
case 3: return BasisFactory::getNodalBasis(MSH_QUA_12);
case 4: return BasisFactory::getNodalBasis(MSH_QUA_16I);
case 5: return BasisFactory::getNodalBasis(MSH_QUA_20);
case 6: return BasisFactory::getNodalBasis(MSH_QUA_24);
case 7: return BasisFactory::getNodalBasis(MSH_QUA_28);
case 8: return BasisFactory::getNodalBasis(MSH_QUA_32);
case 9: return BasisFactory::getNodalBasis(MSH_QUA_36I);
case 10: return BasisFactory::getNodalBasis(MSH_QUA_40);
}
}
switch (order) { switch (order) {
case 0: return BasisFactory::getNodalBasis(MSH_QUA_1); case 0: return BasisFactory::getNodalBasis(MSH_QUA_1);
case 1: return BasisFactory::getNodalBasis(MSH_QUA_4); case 1: return BasisFactory::getNodalBasis(MSH_QUA_4);
...@@ -52,29 +35,12 @@ const nodalBasis* MQuadrangle::getFunctionSpace(int o) const ...@@ -52,29 +35,12 @@ const nodalBasis* MQuadrangle::getFunctionSpace(int o) const
case 10: return BasisFactory::getNodalBasis(MSH_QUA_121); case 10: return BasisFactory::getNodalBasis(MSH_QUA_121);
default: Msg::Error("Order %d quadrangle function space not implemented", order); default: Msg::Error("Order %d quadrangle function space not implemented", order);
} }
return 0; return NULL;
} }
const JacobianBasis* MQuadrangle::getJacobianFuncSpace(int o) const const JacobianBasis* MQuadrangle::getJacobianFuncSpace(int order) const
{ {
int order = (o == -1) ? getPolynomialOrder() : o; if (order == -1) return BasisFactory::getJacobianBasis(getTypeForMSH());
int nf = getNumFaceVertices();
if ((nf == 0) && (o == -1)) {
switch (order) {
case 1: return BasisFactory::getJacobianBasis(MSH_QUA_4);
case 2: return BasisFactory::getJacobianBasis(MSH_QUA_8);
case 3: return BasisFactory::getJacobianBasis(MSH_QUA_12);
case 4: return BasisFactory::getJacobianBasis(MSH_QUA_16I);
case 5: return BasisFactory::getJacobianBasis(MSH_QUA_20);
case 6: return BasisFactory::getJacobianBasis(MSH_QUA_24);
case 7: return BasisFactory::getJacobianBasis(MSH_QUA_28);
case 8: return BasisFactory::getJacobianBasis(MSH_QUA_32);
case 9: return BasisFactory::getJacobianBasis(MSH_QUA_36I);
case 10: return BasisFactory::getJacobianBasis(MSH_QUA_40);
}
}
switch (order) { switch (order) {
case 1: return BasisFactory::getJacobianBasis(MSH_QUA_4); case 1: return BasisFactory::getJacobianBasis(MSH_QUA_4);
case 2: return BasisFactory::getJacobianBasis(MSH_QUA_9); case 2: return BasisFactory::getJacobianBasis(MSH_QUA_9);
...@@ -234,7 +200,7 @@ void MQuadrangle::getIntegrationPoints(int pOrder, int *npts, IntPt **pts) ...@@ -234,7 +200,7 @@ void MQuadrangle::getIntegrationPoints(int pOrder, int *npts, IntPt **pts)
double MQuadrangle::etaShapeMeasure() double MQuadrangle::etaShapeMeasure()
{ {
double AR = 1;//(minEdge()/maxEdge()); double AR = 1;//(minEdge()/maxEdge());
SVector3 v01 (_v[1]->x()-_v[0]->x(),_v[1]->y()-_v[0]->y(),_v[1]->z()-_v[0]->z()); SVector3 v01 (_v[1]->x()-_v[0]->x(),_v[1]->y()-_v[0]->y(),_v[1]->z()-_v[0]->z());
SVector3 v12 (_v[2]->x()-_v[1]->x(),_v[2]->y()-_v[1]->y(),_v[2]->z()-_v[1]->z()); SVector3 v12 (_v[2]->x()-_v[1]->x(),_v[2]->y()-_v[1]->y(),_v[2]->z()-_v[1]->z());
SVector3 v23 (_v[3]->x()-_v[2]->x(),_v[3]->y()-_v[2]->y(),_v[3]->z()-_v[2]->z()); SVector3 v23 (_v[3]->x()-_v[2]->x(),_v[3]->y()-_v[2]->y(),_v[3]->z()-_v[2]->z());
...@@ -268,12 +234,12 @@ double MQuadrangle::etaShapeMeasure() ...@@ -268,12 +234,12 @@ double MQuadrangle::etaShapeMeasure()
} }
/// a shape measure for quadrangles /// a shape measure for quadrangles
/// assume (for now) 2D elements -- /// assume (for now) 2D elements --
/// sf = (1 \pm xi) (1 \pm eta) / 4 /// sf = (1 \pm xi) (1 \pm eta) / 4
/// dsf_xi = \pm (1 \pm eta) / 4 /// dsf_xi = \pm (1 \pm eta) / 4
/// 1 + eta , -(1+eta) , -(1-eta), 1-eta /// 1 + eta , -(1+eta) , -(1-eta), 1-eta
/// dsf_eta = \pm (1 \pm xi) / 4 /// dsf_eta = \pm (1 \pm xi) / 4
/// 1 + xi , 1 - xi , -(1-xi), -(1+xi) /// 1 + xi , 1 - xi , -(1-xi), -(1+xi)
double MQuadrangle::gammaShapeMeasure(){ double MQuadrangle::gammaShapeMeasure(){
return etaShapeMeasure(); return etaShapeMeasure();
} }
...@@ -366,10 +332,10 @@ double MQuadrangle::getInnerRadius() ...@@ -366,10 +332,10 @@ double MQuadrangle::getInnerRadius()
} }
return R; return R;
#else // HAVE_LAPACK #else // HAVE_LAPACK
// Default implementation. Not sure that the following give exactly // Default implementation. Not sure that the following give exactly
// the same value as the HAVE_LAPACK case ! // the same value as the HAVE_LAPACK case !
// but same value for a square // but same value for a square
// Mid-point of each edge of the quadrangle // Mid-point of each edge of the quadrangle
SPoint3 A(_v[0]->x()+_v[1]->x(),_v[0]->y()+_v[1]->y(),_v[0]->z()+_v[1]->z()); SPoint3 A(_v[0]->x()+_v[1]->x(),_v[0]->y()+_v[1]->y(),_v[0]->z()+_v[1]->z());
SPoint3 B(_v[1]->x()+_v[2]->x(),_v[1]->y()+_v[2]->y(),_v[1]->z()+_v[2]->z()); SPoint3 B(_v[1]->x()+_v[2]->x(),_v[1]->y()+_v[2]->y(),_v[1]->z()+_v[2]->z());
......
Geo/MTetrahedron.cpp
+ 8
− 47
View file @ 2344885c
...@@ -100,30 +100,11 @@ void MTetrahedron::xyz2uvw(double xyz[3], double uvw[3]) const ...@@ -100,30 +100,11 @@ void MTetrahedron::xyz2uvw(double xyz[3], double uvw[3]) const
sys3x3(mat, b, uvw, &det); sys3x3(mat, b, uvw, &det);
} }
const nodalBasis* MTetrahedron::getFunctionSpace(int o) const const nodalBasis* MTetrahedron::getFunctionSpace(int order) const
{ {
int order = (o == -1) ? getPolynomialOrder() : o; if (order == -1) return BasisFactory::getNodalBasis(getTypeForMSH());
int nv = getNumVolumeVertices(); switch (order) {
if ((nv == 0) && (o == -1)) {
switch (order) {
case 0: return BasisFactory::getNodalBasis(MSH_TET_1);
case 1: return BasisFactory::getNodalBasis(MSH_TET_4);
case 2: return BasisFactory::getNodalBasis(MSH_TET_10);
case 3: return BasisFactory::getNodalBasis(MSH_TET_16); // not just nv==0
case 4: return BasisFactory::getNodalBasis(MSH_TET_22);
case 5: return BasisFactory::getNodalBasis(MSH_TET_28);
case 6: return BasisFactory::getNodalBasis(MSH_TET_34);
case 7: return BasisFactory::getNodalBasis(MSH_TET_40);
case 8: return BasisFactory::getNodalBasis(MSH_TET_46);
case 9: return BasisFactory::getNodalBasis(MSH_TET_52);
case 10: return BasisFactory::getNodalBasis(MSH_TET_58);
default: Msg::Error("Order %d tetrahedron function space not implemented", order);
}
}
else {
switch (order) {
case 0: return BasisFactory::getNodalBasis(MSH_TET_1); case 0: return BasisFactory::getNodalBasis(MSH_TET_1);
case 1: return BasisFactory::getNodalBasis(MSH_TET_4); case 1: return BasisFactory::getNodalBasis(MSH_TET_4);
case 2: return BasisFactory::getNodalBasis(MSH_TET_10); case 2: return BasisFactory::getNodalBasis(MSH_TET_10);
...@@ -136,34 +117,15 @@ const nodalBasis* MTetrahedron::getFunctionSpace(int o) const ...@@ -136,34 +117,15 @@ const nodalBasis* MTetrahedron::getFunctionSpace(int o) const
case 9: return BasisFactory::getNodalBasis(MSH_TET_220); case 9: return BasisFactory::getNodalBasis(MSH_TET_220);
case 10: return BasisFactory::getNodalBasis(MSH_TET_286); case 10: return BasisFactory::getNodalBasis(MSH_TET_286);
default: Msg::Error("Order %d tetrahedron function space not implemented", order); default: Msg::Error("Order %d tetrahedron function space not implemented", order);
}
} }
return 0; return NULL;
} }
const JacobianBasis* MTetrahedron::getJacobianFuncSpace(int o) const const JacobianBasis* MTetrahedron::getJacobianFuncSpace(int order) const
{ {
int order = (o == -1) ? getPolynomialOrder() : o; if (order == -1) return BasisFactory::getJacobianBasis(getTypeForMSH());
int nv = getNumVolumeVertices(); switch (order) {
if ((nv == 0) && (o == -1)) {
switch (order) {
case 1: return BasisFactory::getJacobianBasis(MSH_TET_4);
case 2: return BasisFactory::getJacobianBasis(MSH_TET_10);
case 3: return BasisFactory::getJacobianBasis(MSH_TET_16); // not just nv==0
case 4: return BasisFactory::getJacobianBasis(MSH_TET_22);
case 5: return BasisFactory::getJacobianBasis(MSH_TET_28);
case 6: return BasisFactory::getJacobianBasis(MSH_TET_34);
case 7: return BasisFactory::getJacobianBasis(MSH_TET_40);
case 8: return BasisFactory::getJacobianBasis(MSH_TET_46);
case 9: return BasisFactory::getJacobianBasis(MSH_TET_52);
case 10: return BasisFactory::getJacobianBasis(MSH_TET_58);
default: Msg::Error("Order %d tetrahedron function space not implemented", order);
}
}
else {
switch (order) {
case 1: return BasisFactory::getJacobianBasis(MSH_TET_4); case 1: return BasisFactory::getJacobianBasis(MSH_TET_4);
case 2: return BasisFactory::getJacobianBasis(MSH_TET_10); case 2: return BasisFactory::getJacobianBasis(MSH_TET_10);
case 3: return BasisFactory::getJacobianBasis(MSH_TET_20); case 3: return BasisFactory::getJacobianBasis(MSH_TET_20);
...@@ -175,9 +137,8 @@ const JacobianBasis* MTetrahedron::getJacobianFuncSpace(int o) const ...@@ -175,9 +137,8 @@ const JacobianBasis* MTetrahedron::getJacobianFuncSpace(int o) const
case 9: return BasisFactory::getJacobianBasis(MSH_TET_220); case 9: return BasisFactory::getJacobianBasis(MSH_TET_220);
case 10: return BasisFactory::getJacobianBasis(MSH_TET_286); case 10: return BasisFactory::getJacobianBasis(MSH_TET_286);
default: Msg::Error("Order %d tetrahedron function space not implemented", order); default: Msg::Error("Order %d tetrahedron function space not implemented", order);
}
} }
return 0; return NULL;
} }
int MTetrahedron10::getNumEdgesRep(){ return 6 * CTX::instance()->mesh.numSubEdges; } int MTetrahedron10::getNumEdgesRep(){ return 6 * CTX::instance()->mesh.numSubEdges; }
......
Geo/MTriangle.cpp
+ 8
− 49
View file @ 2344885c
...@@ -119,30 +119,11 @@ double MTriangle::gammaShapeMeasure() ...@@ -119,30 +119,11 @@ double MTriangle::gammaShapeMeasure()
uvw[2] = 0.; uvw[2] = 0.;
} }
const nodalBasis* MTriangle::getFunctionSpace(int o) const const nodalBasis* MTriangle::getFunctionSpace(int order) const
{ {
int order = (o == -1) ? getPolynomialOrder() : o; if (order == -1) return BasisFactory::getNodalBasis(getTypeForMSH());
int nf = getNumFaceVertices(); switch (order) {
if ((nf == 0) && (o == -1)) {
switch (order) {
case 0: return BasisFactory::getNodalBasis(MSH_TRI_1);
case 1: return BasisFactory::getNodalBasis(MSH_TRI_3);
case 2: return BasisFactory::getNodalBasis(MSH_TRI_6);
case 3: return BasisFactory::getNodalBasis(MSH_TRI_9);
case 4: return BasisFactory::getNodalBasis(MSH_TRI_12);
case 5: return BasisFactory::getNodalBasis(MSH_TRI_15I);
case 6: return BasisFactory::getNodalBasis(MSH_TRI_18);
case 7: return BasisFactory::getNodalBasis(MSH_TRI_21I);
case 8: return BasisFactory::getNodalBasis(MSH_TRI_24);
case 9: return BasisFactory::getNodalBasis(MSH_TRI_27);
case 10: return BasisFactory::getNodalBasis(MSH_TRI_30);
default: Msg::Error("Order %d triangle incomplete function space not implemented", order);
}
}
else {
switch (order) {
case 0: return BasisFactory::getNodalBasis(MSH_TRI_1); case 0: return BasisFactory::getNodalBasis(MSH_TRI_1);
case 1: return BasisFactory::getNodalBasis(MSH_TRI_3); case 1: return BasisFactory::getNodalBasis(MSH_TRI_3);
case 2: return BasisFactory::getNodalBasis(MSH_TRI_6); case 2: return BasisFactory::getNodalBasis(MSH_TRI_6);
...@@ -155,36 +136,15 @@ const nodalBasis* MTriangle::getFunctionSpace(int o) const ...@@ -155,36 +136,15 @@ const nodalBasis* MTriangle::getFunctionSpace(int o) const
case 9: return BasisFactory::getNodalBasis(MSH_TRI_55); case 9: return BasisFactory::getNodalBasis(MSH_TRI_55);
case 10: return BasisFactory::getNodalBasis(MSH_TRI_66); case 10: return BasisFactory::getNodalBasis(MSH_TRI_66);
default: Msg::Error("Order %d triangle function space not implemented", order); default: Msg::Error("Order %d triangle function space not implemented", order);
}
} }
return 0; return NULL;
} }
const JacobianBasis* MTriangle::getJacobianFuncSpace(int o) const const JacobianBasis* MTriangle::getJacobianFuncSpace(int order) const
{ {
if (order == -1) return BasisFactory::getNodalBasis(getTypeForMSH());
int order = (o == -1) ? getPolynomialOrder() : o; switch (order) {
int nf = getNumFaceVertices();
if ((nf == 0) && (o == -1)) {
switch (order) {
case 0: return BasisFactory::getJacobianBasis(MSH_TRI_1);
case 1: return BasisFactory::getJacobianBasis(MSH_TRI_3);
case 2: return BasisFactory::getJacobianBasis(MSH_TRI_6);
case 3: return BasisFactory::getJacobianBasis(MSH_TRI_9);
case 4: return BasisFactory::getJacobianBasis(MSH_TRI_12);
case 5: return BasisFactory::getJacobianBasis(MSH_TRI_15I);
case 6: return BasisFactory::getJacobianBasis(MSH_TRI_18);
case 7: return BasisFactory::getJacobianBasis(MSH_TRI_21I);
case 8: return BasisFactory::getJacobianBasis(MSH_TRI_24);
case 9: return BasisFactory::getJacobianBasis(MSH_TRI_27);
case 10: return BasisFactory::getJacobianBasis(MSH_TRI_30);
default: Msg::Error("Order %d triangle incomplete function space not implemented", order);
}
}
else {
switch (order) {
case 1: return BasisFactory::getJacobianBasis(MSH_TRI_3); case 1: return BasisFactory::getJacobianBasis(MSH_TRI_3);
case 2: return BasisFactory::getJacobianBasis(MSH_TRI_6); case 2: return BasisFactory::getJacobianBasis(MSH_TRI_6);
case 3: return BasisFactory::getJacobianBasis(MSH_TRI_10); case 3: return BasisFactory::getJacobianBasis(MSH_TRI_10);
...@@ -196,9 +156,8 @@ const JacobianBasis* MTriangle::getJacobianFuncSpace(int o) const ...@@ -196,9 +156,8 @@ const JacobianBasis* MTriangle::getJacobianFuncSpace(int o) const
case 9: return BasisFactory::getJacobianBasis(MSH_TRI_55); case 9: return BasisFactory::getJacobianBasis(MSH_TRI_55);
case 10: return BasisFactory::getJacobianBasis(MSH_TRI_66); case 10: return BasisFactory::getJacobianBasis(MSH_TRI_66);
default: Msg::Error("Order %d triangle function space not implemented", order); default: Msg::Error("Order %d triangle function space not implemented", order);
}
} }
return 0; return NULL;
} }
int MTriangleN::getNumEdgesRep(){ return 3 * CTX::instance()->mesh.numSubEdges; } int MTriangleN::getNumEdgesRep(){ return 3 * CTX::instance()->mesh.numSubEdges; }
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment