00001 #ifndef QUADBASISSPECIFIER_H 00002 #define QUADBASISSPECIFIER_H 00003 00004 #include "SundanceDefs.h" 00005 00006 #include "QuadratureFamily.h" 00007 #include "BasisFamily.h" 00008 00009 namespace Sundance 00010 { 00011 00012 using namespace TSF; 00013 using std::string; 00014 00015 using std::ostream; 00016 00017 class WorkSet; 00018 00019 /** \ingroup LowLevelFE QuadBasisSpecifier specifies a combination 00020 * of {quadrature points, basis evaluator, and differentiation 00021 * order}. With such a combination, one can compute and cache a set 00022 * of basis function values (or derivative values) on quadrature 00023 * points. These function values will be stored in a hashtable keyed 00024 * by the QuadBasisSpecifier that created them. 00025 * 00026 * The basis values are stored in a DenseSerialVector indexed as follows: 00027 * values[differentiation direction][cell number][quad pt number] 00028 * [basis function number]. 00029 **/ 00030 00031 class QuadBasisSpecifier 00032 { 00033 public: 00034 /** empty ctor for containers */ 00035 QuadBasisSpecifier() : quad_(), basis_(), derivOrder_(0){} 00036 /** construct given a quad point set, a basis evaluator, and 00037 * the order of differentiation used when evaluating the basis */ 00038 QuadBasisSpecifier(const QuadratureFamily& quad, 00039 const BasisFamily& basis, 00040 const CellTopologyCode& cellType, 00041 int derivOrder); 00042 00043 /** compute the basis function values at all quadrature points, and store 00044 * in a vector */ 00045 TSFSmartPtr<DenseSerialVector> computeValues(const WorkSet& workSet) const ; 00046 /** compute the basis function values at all quadrature points, and store 00047 * in a vector */ 00048 TSFSmartPtr<DenseSerialVector> computeEmbeddedValues(const WorkSet& workSet) const ; 00049 /** test equality with another specifier. Two objects are equal if 00050 * they refer to the same quad points, basis, and deriv order */ 00051 bool operator==(const QuadBasisSpecifier& other) const ; 00052 00053 /** hash code for hashtable lookup */ 00054 int hashCode() const {return hashCode_;} 00055 00056 /** write as a string */ 00057 string toString() const ; 00058 protected: 00059 QuadratureFamily quad_; 00060 BasisFamily basis_; 00061 CellTopologyCode cellType_; 00062 int derivOrder_; 00063 int hashCode_; 00064 }; 00065 00066 } 00067 00068 namespace TSF 00069 { 00070 inline string toString(const Sundance::QuadBasisSpecifier& qbs) 00071 { 00072 return qbs.toString(); 00073 } 00074 00075 inline int hashCode(const Sundance::QuadBasisSpecifier& qbs) 00076 { 00077 return qbs.hashCode(); 00078 } 00079 } 00080 #endif