Chris@49: // Copyright (C) 2012 Ryan Curtin
Chris@49: // Copyright (C) 2012 Conrad Sanderson
Chris@49: //
Chris@49: // This Source Code Form is subject to the terms of the Mozilla Public
Chris@49: // License, v. 2.0. If a copy of the MPL was not distributed with this
Chris@49: // file, You can obtain one at http://mozilla.org/MPL/2.0/.
Chris@49: 
Chris@49: //! \addtogroup SpSubview
Chris@49: //! @{
Chris@49: 
Chris@49: ///////////////////////////////////////////////////////////////////////////////
Chris@49: // SpSubview::iterator_base implementation                                   //
Chris@49: ///////////////////////////////////////////////////////////////////////////////
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: SpSubview<eT>::iterator_base::iterator_base(const SpSubview<eT>& in_M)
Chris@49:   : M(in_M)
Chris@49:   , internal_col(0)
Chris@49:   , internal_pos(0)
Chris@49:   , skip_pos(0)
Chris@49:   {
Chris@49:   // Technically this iterator is invalid (it may not point to a real element).
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: SpSubview<eT>::iterator_base::iterator_base(const SpSubview<eT>& in_M, const uword in_col, const uword in_pos, const uword in_skip_pos)
Chris@49:   : M(in_M)
Chris@49:   , internal_col(in_col)
Chris@49:   , internal_pos(in_pos)
Chris@49:   , skip_pos    (in_skip_pos)
Chris@49:   {
Chris@49:   // Nothing to do.
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: eT
Chris@49: SpSubview<eT>::iterator_base::operator*() const
Chris@49:   {
Chris@49:   return M.m.values[internal_pos + skip_pos];
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: ///////////////////////////////////////////////////////////////////////////////
Chris@49: // SpSubview::const_iterator implementation                                  //
Chris@49: ///////////////////////////////////////////////////////////////////////////////
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: SpSubview<eT>::const_iterator::const_iterator(const SpSubview<eT>& in_M, const uword initial_pos)
Chris@49:   : iterator_base(in_M, 0, initial_pos, 0)
Chris@49:   {
Chris@49:   // Corner case for empty subviews.
Chris@49:   if(in_M.n_nonzero == 0)
Chris@49:     {
Chris@49:     iterator_base::internal_col = in_M.n_cols;
Chris@49:     iterator_base::skip_pos     = in_M.m.n_nonzero;
Chris@49:     return;
Chris@49:     }
Chris@49: 
Chris@49:   // Figure out the row and column of the position.
Chris@49:   // lskip_pos holds the number of values which aren't part of this subview.
Chris@49:   const uword aux_col = iterator_base::M.aux_col1;
Chris@49:   const uword aux_row = iterator_base::M.aux_row1;
Chris@49:   const uword ln_rows = iterator_base::M.n_rows;
Chris@49:   const uword ln_cols = iterator_base::M.n_cols;
Chris@49: 
Chris@49:   uword cur_pos   = 0; // off by one because we might be searching for pos 0
Chris@49:   uword lskip_pos = iterator_base::M.m.col_ptrs[aux_col];
Chris@49:   uword cur_col   = 0;
Chris@49: 
Chris@49:   while(cur_pos < (iterator_base::internal_pos + 1))
Chris@49:     {
Chris@49:     // Have we stepped forward a column (or multiple columns)?
Chris@49:     while(((lskip_pos + cur_pos) >= iterator_base::M.m.col_ptrs[cur_col + aux_col + 1]) && (cur_col < ln_cols))
Chris@49:       {
Chris@49:       ++cur_col;
Chris@49:       }
Chris@49: 
Chris@49:     // See if the current position is in the subview.
Chris@49:     const uword row_index = iterator_base::M.m.row_indices[cur_pos + lskip_pos];
Chris@49:     if(row_index < aux_row)
Chris@49:       {
Chris@49:       ++lskip_pos; // not valid
Chris@49:       }
Chris@49:     else if(row_index < (aux_row + ln_rows))
Chris@49:       {
Chris@49:       ++cur_pos; // valid, in the subview
Chris@49:       }
Chris@49:     else
Chris@49:       {
Chris@49:       // skip to end of column
Chris@49:       const uword next_colptr = iterator_base::M.m.col_ptrs[cur_col + aux_col + 1];
Chris@49:       lskip_pos += (next_colptr - (cur_pos + lskip_pos));
Chris@49:       }
Chris@49:     }
Chris@49: 
Chris@49:   iterator_base::internal_col = cur_col;
Chris@49:   iterator_base::skip_pos     = lskip_pos;
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: SpSubview<eT>::const_iterator::const_iterator(const SpSubview<eT>& in_M, const uword in_row, const uword in_col)
Chris@49:   : iterator_base(in_M, in_col, 0, 0)
Chris@49:   {
Chris@49:   // Corner case for empty subviews.
Chris@49:   if(in_M.n_nonzero == 0)
Chris@49:     {
Chris@49:     // We must be at the last position.
Chris@49:     iterator_base::internal_col = in_M.n_cols;
Chris@49:     iterator_base::skip_pos = in_M.m.n_nonzero;
Chris@49:     return;
Chris@49:     }
Chris@49: 
Chris@49:   // We have a destination we want to be just after, but don't know what position that is.
Chris@49:   // Because we have to count the points in this subview and not in this subview, this becomes a little difficult and slow.
Chris@49:   const uword aux_col = iterator_base::M.aux_col1;
Chris@49:   const uword aux_row = iterator_base::M.aux_row1;
Chris@49:   const uword ln_rows = iterator_base::M.n_rows;
Chris@49:   const uword ln_cols = iterator_base::M.n_cols;
Chris@49: 
Chris@49:   uword cur_pos = 0;
Chris@49:   uword skip_pos = iterator_base::M.m.col_ptrs[aux_col];
Chris@49:   uword cur_col = 0;
Chris@49: 
Chris@49:   while(cur_col < in_col)
Chris@49:     {
Chris@49:     // See if the current position is in the subview.
Chris@49:     const uword row_index = iterator_base::M.m.row_indices[cur_pos + skip_pos];
Chris@49:     if(row_index < aux_row)
Chris@49:       {
Chris@49:       ++skip_pos;
Chris@49:       }
Chris@49:     else if(row_index < (aux_row + ln_rows))
Chris@49:       {
Chris@49:       ++cur_pos;
Chris@49:       }
Chris@49:     else
Chris@49:       {
Chris@49:       // skip to end of column
Chris@49:       const uword next_colptr = iterator_base::M.m.col_ptrs[cur_col + aux_col + 1];
Chris@49:       skip_pos += (next_colptr - (cur_pos + skip_pos));
Chris@49:       }
Chris@49: 
Chris@49:     // Have we stepped forward a column (or multiple columns)?
Chris@49:     while(((skip_pos + cur_pos) >= iterator_base::M.m.col_ptrs[cur_col + aux_col + 1]) && (cur_col < ln_cols))
Chris@49:       {
Chris@49:       ++cur_col;
Chris@49:       }
Chris@49:     }
Chris@49: 
Chris@49:   // Now we are either on the right column or ahead of it.
Chris@49:   if(cur_col == in_col)
Chris@49:     {
Chris@49:     // We have to find the right row index.
Chris@49:     uword row_index = iterator_base::M.m.row_indices[cur_pos + skip_pos];
Chris@49:     while((row_index < (in_row + aux_row)))
Chris@49:       {
Chris@49:       if(row_index < aux_row)
Chris@49:         {
Chris@49:         ++skip_pos;
Chris@49:         }
Chris@49:       else
Chris@49:         {
Chris@49:         ++cur_pos;
Chris@49:         }
Chris@49: 
Chris@49:       // Ensure we didn't step forward a column; if we did, we need to stop.
Chris@49:       while(((skip_pos + cur_pos) >= iterator_base::M.m.col_ptrs[cur_col + aux_col + 1]) && (cur_col < ln_cols))
Chris@49:         {
Chris@49:         ++cur_col;
Chris@49:         }
Chris@49: 
Chris@49:       if(cur_col != in_col)
Chris@49:         {
Chris@49:         break;
Chris@49:         }
Chris@49: 
Chris@49:       row_index = iterator_base::M.m.row_indices[cur_pos + skip_pos];
Chris@49:       }
Chris@49:     }
Chris@49: 
Chris@49:   // Now we need to find the next valid position in the subview.
Chris@49:   uword row_index;
Chris@49:   while(true)
Chris@49:     {
Chris@49:     const uword next_colptr = iterator_base::M.m.col_ptrs[cur_col + aux_col + 1];
Chris@49:     row_index = iterator_base::M.m.row_indices[cur_pos + skip_pos];
Chris@49: 
Chris@49:     // Are we at the last position?
Chris@49:     if(cur_col >= ln_cols)
Chris@49:       {
Chris@49:       cur_col = ln_cols;
Chris@49:       // Make sure we will be pointing at the last element in the parent matrix.
Chris@49:       skip_pos = iterator_base::M.m.n_nonzero - iterator_base::M.n_nonzero;
Chris@49:       break;
Chris@49:       }
Chris@49: 
Chris@49:     if(row_index < aux_row)
Chris@49:       {
Chris@49:       ++skip_pos;
Chris@49:       }
Chris@49:     else if(row_index < (aux_row + ln_rows))
Chris@49:       {
Chris@49:       break; // found
Chris@49:       }
Chris@49:     else
Chris@49:       {
Chris@49:       skip_pos += (next_colptr - (cur_pos + skip_pos));
Chris@49:       }
Chris@49: 
Chris@49:     // Did we move any columns?
Chris@49:     while(((skip_pos + cur_pos) >= iterator_base::M.m.col_ptrs[cur_col + aux_col + 1]) && (cur_col < ln_cols))
Chris@49:       {
Chris@49:       ++cur_col;
Chris@49:       }
Chris@49:     }
Chris@49: 
Chris@49:   // It is possible we have moved another column.
Chris@49:   while(((skip_pos + cur_pos) >= iterator_base::M.m.col_ptrs[cur_col + aux_col + 1]) && (cur_col < ln_cols))
Chris@49:     {
Chris@49:     ++cur_col;
Chris@49:     }
Chris@49: 
Chris@49:   iterator_base::internal_pos = cur_pos;
Chris@49:   iterator_base::skip_pos     = skip_pos;
Chris@49:   iterator_base::internal_col = cur_col;
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: SpSubview<eT>::const_iterator::const_iterator(const SpSubview<eT>& in_M, uword in_row, uword in_col, uword in_pos, uword in_skip_pos)
Chris@49:   : iterator_base(in_M, in_col, in_pos, in_skip_pos)
Chris@49:   {
Chris@49:   // Nothing to do.
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: SpSubview<eT>::const_iterator::const_iterator(const const_iterator& other)
Chris@49:   : iterator_base(other.M, other.internal_col, other.internal_pos, other.skip_pos)
Chris@49:   {
Chris@49:   // Nothing to do.
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: typename SpSubview<eT>::const_iterator&
Chris@49: SpSubview<eT>::const_iterator::operator++()
Chris@49:   {
Chris@49:   const uword aux_col = iterator_base::M.aux_col1;
Chris@49:   const uword aux_row = iterator_base::M.aux_row1;
Chris@49:   const uword ln_rows = iterator_base::M.n_rows;
Chris@49:   const uword ln_cols = iterator_base::M.n_cols;
Chris@49: 
Chris@49:   uword cur_col   = iterator_base::internal_col;
Chris@49:   uword cur_pos   = iterator_base::internal_pos + 1;
Chris@49:   uword lskip_pos = iterator_base::skip_pos;
Chris@49:   uword row_index;
Chris@49: 
Chris@49:   while(true)
Chris@49:     {
Chris@49:     const uword next_colptr = iterator_base::M.m.col_ptrs[cur_col + aux_col + 1];
Chris@49:     row_index = iterator_base::M.m.row_indices[cur_pos + lskip_pos];
Chris@49: 
Chris@49:     // Did we move any columns?
Chris@49:     while((cur_col < ln_cols) && ((lskip_pos + cur_pos) >= iterator_base::M.m.col_ptrs[cur_col + aux_col + 1]))
Chris@49:       {
Chris@49:       ++cur_col;
Chris@49:       }
Chris@49: 
Chris@49:     // Are we at the last position?
Chris@49:     if(cur_col >= ln_cols)
Chris@49:       {
Chris@49:       cur_col = ln_cols;
Chris@49:       // Make sure we will be pointing at the last element in the parent matrix.
Chris@49:       lskip_pos = iterator_base::M.m.n_nonzero - iterator_base::M.n_nonzero;
Chris@49:       break;
Chris@49:       }
Chris@49: 
Chris@49:     if(row_index < aux_row)
Chris@49:       {
Chris@49:       ++lskip_pos;
Chris@49:       }
Chris@49:     else if(row_index < (aux_row + ln_rows))
Chris@49:       {
Chris@49:       break; // found
Chris@49:       }
Chris@49:     else
Chris@49:       {
Chris@49:       lskip_pos += (next_colptr - (cur_pos + lskip_pos));
Chris@49:       }
Chris@49:     }
Chris@49: 
Chris@49:   iterator_base::internal_pos = cur_pos;
Chris@49:   iterator_base::internal_col = cur_col;
Chris@49:   iterator_base::skip_pos     = lskip_pos;
Chris@49: 
Chris@49:   return *this;
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: typename SpSubview<eT>::const_iterator
Chris@49: SpSubview<eT>::const_iterator::operator++(int)
Chris@49:   {
Chris@49:   typename SpSubview<eT>::const_iterator tmp(*this);
Chris@49: 
Chris@49:   ++(*this);
Chris@49: 
Chris@49:   return tmp;
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: typename SpSubview<eT>::const_iterator&
Chris@49: SpSubview<eT>::const_iterator::operator--()
Chris@49:   {
Chris@49:   const uword aux_col = iterator_base::M.aux_col1;
Chris@49:   const uword aux_row = iterator_base::M.aux_row1;
Chris@49:   const uword ln_rows = iterator_base::M.n_rows;
Chris@49: 
Chris@49:   uword cur_col  = iterator_base::internal_col;
Chris@49:   uword cur_pos  = iterator_base::internal_pos - 1;
Chris@49:   uword skip_pos = iterator_base::skip_pos;
Chris@49: 
Chris@49:   // Special condition for end of iterator.
Chris@49:   if((skip_pos + cur_pos + 1) == iterator_base::M.m.n_nonzero)
Chris@49:     {
Chris@49:     // We are at the last element.  So we need to set skip_pos back to what it
Chris@49:     // would be if we didn't manually modify it back in operator++().
Chris@49:     skip_pos = iterator_base::M.m.col_ptrs[cur_col + aux_col] - iterator_base::internal_pos;
Chris@49:     }
Chris@49: 
Chris@49:   uword row_index;
Chris@49: 
Chris@49:   while(true)
Chris@49:     {
Chris@49:     const uword colptr = iterator_base::M.m.col_ptrs[cur_col + aux_col];
Chris@49:     row_index = iterator_base::M.m.row_indices[cur_pos + skip_pos];
Chris@49: 
Chris@49:     // Did we move back any columns?
Chris@49:     while((skip_pos + cur_pos) < iterator_base::M.m.col_ptrs[cur_col + aux_col])
Chris@49:       {
Chris@49:       --cur_col;
Chris@49:       }
Chris@49: 
Chris@49:     if(row_index < aux_row)
Chris@49:       {
Chris@49:       skip_pos -= (colptr - (cur_pos + skip_pos) + 1);
Chris@49:       }
Chris@49:     else if(row_index < (aux_row + ln_rows))
Chris@49:       {
Chris@49:       break; // found
Chris@49:       }
Chris@49:     else
Chris@49:       {
Chris@49:       --skip_pos;
Chris@49:       }
Chris@49:     }
Chris@49: 
Chris@49:   iterator_base::internal_pos = cur_pos;
Chris@49:   iterator_base::skip_pos     = skip_pos;
Chris@49:   iterator_base::internal_col = cur_col;
Chris@49: 
Chris@49:   return *this;
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: typename SpSubview<eT>::const_iterator
Chris@49: SpSubview<eT>::const_iterator::operator--(int)
Chris@49:   {
Chris@49:   typename SpSubview<eT>::const_iterator tmp(*this);
Chris@49: 
Chris@49:   --(*this);
Chris@49: 
Chris@49:   return tmp;
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: bool
Chris@49: SpSubview<eT>::const_iterator::operator==(const const_iterator& rhs) const
Chris@49:   {
Chris@49:   return (rhs.row() == (*this).row()) && (rhs.col() == iterator_base::internal_col);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: bool
Chris@49: SpSubview<eT>::const_iterator::operator!=(const const_iterator& rhs) const
Chris@49:   {
Chris@49:   return (rhs.row() != (*this).row()) || (rhs.col() != iterator_base::internal_col);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: bool
Chris@49: SpSubview<eT>::const_iterator::operator==(const typename SpMat<eT>::const_iterator& rhs) const
Chris@49:   {
Chris@49:   return (rhs.row() == (*this).row()) && (rhs.col() == iterator_base::internal_col);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: bool
Chris@49: SpSubview<eT>::const_iterator::operator!=(const typename SpMat<eT>::const_iterator& rhs) const
Chris@49:   {
Chris@49:   return (rhs.row() != (*this).row()) || (rhs.col() != iterator_base::internal_col);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: bool
Chris@49: SpSubview<eT>::const_iterator::operator==(const const_row_iterator& rhs) const
Chris@49:   {
Chris@49:   return (rhs.row() == (*this).row()) && (rhs.col() == iterator_base::internal_col);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: bool
Chris@49: SpSubview<eT>::const_iterator::operator!=(const const_row_iterator& rhs) const
Chris@49:   {
Chris@49:   return (rhs.row() != (*this).row()) || (rhs.col() != iterator_base::internal_col);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: bool
Chris@49: SpSubview<eT>::const_iterator::operator==(const typename SpMat<eT>::const_row_iterator& rhs) const
Chris@49:   {
Chris@49:   return (rhs.row() == (*this).row()) && (rhs.col() == iterator_base::internal_col);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: bool
Chris@49: SpSubview<eT>::const_iterator::operator!=(const typename SpMat<eT>::const_row_iterator& rhs) const
Chris@49:   {
Chris@49:   return (rhs.row() != (*this).row()) || (rhs.col() != iterator_base::internal_col);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: ///////////////////////////////////////////////////////////////////////////////
Chris@49: // SpSubview<eT>::iterator implementation                                    //
Chris@49: ///////////////////////////////////////////////////////////////////////////////
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: SpValProxy<SpSubview<eT> >
Chris@49: SpSubview<eT>::iterator::operator*()
Chris@49:   {
Chris@49:   return SpValProxy<SpSubview<eT> >(
Chris@49:     iterator_base::row(),
Chris@49:     iterator_base::col(),
Chris@49:     access::rw(iterator_base::M),
Chris@49:     &(access::rw(iterator_base::M.m.values[iterator_base::internal_pos + iterator_base::skip_pos])));
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: typename SpSubview<eT>::iterator&
Chris@49: SpSubview<eT>::iterator::operator++()
Chris@49:   {
Chris@49:   const_iterator::operator++();
Chris@49:   return *this;
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: typename SpSubview<eT>::iterator
Chris@49: SpSubview<eT>::iterator::operator++(int)
Chris@49:   {
Chris@49:   typename SpSubview<eT>::iterator tmp(*this);
Chris@49: 
Chris@49:   const_iterator::operator++();
Chris@49: 
Chris@49:   return tmp;
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: typename SpSubview<eT>::iterator&
Chris@49: SpSubview<eT>::iterator::operator--()
Chris@49:   {
Chris@49:   const_iterator::operator--();
Chris@49:   return *this;
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: typename SpSubview<eT>::iterator
Chris@49: SpSubview<eT>::iterator::operator--(int)
Chris@49:   {
Chris@49:   typename SpSubview<eT>::iterator tmp(*this);
Chris@49: 
Chris@49:   const_iterator::operator--();
Chris@49: 
Chris@49:   return tmp;
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: ///////////////////////////////////////////////////////////////////////////////
Chris@49: // SpSubview<eT>::const_row_iterator implementation                          //
Chris@49: ///////////////////////////////////////////////////////////////////////////////
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: SpSubview<eT>::const_row_iterator::const_row_iterator(const SpSubview<eT>& in_M, uword initial_pos)
Chris@49:   : iterator_base(in_M, 0, initial_pos, 0)
Chris@49:   , internal_row(0)
Chris@49:   , actual_pos(0)
Chris@49:   {
Chris@49:   // Corner case for empty subviews.
Chris@49:   if(in_M.n_nonzero == 0)
Chris@49:     {
Chris@49:     iterator_base::internal_col = 0;
Chris@49:     internal_row = in_M.n_rows;
Chris@49:     iterator_base::skip_pos = in_M.m.n_nonzero;
Chris@49:     return;
Chris@49:     }
Chris@49: 
Chris@49:   const uword aux_col = iterator_base::M.aux_col1;
Chris@49:   const uword aux_row = iterator_base::M.aux_row1;
Chris@49:   const uword ln_cols = iterator_base::M.n_cols;
Chris@49: 
Chris@49:   // We don't know where the elements are in each row.  What we will do is
Chris@49:   // loop across all valid columns looking for elements in row 0 (and add to
Chris@49:   // our sum), then in row 1, and so forth, until we get to the desired
Chris@49:   // position.
Chris@49:   uword cur_pos = -1;
Chris@49:   uword cur_row = 0;
Chris@49:   uword cur_col = 0;
Chris@49: 
Chris@49:   while(true)
Chris@49:     {
Chris@49:     // Is there anything in the column we are looking at?
Chris@49:     const uword colptr = iterator_base::M.m.col_ptrs[cur_col + aux_col];
Chris@49:     const uword next_colptr = iterator_base::M.m.col_ptrs[cur_col + aux_col + 1];
Chris@49: 
Chris@49:     for(uword ind = colptr; (ind < next_colptr) && (iterator_base::M.m.row_indices[ind] <= (cur_row + aux_row)); ++ind)
Chris@49:       {
Chris@49:       // There is something in this column.  Is it in the row we are looking at?
Chris@49:       const uword row_index = iterator_base::M.m.row_indices[ind];
Chris@49:       if(row_index == (cur_row + aux_row))
Chris@49:         {
Chris@49:         // Yes, it is in the right row.
Chris@49:         if(++cur_pos == iterator_base::internal_pos)
Chris@49:           {
Chris@49:           iterator_base::internal_col = cur_col;
Chris@49:           internal_row = cur_row;
Chris@49:           actual_pos = ind;
Chris@49: 
Chris@49:           return;
Chris@49:           }
Chris@49: 
Chris@49:         // We are done with this column.  Break to the column incrementing code (directly below).
Chris@49:         break;
Chris@49:         }
Chris@49:       else if(row_index > (cur_row + aux_row))
Chris@49:         {
Chris@49:         break; // Can't be in this column.
Chris@49:         }
Chris@49:       }
Chris@49: 
Chris@49:     cur_col++; // Done with the column.  Move on.
Chris@49:     if(cur_col == ln_cols)
Chris@49:       {
Chris@49:       // Out of columns.  Loop back to the beginning and look on the next row.
Chris@49:       cur_col = 0;
Chris@49:       cur_row++;
Chris@49:       }
Chris@49:     }
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: SpSubview<eT>::const_row_iterator::const_row_iterator(const SpSubview<eT>& in_M, uword in_row, uword in_col)
Chris@49:   : iterator_base(in_M, in_col, 0, 0)
Chris@49:   , internal_row(0)
Chris@49:   , actual_pos(0)
Chris@49:   {
Chris@49:   // We have a destination we want to be just after, but don't know what that
Chris@49:   // position is.  Because we will have to loop over everything anyway, create
Chris@49:   // another iterator and loop it until it is at the right place, then take its
Chris@49:   // information.
Chris@49:   const_row_iterator it(in_M, 0);
Chris@49:   while((it.row() < in_row) || ((it.row() == in_row) && (it.col() < in_col)))
Chris@49:     {
Chris@49:     ++it;
Chris@49:     }
Chris@49: 
Chris@49:   iterator_base::internal_col = it.col();
Chris@49:   iterator_base::internal_pos = it.pos();
Chris@49:   iterator_base::skip_pos = it.skip_pos;
Chris@49:   internal_row = it.internal_row;
Chris@49:   actual_pos = it.actual_pos;
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: SpSubview<eT>::const_row_iterator::const_row_iterator(const const_row_iterator& other)
Chris@49:   : iterator_base(other.M, other.internal_col, other.internal_pos, other.skip_pos)
Chris@49:   , internal_row(other.internal_row)
Chris@49:   , actual_pos(other.actual_pos)
Chris@49:   {
Chris@49:   // Nothing to do.
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: typename SpSubview<eT>::const_row_iterator&
Chris@49: SpSubview<eT>::const_row_iterator::operator++()
Chris@49:   {
Chris@49:   // We just need to find the next nonzero element.
Chris@49:   ++iterator_base::internal_pos;
Chris@49: 
Chris@49:   // If we have exceeded the bounds, update accordingly.
Chris@49:   if(iterator_base::internal_pos >= iterator_base::M.n_nonzero)
Chris@49:     {
Chris@49:     internal_row = iterator_base::M.n_rows;
Chris@49:     iterator_base::internal_col = 0;
Chris@49:     actual_pos = iterator_base::M.m.n_nonzero;
Chris@49: 
Chris@49:     return *this;
Chris@49:     }
Chris@49: 
Chris@49:   // Otherwise, we need to search.
Chris@49:   uword cur_col = iterator_base::internal_col;
Chris@49:   uword cur_row = internal_row;
Chris@49: 
Chris@49:   const uword aux_col = iterator_base::M.aux_col1;
Chris@49:   const uword aux_row = iterator_base::M.aux_row1;
Chris@49:   const uword ln_cols = iterator_base::M.n_cols;
Chris@49: 
Chris@49:   while(true)
Chris@49:     {
Chris@49:     // Increment the current column and see if we are on a new row.
Chris@49:     if(++cur_col == ln_cols)
Chris@49:       {
Chris@49:       cur_col = 0;
Chris@49:       ++cur_row;
Chris@49:       }
Chris@49: 
Chris@49:     // Is there anything in this new column?
Chris@49:     const uword colptr = iterator_base::M.m.col_ptrs[cur_col + aux_col];
Chris@49:     const uword next_colptr = iterator_base::M.m.col_ptrs[cur_col + aux_col + 1];
Chris@49: 
Chris@49:     for(uword ind = colptr; (ind < next_colptr) && (iterator_base::M.m.row_indices[ind] <= (cur_row + aux_row)); ++ind)
Chris@49:       {
Chris@49:       const uword row_index = iterator_base::M.m.row_indices[ind];
Chris@49: 
Chris@49:       if((row_index - aux_row) == cur_row)
Chris@49:         {
Chris@49:         // We have successfully incremented.
Chris@49:         internal_row = cur_row;
Chris@49:         actual_pos = ind;
Chris@49:         iterator_base::internal_col = cur_col;
Chris@49: 
Chris@49:         return *this;
Chris@49:         }
Chris@49:       }
Chris@49:     }
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: typename SpSubview<eT>::const_row_iterator
Chris@49: SpSubview<eT>::const_row_iterator::operator++(int)
Chris@49:   {
Chris@49:   typename SpSubview<eT>::const_row_iterator tmp(*this);
Chris@49: 
Chris@49:   ++(*this);
Chris@49: 
Chris@49:   return tmp;
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: typename SpSubview<eT>::const_row_iterator&
Chris@49: SpSubview<eT>::const_row_iterator::operator--()
Chris@49:   {
Chris@49:   // We just need to find the previous element.
Chris@49: //  if(iterator_base::pos == 0)
Chris@49: //    {
Chris@49: //    // We cannot decrement.
Chris@49: //    return *this;
Chris@49: //    }
Chris@49: //  else if(iterator_base::pos == iterator_base::M.n_nonzero)
Chris@49: //    {
Chris@49: //    // We will be coming off the last element.  We need to reset the row correctly, because we set row = 0 in the last matrix position.
Chris@49: //    iterator_base::row = iterator_base::M.n_rows - 1;
Chris@49: //    }
Chris@49: //  else if(iterator_base::pos > iterator_base::M.n_nonzero)
Chris@49: //    {
Chris@49: //    // This shouldn't happen...
Chris@49: //    iterator_base::pos--;
Chris@49: //    return *this;
Chris@49: //    }
Chris@49: 
Chris@49:   iterator_base::internal_pos--;
Chris@49: 
Chris@49:   // We have to search backwards.
Chris@49:   uword cur_col = iterator_base::internal_col;
Chris@49:   uword cur_row = internal_row;
Chris@49: 
Chris@49:   const uword aux_col = iterator_base::M.aux_col1;
Chris@49:   const uword aux_row = iterator_base::M.aux_row1;
Chris@49:   const uword ln_cols = iterator_base::M.n_cols;
Chris@49: 
Chris@49:   while(true)
Chris@49:     {
Chris@49:     // Decrement the current column and see if we are on a new row.
Chris@49:     if(--cur_col > ln_cols)
Chris@49:       {
Chris@49:       cur_col = ln_cols - 1;
Chris@49:       cur_row--;
Chris@49:       }
Chris@49: 
Chris@49:     // Is there anything in this new column?
Chris@49:     const uword colptr = iterator_base::M.m.col_ptrs[cur_col + aux_col];
Chris@49:     const uword next_colptr = iterator_base::M.m.col_ptrs[cur_col + aux_col + 1];
Chris@49: 
Chris@49:     for(uword ind = colptr; (ind < next_colptr) && (iterator_base::M.m.row_indices[ind] <= (cur_row + aux_row)); ++ind)
Chris@49:       {
Chris@49:       const uword row_index = iterator_base::M.m.row_indices[ind];
Chris@49: 
Chris@49:       if((row_index - aux_row) == cur_row)
Chris@49:         {
Chris@49:         iterator_base::internal_col = cur_col;
Chris@49:         internal_row = cur_row;
Chris@49:         actual_pos = ind;
Chris@49: 
Chris@49:         return *this;
Chris@49:         }
Chris@49:       }
Chris@49:     }
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: typename SpSubview<eT>::const_row_iterator
Chris@49: SpSubview<eT>::const_row_iterator::operator--(int)
Chris@49:   {
Chris@49:   typename SpSubview<eT>::const_row_iterator tmp(*this);
Chris@49: 
Chris@49:   --(*this);
Chris@49: 
Chris@49:   return tmp;
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: bool
Chris@49: SpSubview<eT>::const_row_iterator::operator==(const const_iterator& rhs) const
Chris@49:   {
Chris@49:   return (rhs.row() == row()) && (rhs.col() == iterator_base::internal_col);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: bool
Chris@49: SpSubview<eT>::const_row_iterator::operator!=(const const_iterator& rhs) const
Chris@49:   {
Chris@49:   return (rhs.row() != row()) || (rhs.col() != iterator_base::internal_col);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: bool
Chris@49: SpSubview<eT>::const_row_iterator::operator==(const typename SpMat<eT>::const_iterator& rhs) const
Chris@49:   {
Chris@49:   return (rhs.row() == row()) && (rhs.col() == iterator_base::internal_col);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: bool
Chris@49: SpSubview<eT>::const_row_iterator::operator!=(const typename SpMat<eT>::const_iterator& rhs) const
Chris@49:   {
Chris@49:   return (rhs.row() != row()) || (rhs.col() != iterator_base::internal_col);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: bool
Chris@49: SpSubview<eT>::const_row_iterator::operator==(const const_row_iterator& rhs) const
Chris@49:   {
Chris@49:   return (rhs.row() == row()) && (rhs.col() == iterator_base::internal_col);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: bool
Chris@49: SpSubview<eT>::const_row_iterator::operator!=(const const_row_iterator& rhs) const
Chris@49:   {
Chris@49:   return (rhs.row() != row()) || (rhs.col() != iterator_base::internal_col);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: bool
Chris@49: SpSubview<eT>::const_row_iterator::operator==(const typename SpMat<eT>::const_row_iterator& rhs) const
Chris@49:   {
Chris@49:   return (rhs.row() == row()) && (rhs.col() == iterator_base::internal_col);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: bool
Chris@49: SpSubview<eT>::const_row_iterator::operator!=(const typename SpMat<eT>::const_row_iterator& rhs) const
Chris@49:   {
Chris@49:   return (rhs.row() != row()) || (rhs.col() != iterator_base::internal_col);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: ///////////////////////////////////////////////////////////////////////////////
Chris@49: // SpSubview<eT>::row_iterator implementation                                //
Chris@49: ///////////////////////////////////////////////////////////////////////////////
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: SpValProxy<SpSubview<eT> >
Chris@49: SpSubview<eT>::row_iterator::operator*()
Chris@49:   {
Chris@49:   return SpValProxy<SpSubview<eT> >(
Chris@49:     const_row_iterator::internal_row,
Chris@49:     iterator_base::internal_col,
Chris@49:     access::rw(iterator_base::M),
Chris@49:     &access::rw(iterator_base::M.m.values[const_row_iterator::actual_pos]));
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: typename SpSubview<eT>::row_iterator&
Chris@49: SpSubview<eT>::row_iterator::operator++()
Chris@49:   {
Chris@49:   const_row_iterator::operator++();
Chris@49:   return *this;
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: typename SpSubview<eT>::row_iterator
Chris@49: SpSubview<eT>::row_iterator::operator++(int)
Chris@49:   {
Chris@49:   typename SpSubview<eT>::row_iterator tmp(*this);
Chris@49: 
Chris@49:   ++(*this);
Chris@49: 
Chris@49:   return tmp;
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: typename SpSubview<eT>::row_iterator&
Chris@49: SpSubview<eT>::row_iterator::operator--()
Chris@49:   {
Chris@49:   const_row_iterator::operator--();
Chris@49:   return *this;
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename eT>
Chris@49: inline
Chris@49: typename SpSubview<eT>::row_iterator
Chris@49: SpSubview<eT>::row_iterator::operator--(int)
Chris@49:   {
Chris@49:   typename SpSubview<eT>::row_iterator tmp(*this);
Chris@49: 
Chris@49:   --(*this);
Chris@49: 
Chris@49:   return tmp;
Chris@49:   }
Chris@49: 
Chris@49: //! @}