operator/submatrix.hh

Go to the documentation of this file.

/* interface for a submatrix and implemention for sparse matrix
 */

#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"

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

    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