cell2D.hh

Go to the documentation of this file.

 
#ifndef cell2D_hh
#define cell2D_hh
 
#include <iostream>
 
#include "basics/outputOperator.hh"
#include "basics/typedefs.hh"
#include "basics/level.hh"
#include "topology.hh"
#include "cell.hh"
#include "cell1D.hh"
#include "elementMaps.hh"
 
namespace concepts {
 
 
  // ************************************************************ Triangle2d **
 
  class Triangle2d : public Cell2 {
    public:
 
      /* Type for indices in the hierarchical definition of the cells
 
      Use a type with a larger range depending on the number of cell
      refinements.
      */
      typedef uint index_type;
 
    // Maximal number of levels in the cell refinement
 
      static uint MAX_LEVEL;
 
      struct Index {
        uchar  l_;
 
        index_type i_;
 
        index_type j_;
 
        Index();
        Index(uchar l, index_type i, index_type j);
        std::ostream& info(std::ostream& os) const;
      };
 
      friend std::ostream& operator<<(std::ostream& os,
          const Triangle2d::Index& i);
 
      Triangle2d(Triangle& cntr, const MappingTriangle2d& map);
      virtual ~Triangle2d();
 
      virtual Triangle2d* child(uint i);
      virtual const Triangle2d* child(uint i) const;
 
      inline Triangle& connector() const { return cntr_; }
 
      ushort level() const { return idx_.l_; }
 
      Real2d chi(Real xi, Real eta) const;
 
 
      Real3d elemMap(const Real2d& coord_local) const
      {
        return chi(coord_local[0], coord_local[1]);
      }
 
      Real jacobianDeterminant(const Real x, const Real y) const;
 
      MapReal2d jacobian(Real x, Real y) const;
 
      Real area() const;
 
      MapReal2d jacobianInverse(const Real x, const Real y) const;
 
      Real3d vertex(uint i) const;
      Real3d center() const {return chi(0.66666666666, 0.5);}
      inline const MappingTriangle2d* map() const {return map_;}
 
      inline Triangle2d* clone(Triangle& cntr, MappingTriangle2d* map) const {
        return new Triangle2d(cntr, map, idx_);}
 
      static Real2d duffyInv(const Real2d& x);
    protected:
      virtual std::ostream& info(std::ostream& os) const;
    private:
      Index idx_;
 
      Triangle& cntr_;
 
      Triangle2d* lnk_;
 
      Triangle2d* chld_;
 
      MappingTriangle2d* map_;
 
      Triangle2d(Triangle& cntr, MappingTriangle2d* map, const Index& idx);
  };
 
  std::ostream& operator<<(std::ostream& os, const Triangle2d::Index& i);
 
  
  // **************************************************************** QuadNd **
  class QuadNd : public Cell2 {
  public:
    /* Type for indices in the hierarchical definition of the cells
       
    Use a type with a larger range depending on the number of cell
    refinements.
    */
    typedef uint index_type;
 
  // Maximal number of levels in the cell refinement
 
    static uint MAX_LEVEL;
 
    struct Index {
      Level<2> level_;
 
      index_type i_;
 
      index_type j_;
 
      Index() : level_(0), i_(0), j_(0) {}
      Index(uchar k, uchar l, index_type i, index_type j) :
        level_(k), i_(i), j_(j) { level_.l_[1] = l; }
      std::ostream& info(std::ostream& os) const;
    };
    
    QuadNd( Quad &cntr ) : cntr_(cntr) { }
    
    QuadNd( Quad &cntr, const Index &idx ) : cntr_(cntr), idx_(idx) { }
 
    friend std::ostream& operator<<(std::ostream& os,
                                    const QuadNd::Index& i);
    
    virtual uint dim() const = 0;
    
    inline Quad& connector() const {
      return cntr_;
    }
    
    virtual Real gramDeterminantRoot(const Real xi, const Real eta) const = 0;
    
    virtual Real lineElement(const Real xi, const uint edge) const = 0;
    
    Level<2> level() const { return idx_.level_; }
    
   protected:
    Quad &cntr_;
    Index idx_;
 
  };
  
  // ***************************************************** Quad2dSubdivision **
 
  class Quad2d;
 
  class Quad2dSubdivision : public OutputOperator {
    public:
      virtual void createChildren(Quad2d& q) const = 0;
      virtual void removeChildren(Quad2d& q) const;
      virtual const QuadSubdivision* topologicalStrategy() const = 0;
  };
 
  // ********************************************************* Quad2dSubdiv4 **
 
  class Quad2dSubdiv4 : public Quad2dSubdivision {
    public:
      virtual void createChildren(Quad2d& q) const;
      virtual const QuadSubdiv4* topologicalStrategy() const;
      static const Quad2dSubdiv4* instance();
      virtual ~Quad2dSubdiv4();
    protected:
      virtual std::ostream& info(std::ostream& os) const;
    private:
      static std::unique_ptr<Quad2dSubdiv4> instance_;
  };
 
  // ******************************************************** Quad2dSubdiv2H **
 
  class Quad2dSubdiv2H : public Quad2dSubdivision {
    public:
      virtual void createChildren(Quad2d& q) const;
      virtual const QuadSubdiv2H* topologicalStrategy() const;
      static const Quad2dSubdiv2H* instance();
      virtual ~Quad2dSubdiv2H();
    protected:
      virtual std::ostream& info(std::ostream& os) const;
    private:
      static std::unique_ptr<Quad2dSubdiv2H> instance_;
  };
 
  // ******************************************************** Quad2dSubdiv2V **
 
  class Quad2dSubdiv2V : public Quad2dSubdivision {
    public:
      virtual void createChildren(Quad2d& q) const;
      virtual const QuadSubdiv2V* topologicalStrategy() const;
      static const Quad2dSubdiv2V* instance();
      virtual ~Quad2dSubdiv2V();
    protected:
      virtual std::ostream& info(std::ostream& os) const;
    private:
      static std::unique_ptr<Quad2dSubdiv2V> instance_;
  };
 
  // **************************************************************** Quad2d **
 
  class Quad2d : public QuadNd {
    friend class Quad2dSubdivision;
    friend class Quad2dSubdiv4;
    friend class Quad2dSubdiv2H;
    friend class Quad2dSubdiv2V;
  public:
    Quad2d(Quad& cntr, const MappingQuad2d& map);
    virtual ~Quad2d();
 
  private:
    Quad2d(const Quad2d& other);
    Quad2d& operator=(const Quad2d& other);
 
  public:
 
    virtual Quad2d* child(uint i);
    virtual const Quad2d* child(uint i) const;
    inline bool hasChildren() const { return chld_ != 0; }
 
    const Quad2d* 
      child(uint i, Real xi, Real eta, Real& xiC, Real& etaC) const {
        const Quad2d* chld = child(i);
        if (chld == 0) return 0;
 
        Level<2> l = chld->level();
        if (idx_.level_.l_[0] == l.l_[0]) xiC = xi;
        else {
          xiC = (xi + idx_.i_) * 2.0 - chld->idx_.i_;
          if (xiC < 0.0 || xiC > 1.0) return 0;
        }
        if (idx_.level_.l_[1] == l.l_[1]) etaC = eta;
        else {
          etaC = (eta + idx_.j_) * 2.0 - chld->idx_.j_;
          if (etaC < 0.0 || etaC > 1.0) return 0;
        }
        return chld;
      }
 
    Real2d chi(Real xi, Real eta) const;
 
    Real3d elemMap(const Real2d& coord_local) const
    {
      return chi(coord_local[0], coord_local[1]);
    }
 
    MapReal2d jacobian(const Real xi, const Real eta) const;
 
    MapReal2d jacobian(const concepts::Real2d& p) const;
 
    MapReal2d jacobianInverse(const Real xi, const Real eta) const {
      return jacobian(xi, eta).inverse();
    }
 
    Real jacobianDeterminant(const Real xi, const Real eta) const {
      return jacobian(xi, eta).determinant();
    }
    
    virtual Real gramDeterminantRoot(const Real xi, const Real eta) const
    {
      return jacobianDeterminant(xi, eta);
    }
 
    MapReal2d inverseLaplace(const Real xi, const Real eta) const;
 
    virtual Real lineElement(const Real xi, const uint edge) const;
 
    Real3d vertex(uint i) const;
    Real3d center() const {return chi(0.5, 0.5);}
 
    inline const MappingQuad2d* map() const {return map_;}
    
    MappingEdge2d* edgeMap(uint edge) const;
    
    inline Quad2d* clone(Quad& cntr, MappingQuad2d* map) const {
      return new Quad2d(cntr, map, idx_);}
 
    virtual void setStrategy(const Quad2dSubdivision* strategy = 0)
      throw(StrategyChange);
 
    const Quad2dSubdivision* getStrategy() const { return subdivStrategy_; }
    
    virtual uint dim() const { return 2; }
  protected:
    virtual std::ostream& info(std::ostream& os) const;
 
    MappingQuad2d* map_;
 
    Quad2d* createChild_(Quad& cntr, const Index& idx, bool flag) const;
 
    Quad2d(Quad& cntr, MappingQuad2d* map, const Index& idx);
    private:
    Quad2d* lnk_;
 
    Quad2d* chld_;
 
    const Quad2dSubdivision* subdivStrategy_;
 
    // Edges as cells
    mutable concepts::Array<Edge2d*> edge_;
 
    virtual Quad2d* createChild_(Quad& cntr, const Index& idx) const;
  };
 
  std::ostream& operator<<(std::ostream& os, const Quad2d::Index& i);
 
  // ******************************************************** InfiniteQuad2d **
 
  class InfiniteQuad2d : public Cell2 {
  public:
    InfiniteQuad2d(InfiniteQuad& cntr, const Real2d& vtx0, const Real2d& vtx1);
 
    virtual Real2d chi(Real xi, Real eta) const = 0;
 
    Real3d vertex(uint i) const;
 
    virtual MapReal2d jacobian(const Real xi, const Real eta) const = 0;
 
    inline InfiniteQuad& connector() const { return cntr_; }
  protected:
    virtual std::ostream& info(std::ostream& os) const;
  private:
    Real2d vtx_[2];
    InfiniteQuad& cntr_;
  };
 
  // ******************************************************** InfiniteRect2d **
 
  class InfiniteRect2d : public InfiniteQuad2d {
    public:
    InfiniteRect2d(InfiniteQuad& cntr, const Real2d& vtx0, const Real2d& vtx1);
    virtual ~InfiniteRect2d();
 
    virtual InfiniteRect2d* child(uint i);
    virtual const InfiniteRect2d* child(uint i) const;
    inline bool hasChildren() const { return chld_ != 0; }
 
    virtual Real2d chi(Real xi, Real eta) const;
 
    Real3d elemMap(const Real2d& coord_local) const
    {
      return chi(coord_local[0], coord_local[1]);
    }
 
    virtual MapReal2d jacobian(const Real xi, const Real eta) const;
 
    Real3d center() const;
 
    inline InfiniteRect2d* clone(InfiniteQuad& cntr,
                                 const Real2d& vtx0, const Real2d& vtx1) const
    {
      return new InfiniteRect2d(cntr, vtx0, vtx1);
    }
 
  protected:
    virtual std::ostream& info(std::ostream& os) const;
  private:
    UnitNd<2> n_;
 
    InfiniteRect2d* lnk_;
 
    InfiniteRect2d* chld_;
  };
 
  // ************************************************************ Triangle3d **
 
  class Triangle3d : public Cell2 {
    public:
 
      /* Type for indices in the hierarchical definition of the cells
 
      Use a type with a larger range depending on the number of cell
      refinements.
      */
      typedef uint index_type;
 
    // Maximal number of levels in the cell refinement
 
      static uint MAX_LEVEL;
 
      struct Index {
        uchar  l_;
 
        index_type i_;
 
        index_type j_;
 
        Index();
        Index(uchar l, index_type i,index_type j);
        std::ostream& info(std::ostream& os) const;
      };
 
      friend std::ostream& operator<<(std::ostream& os,
          const Triangle3d::Index& i);
 
      Triangle3d(Triangle& cntr, const MapTriangle3d& map);
      virtual ~Triangle3d();
 
      virtual Triangle3d* child(uint i);
      virtual const Triangle3d* child(uint i) const;
 
      Triangle& connector() const { return cntr_; }
 
      Real3d chi(Real xi, Real eta) const;
 
      Real3d elemMap(const Real2d& coord_local) const
      {
        return chi(coord_local[0], coord_local[1]);
      }
 
      Real3d vertex(uint i) const;
      Real3d center() const {return chi(0.66666666666, 0.5);}
 
      inline const MapTriangle3d* map() const {return map_;}
 
      inline Triangle3d* clone(Triangle& cntr, MapTriangle3d* map) const {
        return new Triangle3d(cntr, map, idx_);}
    protected:
      virtual std::ostream& info(std::ostream& os) const;
    private:
      Index idx_;
 
      Triangle& cntr_;
 
      Triangle3d* lnk_;
 
      Triangle3d* chld_;
 
      MapTriangle3d* map_;
 
      Triangle3d(Triangle& cntr, MapTriangle3d* map, const Index& key);
  };
 
  std::ostream& operator<<(std::ostream& os, const Triangle3d::Index& i);
 
  // **************************************************************** Quad3d **
 
  // forward declaration
  class MappingQuad3d;
 
  class Quad3d : public QuadNd {
//     friend class Quad3dSubdivision;
//     friend class Quad3dSubdiv4;
//     friend class Quad3dSubdiv2H;
//     friend class Quad3dSubdiv2V;
 
  public:
    Quad3d(Quad& cntr, const MappingQuad3d& map, const Z2 dir = 0);
    virtual ~Quad3d();
 
  private:
    Quad3d(Quad& cntr, MappingQuad3d* map, const Index& idx, const Z2 dir = 0);
    Quad3d(const Quad3d& other);
    Quad3d& operator=(const Quad3d& other);
 
  public:
 
    virtual Quad3d* child(uint i);
    virtual const Quad3d* child(uint i) const;
    inline bool hasChildren() const { return chld_ != 0; }
 
    Real3d chi(Real xi, Real eta) const;
 
    Real3d elemMap(const Real2d& coord_local) const
    {
      return chi(coord_local[0], coord_local[1]);
    }
 
    Real3d jacobianColumn(uint i, const Real xi, const Real eta) const;
 
    Real faceElement(const Real xi, const Real eta) const;
    
    virtual Real gramDeterminantRoot(const Real xi, const Real eta) const
    {
      return faceElement(xi, eta);
    }
    
    virtual Real lineElement(const Real xi, const uint edge) const;
    
    Real3d n0(Real xi, Real eta) const;
 
    Real3d vertex(uint i) const;
    Real3d center() const {return chi(0.5, 0.5);}
 
    inline const MappingQuad3d* map() const {return map_;}
    
    inline Quad3d* clone(Quad& cntr, MappingQuad3d* map) const {
      return new Quad3d(cntr, map, idx_);}
 
//     virtual void setStrategy(const Quad3dSubdivision* strategy = 0)
//       throw(StrategyChange);
 
//     const Quad3dSubdivision* getStrategy() const { return subdivStrategy_; }
    
    virtual uint dim() const { return 3; }
    
    Z2 dir() const { return dir_; }
    
  protected:
    virtual std::ostream& info(std::ostream& os) const;
 
    MappingQuad3d* map_;
    
    Quad3d* createChild_(Quad& cntr, const Index& idx, bool flag) const;
 
  private:
    Quad3d* lnk_;
 
    Quad3d* chld_;
    
    Z2 dir_;
 
//     const Quad3dSubdivision* subdivStrategy_;
 
    virtual Quad3d* createChild_(Quad& cntr, const Index& idx) const;
  };
  
 
} // namespace concepts
 
#endif // cell2D_hh