matrix.h

Go to the documentation of this file.

/*
 *   Copyright (c) 2007 John Weaver
 *
 *   This program is free software; you can redistribute it and/or modify
 *   it under the terms of the GNU General Public License as published by
 *   the Free Software Foundation; either version 2 of the License, or
 *   (at your option) any later version.
 *
 *   This program is distributed in the hope that it will be useful,
 *   but WITHOUT ANY WARRANTY; without even the implied warranty of
 *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *   GNU General Public License for more details.
 *
 *   You should have received a copy of the GNU General Public License
 *   along with this program; if not, write to the Free Software
 *   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
 */

#ifndef _MATRIX_H_
#define _MATRIX_H_

#include <cassert>
#include <cstdlib>
#include <algorithm>

template <class T>
class Matrix {
        public:
                Matrix()
                {
                        m_rows = 0;
                        m_columns = 0;
                        m_matrix = NULL;
                }
                
                Matrix(int rows, int columns)
                {
                        m_matrix = NULL;
                        resize(rows, columns);
                }
                
                Matrix(const Matrix<T> &other)
                {
                        if ( other.m_matrix != NULL ) {
                                // copy arrays
                                m_matrix = NULL;
                                resize(other.m_rows, other.m_columns);
                                for ( int i = 0 ; i < m_rows ; i++ )
                                        for ( int j = 0 ; j < m_columns ; j++ )
                                                m_matrix[i][j] = other.m_matrix[i][j];
                        } else {
                                m_matrix = NULL;
                                m_rows = 0;
                                m_columns = 0;
                        }
                }
                
                Matrix<T> & operator= (const Matrix<T> &other)
                {
                        if ( other.m_matrix != NULL ) {
                                // copy arrays
                                resize(other.m_rows, other.m_columns);
                                for ( int i = 0 ; i < m_rows ; i++ )
                                        for ( int j = 0 ; j < m_columns ; j++ )
                                                m_matrix[i][j] = other.m_matrix[i][j];
                        } else {
                                // free arrays
                                for ( int i = 0 ; i < m_columns ; i++ )
                                        delete [] m_matrix[i];

                                delete [] m_matrix;

                                m_matrix = NULL;
                                m_rows = 0;
                                m_columns = 0;
                        }
                        
                        return *this;
                }
                
                ~Matrix()
                {
                        if ( m_matrix != NULL ) {
                                // free arrays
                                for ( int i = 0 ; i < m_rows ; i++ )
                                        delete [] m_matrix[i];

                                delete [] m_matrix;
                        }
                        m_matrix = NULL;
                }
                
                // all operations except product modify the matrix in-place.
                void resize(int rows, int columns)
                {
                        if ( m_matrix == NULL ) {
                                // alloc arrays
                                m_matrix = new T*[rows]; // rows
                                for ( int i = 0 ; i < rows ; i++ )
                                        m_matrix[i] = new T[columns]; // columns

                                m_rows = rows;
                                m_columns = columns;
                                clear();
                        } else {
                                // save array pointer
                                T **new_matrix;
                                // alloc new arrays
                                new_matrix = new T*[rows]; // rows
                                for ( int i = 0 ; i < rows ; i++ ) {
                                        new_matrix[i] = new T[columns]; // columns
                                        for ( int j = 0 ; j < columns ; j++ )
                                                new_matrix[i][j] = 0;
                                }

                                // copy data from saved pointer to new arrays
                                int minrows = std::min<int>(rows, m_rows);
                                int mincols = std::min<int>(columns, m_columns);
                                for ( int x = 0 ; x < minrows ; x++ )
                                        for ( int y = 0 ; y < mincols ; y++ )
                                                new_matrix[x][y] = m_matrix[x][y];

                                // delete old arrays
                                if ( m_matrix != NULL ) {
                                        for ( int i = 0 ; i < m_rows ; i++ )
                                                delete [] m_matrix[i];

                                        delete [] m_matrix;
                                }

                                m_matrix = new_matrix;
                        }

                        m_rows = rows;
                        m_columns = columns;
                }
                
                void identity(void)
                {
                        assert( m_matrix != NULL );

                        clear();

                        int x = std::min<int>(m_rows, m_columns);
                        for ( int i = 0 ; i < x ; i++ )
                                m_matrix[i][i] = 1;
                }
                
                void clear(void)
                {
                        assert( m_matrix != NULL );

                        for ( int i = 0 ; i < m_rows ; i++ )
                                for ( int j = 0 ; j < m_columns ; j++ )
                                        m_matrix[i][j] = 0;
                }
                
                T& operator () (int x, int y)
                {
                        assert ( x >= 0 );
                        assert ( y >= 0 );
                        assert ( x < m_rows );
                        assert ( y < m_columns );
                        assert ( m_matrix != NULL );
                        return m_matrix[x][y];
                }
                
                T trace(void)
                {
                        assert( m_matrix != NULL );

                        T value = 0;

                        int x = std::min<int>(m_rows, m_columns);
                        for ( int i = 0 ; i < x ; i++ )
                                value += m_matrix[i][i];

                        return value;
                }
                
                
                Matrix<T>& transpose(void)
                {
                        assert( m_rows > 0 );
                        assert( m_columns > 0 );

                        int new_rows = m_columns;
                        int new_columns = m_rows;

                        if ( m_rows != m_columns ) {
                                // expand matrix
                                int m = std::max<int>(m_rows, m_columns);
                                resize(m,m);
                        }

                        for ( int i = 0 ; i < m_rows ; i++ ) {
                                for ( int j = i+1 ; j < m_columns ; j++ ) {
                                        T tmp = m_matrix[i][j];
                                        m_matrix[i][j] = m_matrix[j][i];
                                        m_matrix[j][i] = tmp;
                                }
                        }

                        if ( new_columns != new_rows ) {
                                // trim off excess.
                                resize(new_rows, new_columns);
                        }

                        return *this;
                }
                
                Matrix<T> product(Matrix<T> &other)
                {
                        assert( m_matrix != NULL );
                        assert( other.m_matrix != NULL );
                        assert ( m_columns == other.m_rows );

                        Matrix<T> out(m_rows, other.m_columns);

                        for ( int i = 0 ; i < out.m_rows ; i++ ) {
                                for ( int j = 0 ; j < out.m_columns ; j++ ) {
                                        for ( int x = 0 ; x < m_columns ; x++ ) {
                                                out(i,j) += m_matrix[i][x] * other.m_matrix[x][j];
                                        }
                                }
                        }

                        return out;
                }
                
                int minsize(void)
                {
                        return ((m_rows < m_columns) ? m_rows : m_columns);
                }
                
                int columns(void)
                {
                        return m_columns;
                }
                
                int rows(void)
                {
                        return m_rows;
                }
        private:
                T **m_matrix;
                int m_rows;
                int m_columns;
};

// #ifndef USE_EXPORT_KEYWORD
// #include "matrix.cpp"
//#define export /*export*/
// #endif

#endif /* !defined(_MATRIX_H_) */