diff --git a/FunctionSpace/BasisGenerator.cpp b/FunctionSpace/BasisGenerator.cpp index 089c726903b9231891a151ddf881a7c5cbf3411a..58376a0404fc47b421613e1e141482ea2f84a7fd 100644 --- a/FunctionSpace/BasisGenerator.cpp +++ b/FunctionSpace/BasisGenerator.cpp @@ -2,6 +2,8 @@ #include "GmshDefines.h" #include "Exception.h" +#include "LineNodeBasis.h" + #include "QuadNodeBasis.h" #include "QuadEdgeBasis.h" @@ -26,6 +28,7 @@ Basis* BasisGenerator::generate(int elementType, int basisType, int order){ switch(elementType){ + case TYPE_LIN: return linGen(basisType, order); case TYPE_TRI: return triGen(basisType, order); case TYPE_QUA: return quaGen(basisType, order); case TYPE_TET: return tetGen(basisType, order); @@ -36,6 +39,18 @@ Basis* BasisGenerator::generate(int elementType, } } +Basis* BasisGenerator::linGen(int basisType, + int order){ + switch(basisType){ + case 0: return new LineNodeBasis(order); + case 1: throw Exception("1-form not implemented on Lines"); + case 2: throw Exception("2-form not implemented on Lines"); + case 3: throw Exception("3-form not implemented on Lines"); + + default: throw Exception("There is no %d-form", basisType); + } +} + Basis* BasisGenerator::triGen(int basisType, int order){ switch(basisType){ diff --git a/FunctionSpace/BasisGenerator.h b/FunctionSpace/BasisGenerator.h index 08af96b40d6bb0ce74dee823f9eb93b2f567757c..1f071ea3bcd59ee30a77bb93459acc0a1a022b95 100644 --- a/FunctionSpace/BasisGenerator.h +++ b/FunctionSpace/BasisGenerator.h @@ -24,6 +24,7 @@ class BasisGenerator{ int basisType, int order); + static Basis* linGen(int basisType, int order); static Basis* triGen(int basisType, int order); static Basis* quaGen(int basisType, int order); static Basis* tetGen(int basisType, int order); @@ -56,6 +57,7 @@ class BasisGenerator{ @em instantiated Basis @note Element types are: + @li @c TYPE_LIN for Lines @li @c TYPE_TRI for Triangles @li @c TYPE_QUA for Quadrangles @li @c TYPE_TET for Tetrahedrons @@ -68,6 +70,23 @@ class BasisGenerator{ @li @c 3 for 3-Form ** + @fn BasisGenerator::linGen + @param basisType The Basis type + @param order The order or the requested Basis + + This method will @em instanciate the requested Basis, + with a @em Line for support + + @return Returns a @em pointer to a newly + @em instantiated Basis + + @note Basis types are: + @li @c 0 for 0-Form + @li @c 1 for 1-Form + @li @c 2 for 2-Form + @li @c 3 for 3-Form + ** + @fn BasisGenerator::triGen @param basisType The Basis type @param order The order or the requested Basis diff --git a/FunctionSpace/FunctionSpace.cpp b/FunctionSpace/FunctionSpace.cpp index a51dd7f8755acd977e98fae1006d8a47871a0864..40d9c334336f40d3ecacaf12e30dba0228663ba8 100644 --- a/FunctionSpace/FunctionSpace.cpp +++ b/FunctionSpace/FunctionSpace.cpp @@ -321,34 +321,40 @@ FunctionSpace::locBasis(const MElement& element, // Edge Based // Number of basis function *per* edge // --> should always be an integer ! - const unsigned int nEdge = edgeClosure.size(); - const unsigned int nFPerEdge = nFEdge / nEdge; - unsigned int fEdge = 0; + const unsigned int nEdge = edgeClosure.size(); - for(unsigned int j = 0; j < nFPerEdge; j++){ - for(unsigned int k = 0; k < nEdge; k++){ - fun[i] = - &basis.getEdgeFunction(edgeClosure[k], fEdge); - - fEdge++; - i++; + if(nEdge){ + const unsigned int nFPerEdge = nFEdge / nEdge; + unsigned int fEdge = 0; + + for(unsigned int j = 0; j < nFPerEdge; j++){ + for(unsigned int k = 0; k < nEdge; k++){ + fun[i] = + &basis.getEdgeFunction(edgeClosure[k], fEdge); + + fEdge++; + i++; + } } } // Face Based // Number of basis function *per* face // --> should always be an integer ! - const unsigned int nFace = faceClosure.size(); - const unsigned int nFPerFace = nFFace / nFace; - unsigned int fFace = 0; - - for(unsigned int j = 0; j < nFPerFace; j++){ - for(unsigned int k = 0; k < nFace; k++){ - fun[i] = - &basis.getFaceFunction(faceClosure[k], fFace); - - fFace++; - i++; + const unsigned int nFace = faceClosure.size(); + + if(nFace){ + const unsigned int nFPerFace = nFFace / nFace; + unsigned int fFace = 0; + + for(unsigned int j = 0; j < nFPerFace; j++){ + for(unsigned int k = 0; k < nFace; k++){ + fun[i] = + &basis.getFaceFunction(faceClosure[k], fFace); + + fFace++; + i++; + } } } @@ -390,16 +396,19 @@ FunctionSpace::locBasis(const MElement& element, // Number of basis function *per* edge // --> should always be an integer ! const unsigned int nEdge = edgeClosure.size(); - const unsigned int nFPerEdge = nFEdge / nEdge; - unsigned int fEdge = 0; - for(unsigned int j = 0; j < nFPerEdge; j++){ - for(unsigned int k = 0; k < nEdge; k++){ - fun[i] = - &basis.getEdgeFunction(edgeClosure[k], fEdge); - - fEdge++; - i++; + if(nEdge){ + const unsigned int nFPerEdge = nFEdge / nEdge; + unsigned int fEdge = 0; + + for(unsigned int j = 0; j < nFPerEdge; j++){ + for(unsigned int k = 0; k < nEdge; k++){ + fun[i] = + &basis.getEdgeFunction(edgeClosure[k], fEdge); + + fEdge++; + i++; + } } } @@ -407,16 +416,19 @@ FunctionSpace::locBasis(const MElement& element, // Number of basis function *per* face // --> should always be an integer ! const unsigned int nFace = faceClosure.size(); - const unsigned int nFPerFace = nFFace / nFace; - unsigned int fFace = 0; - for(unsigned int j = 0; j < nFPerFace; j++){ - for(unsigned int k = 0; k < nFace; k++){ - fun[i] = - &basis.getFaceFunction(faceClosure[k], fFace); - - fFace++; - i++; + if(nFace){ + const unsigned int nFPerFace = nFFace / nFace; + unsigned int fFace = 0; + + for(unsigned int j = 0; j < nFPerFace; j++){ + for(unsigned int k = 0; k < nFace; k++){ + fun[i] = + &basis.getFaceFunction(faceClosure[k], fFace); + + fFace++; + i++; + } } } diff --git a/FunctionSpace/LineNodeBasis.cpp b/FunctionSpace/LineNodeBasis.cpp index 94e3afcea41e045a15aecb3f4f5f3eb256a35d0d..762d68b31fb9ce4054d69b2352ceac5af706c9e1 100644 --- a/FunctionSpace/LineNodeBasis.cpp +++ b/FunctionSpace/LineNodeBasis.cpp @@ -36,16 +36,14 @@ LineNodeBasis::LineNodeBasis(int order){ // Vertex Based (Lagrange) // (*node)[0] = - new Polynomial((Polynomial(1, 0, 0, 0) - - Polynomial(1, 1, 0, 0)) * - 0.5); + new Polynomial(Polynomial(0.5, 0, 0, 0) - + Polynomial(0.5, 1, 0, 0)); (*node)[1] = - new Polynomial((Polynomial(1, 0, 0, 0) + - Polynomial(1, 1, 0, 0)) * - 0.5); - - + new Polynomial(Polynomial(0.5, 0, 0, 0) + + Polynomial(0.5, 1, 0, 0)); + + // Edge Based // const int permutation[2] = {0, 1}; const Polynomial x[2] = {