geometry/cell2D.hh

Go to the documentation of this file.

/* 2D Cells, ie. triangles and quads
 * Geometrical information (coordinates), brings element maps and
 * topology together
 */

#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:
      struct Index {
        uchar  l_;

        ushort i_;

        ushort j_;

        Index();
        Index(uchar l, ushort i, ushort 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);

  // ***************************************************** 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::auto_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::auto_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::auto_ptr<Quad2dSubdiv2V> instance_;
  };

  // **************************************************************** Quad2d **

  class Quad2d : public Cell2 {
    friend class Quad2dSubdivision;
    friend class Quad2dSubdiv4;
    friend class Quad2dSubdiv2H;
    friend class Quad2dSubdiv2V;
  public:
    struct Index {
      Level<2> level_;

      ushort i_;

      ushort j_;

      Index() : level_(0), i_(0), j_(0) {}
      Index(uchar k, uchar l, ushort i, ushort j) :
        level_(k), i_(i), j_(j) { level_.l_[1] = l; }
      std::ostream& info(std::ostream& os) const;
    };

    friend std::ostream& operator<<(std::ostream& os,
        const Quad2d::Index& i);

    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;
      }

    inline Quad& connector() const;

    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 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();
    }

    Real lineElement(const Real xi, const uint edge) const;

    Level<2> level() const { return idx_.level_; }

    Real3d vertex(uint i) const;
    Real3d center() const {return chi(0.5, 0.5);};

    inline const MappingQuad2d* map() const {return map_;};
    const 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_; }
    protected:
    virtual std::ostream& info(std::ostream& os) const;

    MappingQuad2d* map_;
    Index idx_;

    Quad2d* createChild_(Quad& cntr, const Index& idx, bool flag) const;

    Quad2d(Quad& cntr, MappingQuad2d* map, const Index& idx);
    private:

    Quad& cntr_;

    Quad2d* lnk_;

    Quad2d* chld_;

    const Quad2dSubdivision* subdivStrategy_;

    // Edges as cells
    mutable concepts::Array<Edge2d*> edge_;

    virtual Quad2d* createChild_(Quad& cntr, const Index& idx) const;
  };

  Quad& Quad2d::connector() const { 
    return cntr_; 
  }

  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:
      struct Index {
        uchar  l_;

        ushort i_;

        ushort j_;

        Index();
        Index(uchar l, ushort i, ushort 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);


} // namespace concepts

#endif // cell2D_hh