submatrix.hh

Go to the documentation of this file.

 
#ifndef submatrix_hh
#define submatrix_hh
 
#include "matrix.hh"
#include "matrixMult.hh"
#include "subMatrixIterator.hh"
#include "CRS.hh"
#include "space/spaceSet.hh"
#include "basics/outputOperator.hh"
 
#include <type_traits>
#include <typeinfo>
 
namespace concepts {
 
  // forward declaration
  template<class F, class G>
  class SubMatrix;
 
  template<class F>
  class SubMatrixN;
 
 
  // ************************************************************* SubMatrix **
 
  template<class F>
  class SubMatrixN : public Matrix<typename F::type>,
                     public CRSConvertable<typename F::type> {
  public:
    typedef typename F::type type;
    typedef typename Realtype<type>::type r_type;
    typedef typename Cmplxtype<type>::type c_type;
    typedef typename std::conditional<std::is_same<typename Realtype<F>::type, F>::value ,
      typename Realtype<F>::type, typename Cmplxtype<F>::type >::type d_type;
 
    typedef _SubMatrix_iterator<F, type&, type*>              iterator;
    typedef _SubMatrix_iterator<F, const type&, const type*>  const_iterator;
 
    SubMatrixN(F& m,
               const Set<IndexRange>& indicesX,
               const Set<IndexRange>& indicesY);
    SubMatrixN();
                  
    virtual ~SubMatrixN() {}
 
    virtual void operator()(const Function<r_type>& fncY,
                            Function<type>& fncX);
    virtual void operator()(const Function<c_type>& fncY,
                            Function<c_type>& fncX);
 
    virtual void transpMult(const Vector<r_type>& fncY,
                            Vector<type>& fncX);
    virtual void transpMult(const Vector<c_type>& fncY,
                            Vector<c_type>& fncX);
 
    template<class H, class I>
    void multiply(const H& fact, I& dest) const {
      matrixMultiplyRowSorting(*this, fact, dest);
    }
 
    virtual type operator()(const uint i, const uint j) const;
    virtual type& operator()(const uint i, const uint j);
 
    iterator begin(uint r = 0);
    const_iterator begin(uint r = 0) const;
    const_iterator end() const;
 
    uint used() const;
 
    virtual void convertCRS(type* a, int* asub, int* xa) const;
    virtual void convertCCS(type* a, int* asub, int* xa) const;
    virtual void convertIJK(type* , int* , int* ) const;
 
    template<class H, class I>
    void addInto(Matrix<H>& dest, const I fact,
                 const uint rowoffset = 0, const uint coloffset = 0) const;
 
    const Set<IndexRange>& indicesX() const { return indicesX_; }
    const Set<IndexRange>& indicesY() const { return indicesY_; }
 
    F& matrix() { return *m_; }
    const F& matrix() const { return *m_; }
  protected:
    virtual std::ostream& info(std::ostream& os) const;
    const Set<IndexRange> indicesX_, indicesY_;
  private:
    F* m_;
  };
 
  template<class F>
  template<class H, class I>
  void SubMatrixN<F>::addInto(Matrix<H>& dest, const I fact,
                              const uint rowoffset,
                              const uint coloffset) const {
    for(const_iterator i = begin(); i != end(); ++i)
      dest(rowoffset + i.row(), coloffset + i.col()) += fact * *i;
  }
 
} // namespace concepts
 
#endif // submatrix_hh