element.hh

Go to the documentation of this file.

 
#ifndef Element1D_h
#define Element1D_h 
 
#include "basics/debug.hh"
#include "basics/output.hh"
#include "basics/typedefs.hh"
#include "basics/vectorsMatrices.hh"
#include "geometry/connector.hh"
#include "geometry/cell1D.hh"
#include "geometry/integral.hh"
#include "space/element.hh"
#include "space/tmatrix.hh"
#include "integration/quadRule.hh"
#include "integration/karniadakis.hh"
#include "integration/legendre.hh"
#include "lineGraphics.hh"
 
// debugging Element
#define ElementConstr_D        0
#define ElementDestr_D         0
 
namespace hp1D {
using concepts::Real;
 
// ********************************************************* IntegrableElm **
 
class IntegrableElm: public concepts::IntegrationCell {
public:
  IntegrableElm(const concepts::EdgeNd& cell);
 
  const concepts::Real3d chi(const Real x) const {
    return cell_.elemMap(x);
  }
 
  inline Real jacobianDeterminant(const Real x) const {
    // factor comes due to integration over [-1,1]^2
    return cell_.jacobianDeterminant(x) * 0.5;
  }
 
  const concepts::QuadratureRule1d* integration() const {
    return int_.get();
  }
 
  static concepts::QuadRuleFactory& rule() {
    return rule_;
  }
 
  virtual bool quadraturePoint(uint i, intPoint& p, intFormType form = ZERO,
      bool localCoord = false) const;
protected:
  const concepts::EdgeNd& cell_;
  std::unique_ptr<concepts::QuadratureRule1d> int_;
  // Factory for the integration rule
  static concepts::QuadRuleFactory rule_;
};
 
// ********************************************************* BaseElement **
 
template<class F>
class BaseElement: public concepts::ElementWithCell<F>, public IntegrableElm {
public:
  typedef F FieldT;
 
  BaseElement(const concepts::EdgeNd& cell, ushort p) :
      IntegrableElm(cell), p_(p) {
  }
  ;
 
  virtual ~BaseElement();
 
  virtual const concepts::Edge& support() const {
    return cell_.connector();
  }
 
  virtual concepts::Real3d vertex(uint i) const {
    return cell_.vertex(i);
  }
 
  virtual const concepts::EdgeNd& cell() const {
    return cell_;
  }
 
  ushort p() const {
    return p_;
  }
 
  virtual const concepts::ElementGraphics<Real>* graphics() const;
 
  virtual void recomputeShapefunctions() = 0;
 
  virtual const concepts::ShapeFunction1D<Real>* shpfct() const = 0;
 
protected:
  virtual std::ostream& info(std::ostream& os) const;
 
  // concepts::TMatrix<Real> T_;
 
  ushort p_;
 
  static std::unique_ptr<LineGraphics> graphics_;
 
};
 
 
// *************************************************************** KarniadakisMixin ***
 
template<class F>
class KarniadakisMixin: public BaseElement<F> {
public:
  KarniadakisMixin(const concepts::EdgeNd& cell, ushort p) :
      BaseElement<F>(cell, p), shpfct_(nullptr), shpfctD_(nullptr), shpfctDD_(
          nullptr), secondDerivativeUsed_(false) {
    recomputeShapefunctions();
  }
 
  void recomputeShapefunctions() override;
 
  const concepts::Karniadakis<1, 0>* shpfct() const override {
    return shpfct_.get();
  }
  const concepts::Karniadakis<1, 1>* shpfctD() const {
    return shpfctD_.get();
  }
 
  const concepts::Karniadakis<1, 2>* shpfctDD() const {
    if (!secondDerivativeUsed_) {
      shpfctDD_.reset(
          new concepts::Karniadakis<1, 2>(this->p_, this->int_->abscissas(),
              this->int_->n()));
      secondDerivativeUsed_ = true;
    }
    return shpfctDD_.get();
  }
 
private:
  void recomputeSecondDerivativeOfShapefunctions_();
 
  std::unique_ptr<concepts::Karniadakis<1, 0> > shpfct_;
  std::unique_ptr<concepts::Karniadakis<1, 1> > shpfctD_;
 
  mutable std::unique_ptr<concepts::Karniadakis<1, 2> > shpfctDD_;
 
  mutable bool secondDerivativeUsed_;
 
};
 
 
 
// ******************************************************** LegendreMixin ***
 
template<class F>
class LegendreMixin: public BaseElement<F> {
public:
  LegendreMixin(const concepts::EdgeNd& cell, ushort p) : BaseElement<F>(cell, p) {
    recomputeShapefunctions();
  }
 
  void recomputeShapefunctions() override;
 
  const concepts::Legendre* shpfct() const override {
    return shpfct_.get();
  }
 
private:
  std::unique_ptr<concepts::Legendre> shpfct_;
 
};
 
 
// ******************************************************* GenericElement **
 
template<class BaseT>
class GenericElement: public BaseT {
public:
  GenericElement(const concepts::EdgeNd& cell, uint p,
      concepts::TColumn<typename BaseT::FieldT>* T0,
      concepts::TColumn<typename BaseT::FieldT>* T1)
    : BaseT(cell,p), T_(T1)
  {
    if (T0)
      T_.append(T0);
    DEBUGL(ElementConstr_D, cell << ", p = " << p << ", T = " << T_);
  }
 
  ~GenericElement() override
  {
    DEBUGL(ElementDestr_D, "done.");
  }
 
  const concepts::TMatrix<typename BaseT::FieldT>& T() const override {
    return T_;
  }
 
  void appendT(concepts::TColumn<typename BaseT::FieldT>* T) {
    T_.append(T);
  }
 
protected:
  virtual std::ostream& info(std::ostream& os) const override
  {
    return os << typeOf(*this);    
  }
 
  concepts::TMatrix<typename BaseT::FieldT> T_;
 
};
 
template<class F>
using Element = GenericElement<KarniadakisMixin<F> >;
 
 
template<class F>
using LegendreElement = GenericElement<hp1D::LegendreMixin<F> >;
 
}// namespace hp1D
 
namespace concepts {
 
// ****************************************************************** Scan **
 
 
template<class F>
class Scan<hp1D::BaseElement<F> > : public Scan<ElementWithCell<F> > {
public:
  hp1D::BaseElement<F>& operator++(int) = 0;
};
 
} // namespace concepts
 
#endif // Element1D_h