diff --git a/FunctionSpace/BasisScalar.h b/FunctionSpace/BasisScalar.h
index 0673f8a1d8f696f3b663698900d65deb1c941382..74c19e667011025236f498fd23326907ff542044 100644
--- a/FunctionSpace/BasisScalar.h
+++ b/FunctionSpace/BasisScalar.h
@@ -25,7 +25,7 @@ class BasisScalar: public Basis{
//! @return Returns the set of @em Polynomial%s
//! defining this (scalar) Basis
- const std::vector<Polynomial>& getBasis(void) const;
+ const std::vector<Polynomial>& getFunctions(void) const;
protected:
//! Instantiate a new BasisScalar
@@ -37,7 +37,10 @@ class BasisScalar: public Basis{
// Inline Functions //
//////////////////////
-inline const std::vector<Polynomial>& BasisScalar::getBasis(void) const{
+inline
+const std::vector<Polynomial>& BasisScalar::
+getFunctions(void) const{
+
return *basis;
}
diff --git a/FunctionSpace/BasisVector.h b/FunctionSpace/BasisVector.h
index 1f60d4ef697a40ebe475756b895d5b180b87b051..77329e2f9901a3ff6c64bd0e35fdfbc581dc10c5 100644
--- a/FunctionSpace/BasisVector.h
+++ b/FunctionSpace/BasisVector.h
@@ -25,7 +25,7 @@ class BasisVector: public Basis{
//! @return Returns the set of @em Polynomial%s
//! defining this (vectorial) Basis
- const std::vector<std::vector<Polynomial> >& getBasis(void) const;
+ const std::vector<std::vector<Polynomial> >& getFunctions(void) const;
protected:
//! Instantiate a new BasisVector
@@ -37,7 +37,10 @@ class BasisVector: public Basis{
// Inline Functions //
//////////////////////
-inline const std::vector<std::vector<Polynomial> >& BasisVector::getBasis(void) const{
+inline
+const std::vector<std::vector<Polynomial> >& BasisVector::
+getFunctions(void) const{
+
return *basis;
}
diff --git a/FunctionSpace/FunctionSpace.h b/FunctionSpace/FunctionSpace.h
index 770d3730d25c32d48af0a3e9f7e2afc1d70f54dd..4b4ce0e00580e26ebd9b92e3d7495c50eefc4acf 100644
--- a/FunctionSpace/FunctionSpace.h
+++ b/FunctionSpace/FunctionSpace.h
@@ -40,7 +40,6 @@ class FunctionSpace{
virtual ~FunctionSpace(void);
const GroupOfElement& getSupport(void) const;
- const Basis& getBasis(const MElement& element) const;
int getType(void) const;
unsigned int getNFunctionPerVertex(const MElement& element) const;
@@ -71,10 +70,6 @@ inline const GroupOfElement& FunctionSpace::getSupport(void) const{
return *goe;
}
-inline const Basis& FunctionSpace::getBasis(const MElement& element) const{
- return *basis;
-}
-
inline int FunctionSpace::getType(void) const{
return type;
}
diff --git a/FunctionSpace/FunctionSpaceEdge.cpp b/FunctionSpace/FunctionSpaceEdge.cpp
index 698822f77297eb8596f8a473749e2677dae5a48c..ee82ae870970e066a491249be4d1075c162278d6 100644
--- a/FunctionSpace/FunctionSpaceEdge.cpp
+++ b/FunctionSpace/FunctionSpaceEdge.cpp
@@ -15,39 +15,39 @@ FunctionSpaceEdge::FunctionSpaceEdge(const GroupOfElement& goe,
FunctionSpaceEdge::~FunctionSpaceEdge(void){
}
-#include <cstdio>
-
fullVector<double> FunctionSpaceEdge::
interpolate(const MElement& element,
const std::vector<double>& coef,
const fullVector<double>& xyz) const{
+ // Const Cast For MElement //
MElement& eelement =
const_cast<MElement&>(element);
- const vector<vector<Polynomial> >& fun =
- static_cast<const BasisVector&>(*basis).getBasis();
+ // Get Basis Functions //
+ const vector<vector<Polynomial> >& fun = getBasis(element).getFunctions();
+ const unsigned int nFun = fun.size();
- const unsigned int nFun = fun.size();
-
- fullVector<double> val(3);
-
+ // Get Jacobian //
fullMatrix<double> invJac(3, 3);
eelement.getJacobian(xyz(0), xyz(1), xyz(2), invJac);
invJac.invertInPlace();
+ // Get First Vertex of Element //
fullVector<double> origin(3);
origin(0) = eelement.getVertex(0)->x();
origin(1) = eelement.getVertex(0)->y();
origin(2) = eelement.getVertex(0)->z();
+ // Map XYZ to Reference Plane (UVW) //
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;
-
-
+
+ // Interpolate (in Reference Place) //
+ fullVector<double> val(3);
+ val(0) = 0;
+ val(1) = 0;
+ 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)),
@@ -57,5 +57,6 @@ interpolate(const MElement& element,
val.axpy(vi, 1);
}
+ // Return Interpolated Value //
return val;
}
diff --git a/FunctionSpace/FunctionSpaceNode.cpp b/FunctionSpace/FunctionSpaceNode.cpp
index 3f0dd595cd9beacf14c2df0a0a5c1f7f5b838482..aadb6aabdee8b3fbba215ead11f9b3b746028734 100644
--- a/FunctionSpace/FunctionSpaceNode.cpp
+++ b/FunctionSpace/FunctionSpaceNode.cpp
@@ -1,5 +1,11 @@
+#include <vector>
+
+#include "BasisScalar.h"
+#include "Mapper.h"
#include "FunctionSpaceNode.h"
+using namespace std;
+
FunctionSpaceNode::FunctionSpaceNode(const GroupOfElement& goe,
int order){
// Build 0Form Basis //
@@ -14,5 +20,36 @@ interpolate(const MElement& element,
const std::vector<double>& coef,
const fullVector<double>& xyz) const{
- return 42;
+ // Const Cast For MElement //
+ MElement& eelement =
+ const_cast<MElement&>(element);
+
+ // Get Basis Functions //
+ const vector<Polynomial>& fun = getBasis(element).getFunctions();
+ const unsigned int nFun = fun.size();
+
+ // Get Jacobian //
+ fullMatrix<double> invJac(3, 3);
+ eelement.getJacobian(xyz(0), xyz(1), xyz(2), invJac);
+
+ invJac.invertInPlace();
+
+ // Get First Vertex of Element (Physical Space) //
+ fullVector<double> origin(3);
+ origin(0) = eelement.getVertex(0)->x();
+ origin(1) = eelement.getVertex(0)->y();
+ origin(2) = eelement.getVertex(0)->z();
+
+ // Map XYZ to Reference Plane (UVW) //
+ fullVector<double> uvw = Mapper::invMap(xyz, origin, invJac);
+
+ // Interpolate (in Reference Place) //
+ double val = 0;
+
+ for(unsigned int i = 0; i < nFun; i++)
+ val +=
+ fun[i].at(uvw(0), uvw(1), uvw(2)) * coef[i];
+
+ // Return Interpolated Value //
+ return val;
}
diff --git a/FunctionSpace/FunctionSpaceScalar.h b/FunctionSpace/FunctionSpaceScalar.h
index e6c942d3ee8c20d4067ee08bbc9428efec8ad6fb..7587672ed9eac6cf6e530f560f0d4a8d4c88130a 100644
--- a/FunctionSpace/FunctionSpaceScalar.h
+++ b/FunctionSpace/FunctionSpaceScalar.h
@@ -2,6 +2,7 @@
#define _FUNCTIONSPACESCALAR_H_
#include "fullMatrix.h"
+#include "BasisScalar.h"
#include "FunctionSpace.h"
class FunctionSpaceScalar : public FunctionSpace{
@@ -12,6 +13,18 @@ class FunctionSpaceScalar : public FunctionSpace{
interpolate(const MElement& element,
const std::vector<double>& coef,
const fullVector<double>& xyz) const = 0;
+
+ const BasisScalar& getBasis(const MElement& element) const;
};
+
+//////////////////////
+// Inline Functions //
+//////////////////////
+
+inline const BasisScalar& FunctionSpaceScalar::
+getBasis(const MElement& element) const{
+ return static_cast<const BasisScalar&>(*basis);
+}
+
#endif
diff --git a/FunctionSpace/FunctionSpaceVector.h b/FunctionSpace/FunctionSpaceVector.h
index f711ea515f825cbab9d85526a97e9a9233715d2f..b1bb914eaa38e5e16da4b0c5efbcf22883a69ecd 100644
--- a/FunctionSpace/FunctionSpaceVector.h
+++ b/FunctionSpace/FunctionSpaceVector.h
@@ -2,6 +2,7 @@
#define _FUNCTIONSPACEVECTOR_H_
#include "fullMatrix.h"
+#include "BasisVector.h"
#include "FunctionSpace.h"
class FunctionSpaceVector : public FunctionSpace{
@@ -12,6 +13,18 @@ class FunctionSpaceVector : public FunctionSpace{
interpolate(const MElement& element,
const std::vector<double>& coef,
const fullVector<double>& xyz) const = 0;
+
+ const BasisVector& getBasis(const MElement& element) const;
};
+
+//////////////////////
+// Inline Functions //
+//////////////////////
+
+inline const BasisVector& FunctionSpaceVector::
+getBasis(const MElement& element) const{
+ return static_cast<const BasisVector&>(*basis);
+}
+
#endif
diff --git a/FunctionSpace/Polynomial.cpp b/FunctionSpace/Polynomial.cpp
index 7d2bf67b2c811a8c972f3eaad69f30acaeda191e..1edda763e19697aaccf4cebdedc9331709790d4a 100644
--- a/FunctionSpace/Polynomial.cpp
+++ b/FunctionSpace/Polynomial.cpp
@@ -1,6 +1,7 @@
#include <cmath>
#include <sstream>
#include <stack>
+
#include "Polynomial.h"
using namespace std;