diff --git a/FunctionSpace/BasisGenerator.cpp b/FunctionSpace/BasisGenerator.cpp index 2b93b9cb4cc80ac0398efd176e81735503c36b27..454c403a1199e524fa09a6f50e14f9658274e106 100644 --- a/FunctionSpace/BasisGenerator.cpp +++ b/FunctionSpace/BasisGenerator.cpp @@ -35,11 +35,11 @@ Basis* BasisGenerator::generate(int elementType, Basis* BasisGenerator::TriGen(int basisType, int order){ switch(basisType){ - case 0: + case 0: return new TriNodeBasis(order); + case 1: if (order == 0) return new TriNedelecBasis(); - else return new TriNodeBasis(order); + else return new TriEdgeBasis(order); - case 1: return new TriEdgeBasis(order); case 2: throw Exception("2-form not implemented on Triangles"); case 3: throw Exception("3-form not implemented on Triangles"); diff --git a/FunctionSpace/CMakeLists.txt b/FunctionSpace/CMakeLists.txt index fbda89a27f384bf3c98ce1b418cdc831c95d1707..d42a0fc7162bdba479050bdabd37bf5ebec4ccd2 100644 --- a/FunctionSpace/CMakeLists.txt +++ b/FunctionSpace/CMakeLists.txt @@ -22,6 +22,10 @@ set(SRC HexEdgeBasis.cpp FunctionSpace.cpp + FunctionSpaceScalar.cpp + FunctionSpaceVector.cpp + FunctionSpaceNode.cpp + FunctionSpaceEdge.cpp ) file(GLOB HDR RELATIVE ${CMAKE_CURRENT_SOURCE_DIR} *.h) diff --git a/FunctionSpace/FunctionSpace.cpp b/FunctionSpace/FunctionSpace.cpp index 44ed9547576d591dfc9d2a095337f2894c111625..e88528b4556ca54567b8f7425977b9a9595d33fd 100644 --- a/FunctionSpace/FunctionSpace.cpp +++ b/FunctionSpace/FunctionSpace.cpp @@ -3,8 +3,15 @@ using namespace std; -FunctionSpace::FunctionSpace(const GroupOfElement& goe, - int basisType, int order){ +FunctionSpace::FunctionSpace(void){ +} + +FunctionSpace::~FunctionSpace(void){ + delete basis; +} + +void FunctionSpace::build(const GroupOfElement& goe, + int basisType, int order){ // Save GroupOfElement & Mesh // this->goe = &goe; this->mesh = &(goe.getMesh()); @@ -41,10 +48,6 @@ FunctionSpace::FunctionSpace(const GroupOfElement& goe, fPerCell = basis->getNCellBased(); // We always got 1 cell } -FunctionSpace::~FunctionSpace(void){ - delete basis; -} - vector<Dof> FunctionSpace::getKeys(const MElement& elem) const{ // Const_Cast // MElement& element = const_cast<MElement&>(elem); @@ -94,9 +97,53 @@ vector<Dof> FunctionSpace::getKeys(const MElement& elem) const{ } } + // Add Edge Based Dof // + for(int i = 0; i < nEdge; i++){ + for(int j = 0; j < nFEdge; j++){ + myDof[it].setDof(mesh->getGlobalId(edge[i]), j); + it++; + } + } + /* + // Add Face Based Dof // + for(int i = 0; i < nFace; i++){ + for(int j = 0; j < nFFace; j++){ + myDof[it].setDof(mesh->getGlobalId(face[i]), j); + it++; + } + } + */ + // Add Cell Based Dof // + for(int j = 0; j < nFCell; j++){ + myDof[it].setDof(mesh->getGlobalId(element), j); + it++; + } + + return myDof; } +int FunctionSpace::getElementType(const Dof& dof) const{ + const unsigned int type = dof.getType(); + + if(type < fPerVertex) // Vertex Based + return 0; + + else if(type < fPerEdge) // Edge Based + return 1; + + else if(type < fPerFace) // Face Based + return 2; + + else // Cell Based + return 3; +} + +int FunctionSpace::getElementGlobalId(const Dof& dof) const{ + return dof.getEntity(); +} + + /* #include "Polynomial.h" #include "BasisScalar.h" diff --git a/FunctionSpace/FunctionSpace.h b/FunctionSpace/FunctionSpace.h index fb0f30dd3e19fe7838766eca860a397ec9c145f4..770d3730d25c32d48af0a3e9f7e2afc1d70f54dd 100644 --- a/FunctionSpace/FunctionSpace.h +++ b/FunctionSpace/FunctionSpace.h @@ -30,37 +30,37 @@ class FunctionSpace{ const GroupOfElement* goe; const Basis* basis; - int fPerVertex; - int fPerEdge; - int fPerFace; - int fPerCell; - int type; + unsigned int fPerVertex; + unsigned int fPerEdge; + unsigned int fPerFace; + unsigned int fPerCell; + unsigned int type; public: - FunctionSpace(const GroupOfElement& goe, - int basisType, int order); - virtual ~FunctionSpace(void); const GroupOfElement& getSupport(void) const; const Basis& getBasis(const MElement& element) const; int getType(void) const; - int getNFunctionPerVertex(const MElement& element) const; - int getNFunctionPerEdge(const MElement& element) const; - int getNFunctionPerFace(const MElement& element) const; - int getNFunctionPerCell(const MElement& element) const; + unsigned int getNFunctionPerVertex(const MElement& element) const; + unsigned int getNFunctionPerEdge(const MElement& element) const; + unsigned int getNFunctionPerFace(const MElement& element) const; + unsigned int getNFunctionPerCell(const MElement& element) const; std::vector<Dof> getKeys(const MElement& element) const; std::vector<Dof> getKeys(const MVertex& element) const; std::vector<Dof> getKeys(const MEdge& element) const; std::vector<Dof> getKeys(const MFace& element) const; - /* - void interpolateAtNodes(const MElement& element, - const std::vector<double>& coef, - std::vector<double>& nodeValue) const; - */ + int getElementType(const Dof& dof) const; + int getElementGlobalId(const Dof& dof) const; + + protected: + FunctionSpace(); + + void build(const GroupOfElement& goe, + int basisType, int order); }; ////////////////////// @@ -79,19 +79,19 @@ inline int FunctionSpace::getType(void) const{ return type; } -inline int FunctionSpace::getNFunctionPerVertex(const MElement& element) const{ +inline unsigned int FunctionSpace::getNFunctionPerVertex(const MElement& element) const{ return fPerVertex; } -inline int FunctionSpace::getNFunctionPerEdge(const MElement& element) const{ +inline unsigned int FunctionSpace::getNFunctionPerEdge(const MElement& element) const{ return fPerEdge; } -inline int FunctionSpace::getNFunctionPerFace(const MElement& element) const{ +inline unsigned int FunctionSpace::getNFunctionPerFace(const MElement& element) const{ return fPerFace; } -inline int FunctionSpace::getNFunctionPerCell(const MElement& element) const{ +inline unsigned int FunctionSpace::getNFunctionPerCell(const MElement& element) const{ return fPerCell; } diff --git a/FunctionSpace/FunctionSpaceEdge.cpp b/FunctionSpace/FunctionSpaceEdge.cpp new file mode 100644 index 0000000000000000000000000000000000000000..698822f77297eb8596f8a473749e2677dae5a48c --- /dev/null +++ b/FunctionSpace/FunctionSpaceEdge.cpp @@ -0,0 +1,61 @@ +#include <vector> + +#include "BasisVector.h" +#include "Mapper.h" +#include "FunctionSpaceEdge.h" + +using namespace std; + +FunctionSpaceEdge::FunctionSpaceEdge(const GroupOfElement& goe, + int order){ + // Build 1Form Basis // + build(goe, 1, order); +} + +FunctionSpaceEdge::~FunctionSpaceEdge(void){ +} + +#include <cstdio> + +fullVector<double> FunctionSpaceEdge:: +interpolate(const MElement& element, + const std::vector<double>& coef, + const fullVector<double>& xyz) const{ + + MElement& eelement = + const_cast<MElement&>(element); + + const vector<vector<Polynomial> >& fun = + static_cast<const BasisVector&>(*basis).getBasis(); + + const unsigned int nFun = fun.size(); + + fullVector<double> val(3); + + fullMatrix<double> invJac(3, 3); + eelement.getJacobian(xyz(0), xyz(1), xyz(2), invJac); + invJac.invertInPlace(); + + fullVector<double> origin(3); + origin(0) = eelement.getVertex(0)->x(); + origin(1) = eelement.getVertex(0)->y(); + origin(2) = eelement.getVertex(0)->z(); + + fullVector<double> uvw = Mapper::invMap(xyz, origin, invJac); + + printf("%g, %g, %g\n", uvw(0), uvw(1), uvw(2)); + + //val(0) = 1; val(1) = 1; val(2) = 0; + + + for(unsigned int i = 0; i < nFun; i++){ + fullVector<double> vi = + Mapper::grad(Polynomial::at(fun[i], uvw(0), uvw(1), uvw(2)), + invJac); + + vi.scale(coef[i]); + val.axpy(vi, 1); + } + + return val; +} diff --git a/FunctionSpace/FunctionSpaceEdge.h b/FunctionSpace/FunctionSpaceEdge.h new file mode 100644 index 0000000000000000000000000000000000000000..877289a27cdf81fe3e1ee33dcd25cf7b315b0d08 --- /dev/null +++ b/FunctionSpace/FunctionSpaceEdge.h @@ -0,0 +1,18 @@ +#ifndef _FUNCTIONSPACEEDGE_H_ +#define _FUNCTIONSPACEEDGE_H_ + +#include "FunctionSpaceVector.h" + +class FunctionSpaceEdge : public FunctionSpaceVector{ + public: + FunctionSpaceEdge(const GroupOfElement& goe, int order); + + virtual ~FunctionSpaceEdge(void); + + virtual fullVector<double> + interpolate(const MElement& element, + const std::vector<double>& coef, + const fullVector<double>& xyz) const; +}; + +#endif diff --git a/FunctionSpace/FunctionSpaceNode.cpp b/FunctionSpace/FunctionSpaceNode.cpp new file mode 100644 index 0000000000000000000000000000000000000000..3f0dd595cd9beacf14c2df0a0a5c1f7f5b838482 --- /dev/null +++ b/FunctionSpace/FunctionSpaceNode.cpp @@ -0,0 +1,18 @@ +#include "FunctionSpaceNode.h" + +FunctionSpaceNode::FunctionSpaceNode(const GroupOfElement& goe, + int order){ + // Build 0Form Basis // + build(goe, 0, order); +} + +FunctionSpaceNode::~FunctionSpaceNode(void){ +} + +double FunctionSpaceNode:: +interpolate(const MElement& element, + const std::vector<double>& coef, + const fullVector<double>& xyz) const{ + + return 42; +} diff --git a/FunctionSpace/FunctionSpaceNode.h b/FunctionSpace/FunctionSpaceNode.h new file mode 100644 index 0000000000000000000000000000000000000000..5b1cfd8ac1798e412479f971b76bfad982a5f1b4 --- /dev/null +++ b/FunctionSpace/FunctionSpaceNode.h @@ -0,0 +1,18 @@ +#ifndef _FUNCTIONSPACENODE_H_ +#define _FUNCTIONSPACENODE_H_ + +#include "FunctionSpaceScalar.h" + +class FunctionSpaceNode : public FunctionSpaceScalar{ + public: + FunctionSpaceNode(const GroupOfElement& goe, int order); + + virtual ~FunctionSpaceNode(void); + + virtual double + interpolate(const MElement& element, + const std::vector<double>& coef, + const fullVector<double>& xyz) const; +}; + +#endif diff --git a/FunctionSpace/FunctionSpaceScalar.cpp b/FunctionSpace/FunctionSpaceScalar.cpp new file mode 100644 index 0000000000000000000000000000000000000000..bda936c12925b8ec72e49d9ec5c878d2402b8f14 --- /dev/null +++ b/FunctionSpace/FunctionSpaceScalar.cpp @@ -0,0 +1,4 @@ +#include "FunctionSpaceScalar.h" + +FunctionSpaceScalar::~FunctionSpaceScalar(void){ +} diff --git a/FunctionSpace/FunctionSpaceScalar.h b/FunctionSpace/FunctionSpaceScalar.h new file mode 100644 index 0000000000000000000000000000000000000000..e6c942d3ee8c20d4067ee08bbc9428efec8ad6fb --- /dev/null +++ b/FunctionSpace/FunctionSpaceScalar.h @@ -0,0 +1,17 @@ +#ifndef _FUNCTIONSPACESCALAR_H_ +#define _FUNCTIONSPACESCALAR_H_ + +#include "fullMatrix.h" +#include "FunctionSpace.h" + +class FunctionSpaceScalar : public FunctionSpace{ + public: + virtual ~FunctionSpaceScalar(void); + + virtual double + interpolate(const MElement& element, + const std::vector<double>& coef, + const fullVector<double>& xyz) const = 0; +}; + +#endif diff --git a/FunctionSpace/FunctionSpaceVector.cpp b/FunctionSpace/FunctionSpaceVector.cpp new file mode 100644 index 0000000000000000000000000000000000000000..821e8505ce8dd95d94e9b51cb42f19143048dfa9 --- /dev/null +++ b/FunctionSpace/FunctionSpaceVector.cpp @@ -0,0 +1,4 @@ +#include "FunctionSpaceVector.h" + +FunctionSpaceVector::~FunctionSpaceVector(void){ +} diff --git a/FunctionSpace/FunctionSpaceVector.h b/FunctionSpace/FunctionSpaceVector.h new file mode 100644 index 0000000000000000000000000000000000000000..f711ea515f825cbab9d85526a97e9a9233715d2f --- /dev/null +++ b/FunctionSpace/FunctionSpaceVector.h @@ -0,0 +1,17 @@ +#ifndef _FUNCTIONSPACEVECTOR_H_ +#define _FUNCTIONSPACEVECTOR_H_ + +#include "fullMatrix.h" +#include "FunctionSpace.h" + +class FunctionSpaceVector : public FunctionSpace{ + public: + virtual ~FunctionSpaceVector(void); + + virtual fullVector<double> + interpolate(const MElement& element, + const std::vector<double>& coef, + const fullVector<double>& xyz) const = 0; +}; + +#endif