---

SparseRowMatrix.h

/****************************************************************************
 *                                                                          *
 *  Program: FDDlib                                                         *
 *  Version: 1.1                                                            *
 *                                                                          *
 *  Copyright (C) 2002 - 2004 by Eric Miller and Dana Brooks                *
 *  All rights reserved.                                                    *
 *                                                                          *
 *  This software is Version 1.1 of the fddlib tomography toolbox.          *
 *  It is not to be redistributed or used for any commercial purpose        *
 *  without the prior written consent of the authors and Northeastern       *
 *  University.                                                             *
 *                                                                          *
 *  This software is provided as is, and any express or implied warranties, *
 *  including but not limited to the implied warranty of merchantability    *
 *  and the implied warranty of fitness for a particular purpose, are dis-  *
 *  claimed.  In no event shall the authors or Northeastern University be   *
 *  liable for any direct, indirect, incidental, special, exemplary, or     *
 *  consequential damages (including but not limited to procurement of      * 
 *  substitute goods or services; loss of use, data, or profits; or busi-   *
 *  ness interruption) however caused and on any theory of liability,       *
 *  whether in contract, strict liability, or tort (including negligence or *
 *  otherwise) arising in any way from the use of this software, even if    *
 *  advised of the possibility of such damage.                              *
 *                                                                          *
 *  Portions of this code benefit from ideas from Kyle Guilbert, Greg       *
 *  Boverman, Derek Uluski, David Kaeli, and Jennifer Black                 *
 *                                                                          *
 ****************************************************************************/

#ifndef _SparseRowMatrix_H_
#define _SparseRowMatrix_H_

#include <stdlib.h>
#include <stdio.h>
#include <string>
#include "Matrix.h"
#include "RealVector.h"

namespace FDDlib {

template <class T>
struct SparseRowEntry
{
  int column;
  T   val;
};

template<class T>
class SparseRowMatrix : public Matrix<T>
{
protected:
  int              MaxRows_;
  int              EntriesPerRow_;
  int*             NumCols_;
  SparseRowEntry<T>*        Entries_;

public:
  SparseRowMatrix();
  
  SparseRowMatrix(int maxrows, int max_per_row);
  
  SparseRowMatrix(const SparseRowMatrix<T>& SRM);
  
  virtual ~SparseRowMatrix();

  virtual SparseRowMatrix<T>& operator=(const SparseRowMatrix<T>& SRM);

  virtual void init(int maxrows, int max_per_row);

  virtual void deallocate();

  virtual T val(int row, int col) const;

  virtual int getEntriesPerRow() const;

  virtual SparseRowEntry<T>* getRowEntry(int row) const;

  virtual int getCol(int row, int n) const;

  virtual int getCol(int row, int n, T& v) const;

  virtual int getNumCols(int row) const;

  virtual int getNumCols() const;

  virtual int getNumRows() const;

  virtual void set(int row, int col, T val) throw(std::string);

  template<class Vector>
  Vector operator*(const Vector& m) const;
  
  RealVector<T> getRow(int row) const;

  RealVector<T> getColumn(int col) const;

};


template <class T>
SparseRowMatrix<T>::SparseRowMatrix()
{
  MaxRows_ = 0;
  EntriesPerRow_ = 0;
  NumCols_ = (int*) NULL;
  Entries_ = (SparseRowEntry<T>*) NULL;
}

template <class T>
SparseRowMatrix<T>::SparseRowMatrix(int maxrows,
                                            int max_per_row)
{
  MaxRows_ = maxrows;
  EntriesPerRow_ = max_per_row;
  NumCols_ = new int[MaxRows_];
  Entries_ = new SparseRowEntry<T>[MaxRows_*EntriesPerRow_];
  for (int i = 0; i < MaxRows_; i++)
  {
    NumCols_[i] = 0;
  }
}

template<class T>
SparseRowMatrix<T>::SparseRowMatrix(const SparseRowMatrix<T>& ESR)
{
  int       i;
  SparseRowEntry<T>* s, *d;

  MaxRows_ =  ESR.getNumRows();
  EntriesPerRow_ = ESR.getEntriesPerRow();
  NumCols_ = new int[MaxRows_];
  Entries_ = new SparseRowEntry<T>[MaxRows_*EntriesPerRow_];
  for (i = 0; i < MaxRows_; i++)
  {
    s = ESR.getRowEntry(i);
    d = getRowEntry(i);
    NumCols_[i] = ESR.getNumCols(i);
    for (j = 0; j < EntriesPerRow_; j++)
    {
      d[j].column = s[j].column;
      d[j].val = s[j].val;
    }
  }
}

template<class T>
SparseRowMatrix<T>& SparseRowMatrix<T>::operator=(const SparseRowMatrix<T>& ESR)
{
  int       i, j;
  SparseRowEntry<T>* s, *d;

  if (getNumRows() > 0)
  {
    delete[] NumCols_;
    delete[] Entries_;
  }
  
  MaxRows_ = ESR.getNumRows();
  EntriesPerRow_ = ESR.getEntriesPerRow();

  NumCols_ = new int[ESR.getNumRows()];
  Entries_ = new SparseRowEntry<T>[MaxRows_*EntriesPerRow_];
  for (i = 0; i < MaxRows_; i++)
  {
    s = ESR.getRowEntry(i);
    d = getRowEntry(i);
    NumCols_[i] = ESR.getNumCols(i);
    for (j = 0; j < NumCols_[i]; j++)
    {
      d[j].column = s[j].column;
      d[j].val = s[j].val;
    }
  }

  return *this;
}

template <class T>
void SparseRowMatrix<T>::init(int maxrows, int max_per_row)
{
  MaxRows_ = maxrows;
  EntriesPerRow_ = max_per_row;

  if (NumCols_ != (int*) NULL)
    delete[] NumCols_;

  if (Entries_ != (SparseRowEntry<T>*) NULL)
      delete[] Entries_;

  NumCols_ = new int[MaxRows_];
  Entries_ = new SparseRowEntry<T>[MaxRows_*EntriesPerRow_];
  for (int i = 0; i < MaxRows_; i++)
  {
    NumCols_[i] = 0;
  }
}

template<class T>
SparseRowMatrix<T>::~SparseRowMatrix()
{
  if (Entries_ != (SparseRowEntry<T>*) NULL)
    delete[] Entries_;

  if (NumCols_ != (int*) NULL)
    delete[] NumCols_;
}

template<class T>
void SparseRowMatrix<T>::deallocate()
{
  if (Entries_ != (SparseRowEntry<T>*) NULL)
    delete[] Entries_;

  if (NumCols_ != (int*) NULL)
    delete[] NumCols_;

  Entries_ = (SparseRowEntry<T>*) NULL;
  NumCols_ = (int*) NULL;
  MaxRows_ = 0;
  EntriesPerRow_ = 0;
}


template<class T>
inline T SparseRowMatrix<T>::val(int row, int col) const
{
  int i;

  for (i = 0; i < NumCols_[row]; i++)
  {
    if (Entries_[row*EntriesPerRow_ + i].column == col)
      return Entries_[row*EntriesPerRow_ + i].val;
  }
  return 0;
}

template<class T>
inline int SparseRowMatrix<T>::getEntriesPerRow() const
{
  return EntriesPerRow_;
}


template<class T>
inline SparseRowEntry<T>* SparseRowMatrix<T>::getRowEntry(int r) const
{
  return (Entries_ + r*EntriesPerRow_);
}


template<class T>
inline int SparseRowMatrix<T>::getCol(int row, int n) const
{
  return (Entries_[row*EntriesPerRow_ + n].column);
}


template<class T>
inline int SparseRowMatrix<T>::getCol(int row, int n, T& v) const
{
  v = Entries_[row*EntriesPerRow_ + n].val;
  return (Entries_[row*EntriesPerRow_ + n].column);
}


template<class T>
inline int SparseRowMatrix<T>::getNumCols(int row) const
{
  // if we used the default constructor NumCols will be null.
  if (NumCols_ == (int*)NULL)
    return 0;
    
  return  NumCols_[row];
}

template<class T>
inline int SparseRowMatrix<T>::getNumCols() const
{
  int i, j;

  int col;

  // Find the maximum column index of all the rows
  int maxcol = 0;

  for (i = 0; i < getNumRows(); i++)
  {
    for (j = 0; j < getNumCols(i); j++)
    {
      col = getCol(i, j);
      if (col > maxcol)
      maxcol = col;
    }
  }

  return maxcol+1;
}

template<class T>
inline int SparseRowMatrix<T>::getNumRows() const
{
  return MaxRows_;
}


template<class T>
inline void SparseRowMatrix<T>::set(int row, int col, T val) throw(std::string)
{
  int       i;
  SparseRowEntry<T>* ep;
  
  // this will likely prevent a segfault
  if (row < 0 || col < 0)
    throw std::string("SRM::set: negative row or column parameter");
  
  ep = getRowEntry(row);
  for (i = 0; i < NumCols_[row]; i++)
  {
    if (ep->column == col)
      break;
    ep++;
  }
  if (ep - getRowEntry(row) >= EntriesPerRow_)
    throw std::string("SRM::set: hit max number of entries per row");
  if (ep - getRowEntry(row) == NumCols_[row])
    (NumCols_[row])++;

  ep->val = val;
  ep->column = col;
}


template<class T>
template<class Vector>
inline Vector SparseRowMatrix<T>::operator*(const Vector& m) const
{
  int       i;
  int       j;
  T         rv;
  int       cols;
  SparseRowEntry<T>* ep;

  Vector ret(MaxRows_);
  for (i = 0; i < MaxRows_; i++)
  {
    rv = (T) 0;
    cols = NumCols_[i];
    ep = getRowEntry(i);
    for (j = 0; j < cols; j++)
    {
      rv += m(ep->column) * (ep->val);
      ep++;
    }
    ret(i) = rv;
  }

  return ret;
}

template<class T>
RealVector<T> SparseRowMatrix<T>::getRow(int row) const
{
  int i;
  const int numcols = getNumCols();
  RealVector<T> ret(numcols);
  
  for(i=0;i<numcols;i++)
    ret.set(i,val(row,i));
  
  return ret;
}

template<class T>
RealVector<T> SparseRowMatrix<T>::getColumn(int col) const
{
  int i;
  RealVector<T> ret(MaxRows_);
  
  for(i=0;i<MaxRows_;i++)
    ret.set(i,val(i,col));
    
  return ret;
}

} // end namespace
#endif

---