source: trunk/src/matrix.h @ 40

#ifndef _thep_gsl_api_matrix_ 
#define _thep_gsl_api_matrix_

#include <iostream>
#include <fstream>
#include <iomanip>
#include <cmath>
#include <cstdlib>
#include <cassert>

#include <gsl/gsl_math.h>
#include <gsl/gsl_matrix.h>
#include <gsl/gsl_linalg.h>
#include <gsl/gsl_blas.h>

#include "vector.h"

namespace thep_gsl_api
{

  /**
     This is a matrix interface to GSL for c++.
     Note that the namespace is "thep_gsl_api".
  */

  class matrix
  {
  public:
    typedef double Type;
  
    /**
       Default constructor
       Does not allocate any memory
    */
    matrix();

    /**
       Constructor taking three argumens:
       N = Number of rows
       M = Number of columns
       init_to_zero, if true matrix will be initialized with zeros.
    */
    matrix( const size_t& rows, const size_t& cols, 
      bool init_to_zero = true );




 
    /** 
  Copy constructor
    */
    matrix( const matrix& other );


    /**
       Constructor taking as only argument a matrix of gsl type.
       Be causious here! Do not delete this matrix outside of
       class: NO COPY IS MADE!!!
    */
    matrix( gsl_matrix* );
    

    /**
       Constructor taking as only argument a an istream that should
       contain the matrix where columns are separated with whitespace
       (not '\\n') and rows by '\\n'.
    */
    matrix(std::istream &);

    /**
       Destructor will clear matrix's allocated memory
    */
    ~matrix();


    /**
       The most important function when adding new classes to 
       c++_tools. This function will return a pointer to the 
       inner data of the class. Should be forbidden in a strictly 
       object-orientated environment.
       Motivation: Effeciency! ( CAN IT BE MADE CONST? )
       To be used when writting new classes using gsl-functions. 
    */
    inline gsl_matrix* get_gsl_matrix() const { return m_; }  


    /**
       Returns the number of rows in the matrix.
    */
    size_t rows() const { return m_->size1; }


    /**
       Returns the number of rows in the matrix.
    */
    size_t cols() const { return m_->size2; }

    /**
       get( i, j ) will return the Aij element in matrix A
       where i = row number and j = column number
    */
    inline Type get( const size_t& i, const size_t& j ) const
      { return gsl_matrix_get(m_,i_row,j_col); }


    /**
       set( i, j, val ) will make Aij = val in matrix A
       where i = row number and j = column number
    */
    inline void set( const size_t& i, const size_t& j, const Type& val )
      { gsl_matrix_set( m_, i_row, j_col, val ); }

    /**
       set_all( val ) will make Aij = val in matrix A
       for all "i"
    */
    inline void set_all( const Type& val ) { gsl_matrix_set_all(m_,val); }



    /**
       Assignment-operator. No requirements on dimensions,
       i.e. if A = B, then A will be a copy of B with B:s
       dimensions.
     */
    matrix& operator=( const matrix& other );
    
    /**
       operators for comparison are:
       == equal
       != not equal
    */
    bool operator==( const matrix &other ) const;
    bool operator!=( const matrix &other ) const; 

    /**
       operators for basic arithmetics are: \n
       \f$+\f$ matrix-matrix addition \n
       \f$-\f$ matrix-matrix subtraction \n 
       \f$*\f$ matrix-matrix multiplication
    */
    matrix operator+( const matrix &other ) const;
    matrix operator-( const matrix &other ) const;
    matrix operator*( const matrix &other ) const;
    vector operator*( const vector &other ) const;

    /** 
  multiplying each element of two matrices
  dimension MUST AGREE
    */
    matrix mul_elements( const matrix& a, const matrix& b );

    /**
       Method for transposing a matrix
     */
    matrix transpose() const;
    
    /**
       Method for calculating the n:th sqrt of the sum of 
       the n:th power of all elements
    */
    double norm( double n ) const;

    /**
       Calculates sum of all elements in matrix
    */    
    double matrix::sum() const;
      
    /**
       Exchange row i with row j (vice versa)
     */
    void swap_rows( const size_t& i, const size_t& j );

    /**
       Exchange column i with column j (vice versa)
     */
    void swap_cols( const size_t& i, const size_t& j );

    /**
       Returns row i in matrix
    */    
    matrix row( const size_t& i ) const;


    /**
       Returns row i in matrix (vector form)
    */
    vector row_vector( const size_t& i ) const;


    /**
       Returns col i in matrix
    */    
    matrix col( const size_t& i ) const;


    /**
       Returns col i in matrix (vector form)
    */
    vector col_vector( const size_t& i ) const;


    /**
       Method for calculating the row sum
    */    
    vector row_sum() const;

    /**
       Method for calculating the mean row sum
    */    
    vector mean_row_sum() const;

    /** 
  standard output operator is defined
    */
    friend std::ostream& operator<< ( std::ostream& s_out, const matrix& );

    /**
       Getting vector operator
    */
    vector operator[]( const size_t& i ) const;
  
  private:
    gsl_matrix* new_copy( const gsl_matrix* );

    gsl_matrix* m_;
  };


};




#endif