Commit fca0e1c4 authored by Joe Frazier's avatar Joe Frazier
Browse files

Now works with NekManager.


git-svn-id: https://gforge.sci.utah.edu/svn/nektar/trunk@245 305cdda6-5ce1-45b3-a98d-dfc68c8b3305
parent 17368e5e
......@@ -36,11 +36,9 @@
#ifndef FOUNDATIONS_H
#define FOUNDATIONS_H
#include <string>
namespace Nektar
{
namespace LibUtiltities
namespace LibUtilities
{
enum BasisType
{
......
......@@ -45,18 +45,26 @@ namespace Nektar
{
namespace LibUtilities
{
class PointsKey;
// Use for looking up the creator. The creator for number of points
// can generate for any number, so we want the same creator called
// for all number.
struct opLess
{
bool operator()(const PointsKey &lhs, const PointsKey &rhs);
};
class PointsKey
{
public:
PointsKey()
{
NEKERROR(ErrorUtil::efatal,"Default Constructor for PointsKey should not be called");
}
PointsKey(const int &pointsorder, const PointsType &pointstype,
const PointsIdentifier &pointsid): m_pointsorder(pointsorder),
m_pointstype(pointstype), m_pointsid(pointsid)
PointsKey(const unsigned int &pointsdim, const int &numpoints,
const PointsType &pointstype, const PointsIdentifier &pointsid = eWildcard):
m_pointsdim(pointsdim),
m_numpoints(numpoints),
m_pointstype(pointstype),
m_pointsid(pointsid)
{
}
......@@ -66,20 +74,27 @@ namespace Nektar
PointsKey(const PointsKey &key)
{
*this = key;
*this = key; // defer to assignment operator
}
PointsKey& operator=(const PointsKey &key)
{
m_pointsorder = key.m_pointsorder;
m_pointsdim = key.m_pointsdim;
m_numpoints = key.m_numpoints;
m_pointstype = key.m_pointstype;
m_pointsid = key.m_pointsid;
return *this;
}
inline int GetPointsOrder() const
inline unsigned int GetPointsDim() const
{
return m_pointsorder;
return m_pointsdim;
}
inline int GetNumPoints() const
{
return m_numpoints;
}
inline PointsType GetPointsType() const
......@@ -92,215 +107,163 @@ namespace Nektar
return m_pointsid;
}
bool operator==(const PointsKey &key)
{
return (m_pointsdim == key.m_pointsdim &&
m_numpoints == key.m_numpoints &&
m_pointstype == key.m_pointstype &&
m_pointsid == key.m_pointsid);
}
bool operator == (const PointsKey *y)
{
return (*this == *y);
}
bool operator != (const PointsKey& y)
{
return (!(*this == y));
}
bool operator != (const PointsKey *y)
{
return (!(*this == *y));
}
friend bool operator<(const PointsKey &lhs, const PointsKey &rhs);
friend bool opLess::operator()(const PointsKey &lhs, const PointsKey &rhs);
protected:
int m_pointsorder; //!< "Order" of the points (as appropriately defined for PointsType)
unsigned int m_pointsdim; //!< dimension of the points
int m_numpoints; //!< number of the points (as appropriately defined for PointsType)
PointsType m_pointstype; //!< Type of Points
PointsIdentifier m_pointsid; //!< Unique indentifier (when needed)
PointsIdentifier m_pointsid; //!< Unique indentifier (when needed)
private:
PointsKey()
{
NEKERROR(ErrorUtil::efatal,"Default Constructor for PointsKey should not be called");
}
};
bool operator<(const PointsKey &lhs, const PointsKey &rhs);
std::ostream& operator<<(std::ostream& os, const PointsKey& rhs);
template<typename DataType, unsigned int dim>
template<typename DataType>
class Points
{
public:
Points()
Points(const PointsKey &key): m_pkey(key)
{
NEKERROR(ErrorUtil::efatal,"Default Constructor for Points should not be called");
}
Points(const PointsKey &key): m_pkey(key)
{
CalculateNumPoints(); //populate m_numpoints
virtual ~Points()
// Allocate Memory
for(unsigned int i = 0; i < dim; ++i)
{
if(m_pkey.GetNumPoints())
{
m_points[i] = new DataType[m_numpoints];
//unsigned int dim = m_pkey.GetPointsDim();
//for(unsigned int i = 0; i < dim; ++i)
//{
// delete[] m_points[i];
//}
//delete[] m_points;
//delete[] m_weights;
}
}
m_weights = new DataType[m_numpoints];
m_derivmatrix.reset( new NekMatrix<DataType>(m_numpoints) );
CalculatePoints();
inline int GetPointsDim() const
{
return m_pkey.GetPointsDim();
}
CalculateWeights();
inline int GetNumPoints() const
{
return m_pkey.GetNumPoints();
}
CalculateDerivMatrix();
inline PointsType GetPointsType() const
{
return m_pkey.GetPointsType();
}
virtual ~Points()
{
if(m_numpoints)
{
for(unsigned int i = 0; i < dim; ++i)
{
delete[] m_points[i];
}
delete[] m_weights;
key.m_numpoints = 0;
}
}
inline int GetPointsOrder() const
{
return m_pkey.GetNumPoints();
}
inline int GetNumPoints() const
{
return m_numpoints;
}
inline PointsType GetPointsType() const
{
return m_pkey.GetPointsType();
}
inline PointsIdentifier GetPointsId() const
{
return m_pkey.GetPointsId();
}
inline double *GetZ() const
{
BOOST_STATIC_ASSERT(dim == 1);
return m_points[0];
}
inline double *GetW() const
{
return m_weights;
}
inline void GetZW(const double *&z, const double *&w) const
{
BOOST_STATIC_ASSERT(dim == 1);
z = m_points[0];
w = m_weights;
}
inline void GetPoints(const double *&x) const
{
BOOST_STATIC_ASSERT(dim == 1);
x = m_points[0];
}
inline void GetPoints(const double *&x, const double *&y) const
{
BOOST_STATIC_ASSERT(dim == 2);
x = m_points[0];
y = m_points[1];
}
inline void GetPoints(const double *&x, const double *&y, const double *&z) const
{
BOOST_STATIC_ASSERT(dim == 3);
x = m_points[0];
y = m_points[1];
z = m_points[2];
}
inline const boost::shared_ptr<NekMatrix<DataType> > GetD() const
{
return m_derivmatrix;
inline PointsIdentifier GetPointsId() const
{
return m_pkey.GetPointsId();
}
inline double *GetZ() const
{
BOOST_STATIC_ASSERT(dim == 1);
return m_points[0];
}
inline double *GetW() const
{
return m_weights;
}
inline void GetZW(const double *&z, const double *&w) const
{
BOOST_STATIC_ASSERT(dim == 1);
z = m_points[0];
w = m_weights;
}
inline void GetPoints(const double *&x) const
{
BOOST_STATIC_ASSERT(dim == 1);
x = m_points[0];
}
inline void GetPoints(const double *&x, const double *&y) const
{
BOOST_STATIC_ASSERT(dim == 2);
x = m_points[0];
y = m_points[1];
}
inline void GetPoints(const double *&x, const double *&y, const double *&z) const
{
BOOST_STATIC_ASSERT(dim == 3);
x = m_points[0];
y = m_points[1];
z = m_points[2];
}
inline const boost::shared_ptr<NekMatrix<DataType> > GetD() const
{
return m_derivmatrix;
}
protected:
PointsKey m_pkey;
int m_numpoints;
DataType *m_points[dim];
DataType **m_points;
DataType *m_weights;
boost::shared_ptr<NekMatrix<DataType> > m_derivmatrix;
private:
virtual void CalculateNumPoints() = 0;
virtual void CalculatePoints() = 0;
virtual void CalculateWeights() = 0;
virtual void CalculateDerivMatrix() = 0;
// This should never be called.
Points()
{
NEKERROR(ErrorUtil::efatal,"Default Constructor for Points should not be called");
}
};
bool operator == (const PointsKey& x, const PointsKey& y)
{
if( (x.m_pointsorder == y.m_pointsorder) &&
(x.m_pointstype == y.m_pointstype) &&
(x.m_pointsid == y.m_pointsid) )
{
return true;
}
return false;
}
bool operator == (const PointsKey* x, const PointsKey& y)
{
if( ((*x).m_pointsorder == y.m_pointsorder) &&
((*x).m_pointstype == y.m_pointstype) &&
((*x).m_pointsid == y.m_pointsid) )
{
return true;
}
return false;
}
bool operator == (const PointsKey& x, const PointsKey *y)
{
if( (x.m_pointsorder == (*y).m_pointsorder) &&
(x.m_pointstype == (*y).m_pointstype) &&
(x.m_pointsid == (*y).m_pointsid) )
{
return true;
}
return false;
}
bool operator != (const PointsKey& x, const PointsKey& y)
{
if( (x.m_pointsorder == y.m_pointsorder) &&
(x.m_pointstype == y.m_pointstype) &&
(x.m_pointsid == y.m_pointsid) )
{
return false;
}
return true;
}
bool operator != (const PointsKey* x, const PointsKey& y)
{
if( ((*x).m_pointsorder == y.m_pointsorder) &&
((*x).m_pointstype == y.m_pointstype) &&
((*x).m_pointsid == y.m_pointsid) )
{
return false;
}
return true;
}
bool operator != (const PoinstKey& x, const PointsKey *y)
{
if( (x.m_pointsorder == (*y).m_pointsorder) &&
(x.m_pointstype == (*y).m_pointstype) &&
(x.m_pointsid == (*y).m_pointsid) )
{
return false;
}
return true;
}
} // end of namespace
}; // end of namespace
} // end of namespace
#endif //POINTS_HPP
......@@ -34,10 +34,12 @@
///////////////////////////////////////////////////////////////////////////////
#include <iostream>
#include <StdRegions/StdPoints.h>
#include <LibUtilities/Foundations/Points.hpp>
#include <LibUtilities/Foundations/Foundations.hpp>
#include <LibUtilities/BasicUtils/ErrorUtil.hpp>
#include <LibUtilities/Polylib/Polylib.h>
#include <LibUtilities/Foundations/Points1d.h>
namespace Nektar
{
......@@ -55,78 +57,20 @@ namespace Nektar
{
}
GaussPolyPoints::GaussPolyPoints(const PointsKey &key,
const double alpha, const double beta):
Points<double>(key),
m_alpha(alpha),
m_beta(beta)
{
}
boost::shared_ptr< Points<double> > GaussPolyPoints::Create(const PointsKey &key)
{
boost::shared_ptr< Points<double> > returnval(new GaussPolyPoints(key, 1.0, 1.0));
return returnval;
}
GaussPolyPoints::GaussPolyPoints(const int npts, PointsType ptype,
const double alpha, const double beta):
PolyPoints(npts),
m_gausspointstype(ptype),
m_alpha(alpha),
m_beta(beta)
{
switch(m_gausspointstype)
{
case eGauss:
Polylib::zwgj(m_zeros,m_weights,m_pointsorder,m_alpha,
m_beta);
Polylib::Dgj(m_derivmatrix,m_zeros,m_pointsorder,m_alpha,
m_beta);
break;
case eLobatto:
Polylib::zwglj(m_zeros,m_weights,m_pointsorder,m_alpha,
m_beta);
Polylib::Dglj(m_derivmatrix,m_zeros,m_pointsorder,m_alpha,
m_beta);
break;
case eRadauM:
Polylib::zwgrjm(m_zeros,m_weights,m_pointsorder,m_alpha,
m_beta);
Polylib::Dgrjm(m_derivmatrix,m_zeros,m_pointsorder,m_alpha,
m_beta);
break;
case eRadauP:
Polylib::zwgrjp(m_zeros,m_weights,m_pointsorder,m_alpha,
m_beta);
Polylib::Dgrjp(m_derivmatrix,m_zeros,m_pointsorder,m_alpha,
m_beta);
break;
default:
ErrorUtil::Error(ErrorUtil::efatal,
"GaussPolyPoints::GaussPolyPoints",
"Unknown points type");
break;
}
}
FourierPoints::FourierPoints(const int npts):
Points(npts)
{
if(m_pointsorder%2)
{
ErrorUtil::Error(ErrorUtil::efatal,
"FourierPoints::FourierPoints",
"Fourier points need to be of even order");
}
// define points in the region [-1,1]
for(int i = 0; i < m_pointsorder; ++i)
{
m_zeros [i] = -1.0 + i*2.0/(double)m_pointsorder;
m_weights[i] = 2.0/(double)m_pointsorder;
}
ErrorUtil::Error(ErrorUtil::ewarning,
"FourierPoints::FourierPoints",
"Fourier points Derivative matrix "
"need defining for Fourier spacing");
}
} // end of namespace stdregion
} // end of namespace stdregion
......
......@@ -38,49 +38,40 @@
#include <math.h>
#include <LibUtilities/Foundations/Foundations.hpp>
#include <LibUtilities/ErrorUtil.hpp>
#include <LibUtilities/BasicUtils/ErrorUtil.hpp>
#include <LibUtilities/Foundations/Points.hpp>
namespace Nektar
{
namespace LibUtilities
{
class GaussPolyPoints: public Points<double,1>
class GaussPolyPoints: public Points<double>
{
public:
GaussPolyPoints()
{
}
GaussPolyPoints(const GaussPolyPoints &g) const
{
NEKERROR(efatal, "This copy constructor should not be called");
}
GaussPolyPoints::GaussPolyPoints(const PointsKey &key,
const double alpha, const double beta);
virtual ~GaussPolyPoints()
{
}
// Not sure what to return here for the create function.
// Since this function will be registered in a manager with
// other points types, it needs to return a base type.
// I am not sure if Points<> will work because of its type
// obtained from the template parameters. That would make
// it so the given manager can only hold, say, Points<double>.
static boost::shared_ptr< Points<double> > Create(const PointsKey &key);
protected:
double m_alpha;
double m_beta;
private:
void CalculatePoints();
void CalculateWeights();
void CalculateDerivMatrix();
};
class FourierPoints: public Points
{
public:
protected:
private:
};
} // end of namespace
} // end of namespace
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment