hp2D/spacePreBuilder.hh

Go to the documentation of this file.

/* space pre builder classes for 2D spaces
*/

#ifndef hpSpacePreBuilder2D_hh
#define hpSpacePreBuilder2D_hh

#include <set>
#include <queue>
#include "basics/outputOperator.hh"
#include "geometry/cell2D.hh"
#include "space/spacePreBuilder.hh"
#include "space/space.hh"
#include "space/hpMethod.hh"          // for AdaptiveControlP
#include "toolbox/stiffArray.hh"
#include "toolbox/multiArray.hh"
#include "toolbox/hashMap.hh"
#include "toolbox/sequence.hh"
#include "toolbox/set.hh"
#include "hp2D/edgeMinMax.hh"

// If flag is set to 1, edgelist_ and vertexlist_ collect keys, not
// local indices. That is easier to debug, but lightly slower.
#define KeyLists 1

namespace hp2D {

  // *********************************************************** NotValidDof **

  //   class NotValidDof : public concepts::ExceptionBase {};

  // ******************************************************* SpacePreBuilder **

  class SpacePreBuilder : public concepts::SpacePreBuilder {
  public:
    SpacePreBuilder(concepts::Mesh2& msh) : msh_(msh) {}
    virtual ~SpacePreBuilder() {}

    virtual concepts::Mesh2& mesh() const { return msh_; }
    virtual bool vtxDof(const concepts::Connector2& cntr, uint i) const = 0;
    virtual ushort edgeDof(const concepts::Connector2& cntr, uint i) const = 0;
    virtual const ushort* innerDof(const concepts::Connector2& cntr) const = 0;
    virtual const ushort* pMax(const concepts::Connector2& cntr,
                               ushort* p = 0) const = 0;
    virtual ushort edgeP(const concepts::Connector2& cntr, uint i) const = 0;
  protected:
    concepts::Mesh2& msh_;
  };
  
  typedef concepts::Adaptivity<concepts::Connector,
                               concepts::AdaptiveAdjustP<2> > Adaptivity;

  // *********************************************************** DegreeDim **

  template<typename F>
  struct DegreeDim {};

  template<>
  struct DegreeDim<concepts::Quad>         { static const ushort dim = 2; };

  template<>
  struct DegreeDim<concepts::InfiniteQuad> { static const ushort dim = 2; };

  template<>
  struct DegreeDim<concepts::Triangle>     { static const ushort dim = 1; };
  

  // **************************************************************** hpFull **

  class hpFull : public SpacePreBuilder, public Adaptivity {
  public:
    hpFull(concepts::Mesh2& msh, uint l = 0, uint p = 1, bool build = false);
    hpFull(const hpFull& preBuild);
    virtual ~hpFull() {}

    void setPolyEdgeRule(const PolyEdgeRule& rule);

    virtual void adjust(const concepts::Connector& cell,
                        const concepts::AdaptiveAdjustP<2>& a);
    void adjust(const uint nrCell, const concepts::AdaptiveAdjustP<2>& a);
    void adjust(const concepts::AdaptiveAdjustP<2>& a);

    uint available() const { return rebuild_ ? 0 : build_; }
    virtual void rebuildMesh();
    virtual void rebuildDof();
    virtual bool vtxDof(const concepts::Connector2& cntr, uint i) const;
    virtual ushort edgeDof(const concepts::Connector2& cntr, uint i) const;
    virtual const ushort* innerDof(const concepts::Connector2& cntr) const;
    virtual const ushort* pMax(const concepts::Connector2& cntr,
                               ushort* p = 0) const;
    virtual ushort edgeP(const concepts::Connector2& cntr, uint i) const;
    const concepts::IndexRange& index(uint dim,
                                      const concepts::Connector2& cntr,
                                      uint i = 0, uint spcNo = 0) const;
    virtual concepts::IndexRange&
      setIndex(uint& firstIdx, uint noIdx, uint dim,
         const concepts::Connector& cntr, uint i = 0,
         uint spcNo = 0) throw(concepts::NotValidDof);

    virtual concepts::Set<concepts::IndexRange> indices
      (uint dim, const concepts::Connector& cntr, uint spcNo = 0) const
      throw(concepts::MissingFeature);

    void clearAllIndices(uint spcNo = 0);

    uint ncells() { return cells_.size(); }

    bool isHangEdge(const concepts::Connector1& cntr) const;
    void storeMatlab(const std::string name);

    concepts::Set<concepts::Cell2*> cells() const { return cells_; }
    std::set<concepts::Cell2*> allCells() const { return allCells_; }
  protected:
    virtual std::ostream& info(std::ostream& os) const;
  private:
    bool rebuild_;
    uint build_;
    std::auto_ptr<const PolyEdgeRule> polyEdgeRule_;

    concepts::HashMap<concepts::AdaptiveAdjustP<2> > adj_;

    // mapping from space number to an index range
    typedef std::map<uint, concepts::IndexRange> IndexType;
    // index range per space number for all 'member of space' vertices dof
    typedef concepts::AdaptiveControl<IndexType> CtrlType_0;
    // index range per space number for all 'member of space' edge dof
    typedef concepts::AdaptiveControlP<1, IndexType> CtrlType_1;
    // index range per space number for all inner dof
    typedef concepts::AdaptiveControlP<2, IndexType> CtrlType_2;

    concepts::HashMap<CtrlType_0> ctrl0_;
    concepts::HashMap<CtrlType_1> ctrl1_;
    concepts::HashMap<CtrlType_2> ctrl2_;

    concepts::HashMap<ushort> polyDep_;
    concepts::HashMap<ushort> polyHang_;
    mutable concepts::HashMap<concepts::StiffArray<2,ushort> > pmax_;
    mutable concepts::HashMap<concepts::StiffArray<2,ushort> > pminInner_;

    concepts::Set<concepts::Cell2*> cells_;
    std::set<concepts::Cell2*> allCells_;
    concepts::HashMap<std::set<uint> > vertexList_;
    concepts::HashMap<std::set<uint> > edgeList_;
    concepts::Set<uint> hangEdges_;

    bool isMember_(const concepts::Connector2& cntr) const;
    bool isMember_(const concepts::Connector1& cntr) const;

    const ushort* pMinInner_(const concepts::Quad& quad);

    void refine_(concepts::Quad2d& cell, concepts::Level<2> L, short* P,
                 concepts::HashMap<CtrlType_2>::iterator& i2);
    void refine_(concepts::InfiniteQuad2d& cell, concepts::Level<2> L,
                 short* P, concepts::HashMap<CtrlType_2>::iterator& i2);
    void refine_(concepts::Cell2& cell, ushort* P,
                 concepts::HashMap<CtrlType_2>::iterator& i2);

    void memberEdges_(const concepts::Connector2& cntr,
                      std::set<const concepts::Connector1*>& edges);

    ushort takeOutchildrenEdges_
    (const concepts::Connector1& edge, const concepts::Connector1& cntr1,
     std::set<const concepts::Connector1*>& edges,
     std::set<const concepts::Connector0*>& vertices);

    template<class F>
    bool getEdgeSupport_(F cntr);

    template<class F>
    bool getPretendedVertexSupport_
    (F cntr,
     concepts::HashMap<std::queue<concepts::StiffArray<2, uint> > >& nonLocal);

    void getRootCells_
    (const concepts::Connector2& cntr, const std::set<uint>& vertices,
     std::map<uint, std::vector<const concepts::Connector2*> >& rootCells);

    void getNonLocalVertexSupport_
    (concepts::HashMap<std::queue<concepts::StiffArray<2, uint> > >& nonLocal,
     const std::map<uint, std::vector<const concepts::Connector2*> >&
     rootCells);

    void getEdgeMappings_
    (uint key, const std::vector<const concepts::Connector2*>& root,
     concepts::MultiArray<2, const concepts::Connector2*>& edgesToCell, 
     std::map<uint,
       concepts::Sequence<const concepts::Connector2*> >& edgeToCell,
     std::map<uint, uint>& fatherEdges, std::map<uint, uint>& mosEdges) const;

    void dependentEdges_(const concepts::Quad& quad);
    void polyDependentEdges_(const concepts::Quad& quad,
           const concepts::Edge& edge,
           const concepts::Edge& oppEdge, const uint p);

    ushort* pMax_(const concepts::Quad& cntr) const;
    ushort* pMax_(const concepts::InfiniteQuad& cntr) const;

    void setPmax_(const concepts::Connector2& cntr, ushort* Pmax);

    void polyHangingEdges_(const concepts::Connector1& edge, const uint& p);

    template<class F>
    void findControl_(F*& ctrl, const concepts::Connector2& cntr, uint i);

          IndexType& findIndex_(uint dim, const concepts::Connector2& cntr, 
                                uint i);
    const IndexType& findIndex_(uint dim, const concepts::Connector2& cntr,
                                uint i) const;
  };

} // namespace hp2D

#endif // hpSpacePreBuilder2D_hh