Provided by: librheolef-dev_6.7-6_amd64 bug

NAME

       csr - compressed sparse row matrix (rheolef-6.7)

SYNOPSYS

       Distributed compressed sparse matrix container stored row by row.

DESCRIPTION

       Sparse matrix are compressed by rows. In distributed environment, the distribution follows
       the row distributor (see distributor(2)).

ALGEBRA

       Adding or substracting two matrices writes a+b and a-b, respectively,  and  multiplying  a
       matrix  by  a  scalar writes lambda*x.  Thus, any linear combination of sparse matrices is
       available.

       Matrix-vector product writes a*x where x is a vector (see vec(2)).

LIMITATIONS

       Some basic linear algebra is still under  development:  a.trans_mult(x)  matrix  transpose
       vector product, trans(a) matrix transpose, a*b matrix product.

IMPLEMENTATION

       template<class T>
       class csr<T,sequential> : public smart_pointer<csr_rep<T,sequential> > {
       public:

       // typedefs:

           typedef csr_rep<T,sequential>             rep;
           typedef smart_pointer<rep>                base;
           typedef typename rep::memory_type         memory_type;
           typedef typename rep::size_type           size_type;
           typedef typename rep::element_type        element_type;
           typedef typename rep::iterator            iterator;
           typedef typename rep::const_iterator      const_iterator;
           typedef typename rep::data_iterator       data_iterator;
           typedef typename rep::const_data_iterator const_data_iterator;

       // allocators/deallocators:

           csr() : base(new_macro(rep())) {}
           template<class A>
           explicit csr(const asr<T,sequential,A>& a) : base(new_macro(rep(a))) {}
           void resize (size_type loc_nrow1 = 0, size_type loc_ncol1 = 0, size_type loc_nnz1 = 0)
               { base::data().resize(loc_nrow1, loc_ncol1, loc_nnz1); }
           void resize (const distributor& row_ownership, const distributor& col_ownership, size_type nnz1 = 0)
               { base::data().resize(row_ownership, col_ownership, nnz1); }

       // allocators from initializer list (c++ 2011):

       #ifdef _RHEOLEF_HAVE_STD_INITIALIZER_LIST
           csr (const std::initializer_list<csr_concat_value<T,sequential> >& init_list);
           csr (const std::initializer_list<csr_concat_line<T,sequential> >&  init_list);
       #endif // _RHEOLEF_HAVE_STD_INITIALIZER_LIST

       // accessors:

           // global sizes
           const distributor& row_ownership() const { return base::data().row_ownership(); }
           const distributor& col_ownership() const { return base::data().col_ownership(); }
           size_type dis_nrow () const              { return row_ownership().dis_size(); }
           size_type dis_ncol () const              { return col_ownership().dis_size(); }
           size_type dis_nnz () const               { return base::data().nnz(); }
           size_type dis_ext_nnz () const           { return 0; }
           bool is_symmetric() const                { return base::data().is_symmetric(); }
           void set_symmetry (bool is_symm) const   { base::data().set_symmetry(is_symm); }
           void set_symmetry_by_check (const T& tol = std::numeric_limits<T>::epsilon()) const
                                                    { base::data().set_symmetry_by_check(); }
           bool is_definite_positive() const        { return base::data().is_definite_positive(); }
           void set_definite_positive (bool is_defpos) const { base::data().set_definite_positive(is_defpos); }
           size_type pattern_dimension() const      { return base::data().pattern_dimension(); }
           void set_pattern_dimension(size_type dim) const { base::data().set_pattern_dimension(dim); }
           T max_abs () const                       { return base::data().max_abs(); }

           // local sizes
           size_type nrow () const                  { return base::data().nrow(); }
           size_type ncol () const                  { return base::data().ncol(); }
           size_type nnz () const                   { return base::data().nnz(); }

           // range on local memory
           size_type row_first_index () const       { return base::data().row_first_index(); }
           size_type row_last_index () const        { return base::data().row_last_index(); }
           size_type col_first_index () const       { return base::data().col_first_index(); }
           size_type col_last_index () const        { return base::data().col_last_index(); }

           const_iterator begin() const             { return base::data().begin(); }
           const_iterator end()   const             { return base::data().end(); }
           iterator begin_nonconst()                { return base::data().begin(); }
           iterator end_nonconst()                  { return base::data().end(); }

           // accessors, only for distributed (for interface compatibility)
           size_type ext_nnz() const                { return 0; }
           const_iterator ext_begin() const         { return const_iterator(); }
           const_iterator ext_end()   const         { return const_iterator(); }
                 iterator ext_begin_nonconst()      { return iterator(); }
                 iterator ext_end_nonconst()        { return iterator(); }
           size_type jext2dis_j (size_type jext) const { return 0; }

       // algebra:

           // y := a*x
           void mult (const vec<element_type,sequential>& x, vec<element_type,sequential>& y) const {
             base::data().mult (x,y);
           }
           vec<element_type,sequential> operator* (const vec<element_type,sequential>& x) const {
             vec<element_type,sequential> y (row_ownership(), element_type());
             mult (x, y);
             return y;
           }
           void trans_mult (const vec<element_type,sequential>& x, vec<element_type,sequential>& y) const {
             base::data().trans_mult (x,y);
           }
           vec<element_type,sequential> trans_mult (const vec<element_type,sequential>& x) const {
             vec<element_type,sequential> y (col_ownership(), element_type());
             trans_mult (x, y);
             return y;
           }
           // a+b, a-b, a*b
           csr<T,sequential> operator+ (const csr<T,sequential>& b) const;
           csr<T,sequential> operator- (const csr<T,sequential>& b) const;
           csr<T,sequential> operator* (const csr<T,sequential>& b) const;

           // lambda*a
           csr<T,sequential>& operator*= (const T& lambda) {
             base::data().operator*= (lambda);
             return *this;
           }
       // output:

           void dump (const std::string& name) const { base::data().dump(name); }
       };
       // lambda*a
       template<class T>
       inline
       csr<T,sequential>
       operator* (const T& lambda, const csr<T,sequential>& a)
       {
         csr<T,sequential> b = a;
         b.operator*= (lambda);
         return b;
       }
       // -a
       template<class T>
       inline
       csr<T,sequential>
       operator- (const csr<T,sequential>& a)
       {
         return T(-1)*a;
       }
       // trans(a)
       template<class T>
       inline
       csr<T,sequential>
       trans (const csr<T,sequential>& a)
       {
         csr<T,sequential> b;
         a.data().build_transpose (b.data());
         return b;
       }

IMPLEMENTATION

       template<class T>
       class csr<T,distributed> : public smart_pointer<csr_rep<T,distributed> > {
       public:

       // typedefs:

           typedef csr_rep<T,distributed>                    rep;
           typedef smart_pointer<rep>                base;
           typedef typename rep::memory_type         memory_type;
           typedef typename rep::size_type           size_type;
           typedef typename rep::element_type        element_type;
           typedef typename rep::iterator            iterator;
           typedef typename rep::const_iterator      const_iterator;
           typedef typename rep::data_iterator       data_iterator;
           typedef typename rep::const_data_iterator const_data_iterator;

       // allocators/deallocators:

           csr() : base(new_macro(rep())) {}
           template<class A>
           explicit csr(const asr<T,memory_type,A>& a) : base(new_macro(rep(a))) {}
           void resize (const distributor& row_ownership, const distributor& col_ownership, size_type nnz1 = 0)
               { base::data().resize(row_ownership, col_ownership, nnz1); }

       // allocators from initializer list (c++ 2011):

       #ifdef _RHEOLEF_HAVE_STD_INITIALIZER_LIST
           csr (const std::initializer_list<csr_concat_value<T,distributed> >& init_list);
           csr (const std::initializer_list<csr_concat_line<T,distributed> >&  init_list);
       #endif // _RHEOLEF_HAVE_STD_INITIALIZER_LIST

       // accessors:

           // global sizes
           const distributor& row_ownership() const { return base::data().row_ownership(); }
           const distributor& col_ownership() const { return base::data().col_ownership(); }
           size_type dis_nrow () const              { return row_ownership().dis_size(); }
           size_type dis_ncol () const              { return col_ownership().dis_size(); }
           size_type dis_nnz () const               { return base::data().dis_nnz(); }
           size_type dis_ext_nnz () const           { return base::data().dis_ext_nnz(); }
           bool is_symmetric() const                { return base::data().is_symmetric(); }
           void set_symmetry (bool is_symm) const   { base::data().set_symmetry(is_symm); }
           void set_symmetry_by_check (const T& tol = std::numeric_limits<T>::epsilon()) const
                                                    { base::data().set_symmetry_by_check(); }
           bool is_definite_positive() const        { return base::data().is_definite_positive(); }
           void set_definite_positive (bool is_defpos) const { base::data().set_definite_positive(is_defpos); }
           size_type pattern_dimension() const      { return base::data().pattern_dimension(); }
           void set_pattern_dimension(size_type dim) const { base::data().set_pattern_dimension(dim); }
           T max_abs () const                       { return base::data().max_abs(); }

           // local sizes
           size_type nrow () const                  { return base::data().nrow(); }
           size_type ncol () const                  { return base::data().ncol(); }
           size_type nnz () const                   { return base::data().nnz(); }

           // range on local memory
           size_type row_first_index () const       { return base::data().row_first_index(); }
           size_type row_last_index () const        { return base::data().row_last_index(); }
           size_type col_first_index () const       { return base::data().col_first_index(); }
           size_type col_last_index () const        { return base::data().col_last_index(); }

           const_iterator begin() const             { return base::data().begin(); }
           const_iterator end()   const             { return base::data().end(); }
           iterator begin_nonconst()                { return base::data().begin(); }
           iterator end_nonconst()                  { return base::data().end(); }

           // accessors, only for distributed
           size_type ext_nnz() const                { return base::data().ext_nnz(); }
           const_iterator ext_begin() const         { return base::data().ext_begin(); }
           const_iterator ext_end()   const         { return base::data().ext_end(); }
                 iterator ext_begin_nonconst()      { return base::data().ext_begin(); }
                 iterator ext_end_nonconst()        { return base::data().ext_end(); }
           size_type jext2dis_j (size_type jext) const { return base::data().jext2dis_j(jext); }

       // algebra:

           // y := a*x
           void mult (const vec<element_type,distributed>& x, vec<element_type,distributed>& y) const {
             base::data().mult (x,y);
           }
           vec<element_type,distributed> operator* (const vec<element_type,distributed>& x) const {
             vec<element_type,distributed> y (row_ownership(), element_type());
             mult (x, y);
             return y;
           }
           void trans_mult (const vec<element_type,distributed>& x, vec<element_type,distributed>& y) const {
             base::data().trans_mult (x,y);
           }
           vec<element_type,distributed> trans_mult (const vec<element_type,distributed>& x) const {
             vec<element_type,distributed> y (col_ownership(), element_type());
             trans_mult (x, y);
             return y;
           }
           // a+b, a-b, a*b
           csr<T,distributed> operator+ (const csr<T,distributed>& b) const;
           csr<T,distributed> operator- (const csr<T,distributed>& b) const;
           csr<T,distributed> operator* (const csr<T,distributed>& b) const;

           // lambda*a
           csr<T,distributed>& operator*= (const T& lambda) {
             base::data().operator*= (lambda);
             return *this;
           }
       // output:

           void dump (const std::string& name) const { base::data().dump(name); }
       };
       // lambda*a
       template<class T>
       inline
       csr<T,distributed>
       operator* (const T& lambda, const csr<T,distributed>& a)
       {
         csr<T,distributed> b = a;
         b.operator*= (lambda);
         return b;
       }
       // -a
       template<class T>
       inline
       csr<T,distributed>
       operator- (const csr<T,distributed>& a)
       {
         return T(-1)*a;
       }
       // trans(a)
       template<class T>
       inline
       csr<T,distributed>
       trans (const csr<T,distributed>& a)
       {
         csr<T,distributed> b;
         a.data().build_transpose (b.data());
         return b;
       }
       #endif // _RHEOLEF_HAVE_MPI

       // b = f(a); f as a class-function or usual fct
       template<class T, class M, class Function>
       csr<T,M>
       apply (Function f, const csr<T,M>& a)
       {
         csr<T,M> b = a;
         typename csr<T,M>::size_type n = a.nrow();
         typename csr<T,M>::const_iterator dia_ia = a.begin();
         typename csr<T,M>::iterator       dia_ib = b.begin_nonconst();
         pair_transform_second (dia_ia[0], dia_ia[n], dia_ib[0], f);
         if (a.ext_nnz() != 0) {
           typename csr<T,M>::const_iterator ext_ia = a.ext_begin();
           typename csr<T,M>::iterator       ext_ib = b.ext_begin_nonconst();
           pair_transform_second (ext_ia[0], ext_ia[n], ext_ib[0], f);
         }
         return b;
       }
       template<class T, class M, class Function>
       csr<T,M>
       apply (T (*f)(const T&), const csr<T,M>& a)
       {
         return apply (std::ptr_fun(f), a);
       }

SEE ALSO

       distributor(2), vec(2)