Matrix.h

Go to the documentation of this file.

//Copyright (c) 2004-2005, Baris Sumengen
//All rights reserved.
//
// CIMPL Matrix Performance Library
//
//Redistribution and use in source and binary
//forms, with or without modification, are
//permitted provided that the following
//conditions are met:
//
//    * No commercial use is allowed. 
//    This software can only be used
//    for non-commercial purposes. This 
//    distribution is mainly intended for
//    academic research and teaching.
//    * Redistributions of source code must
//    retain the above copyright notice, this
//    list of conditions and the following
//    disclaimer.
//    * Redistributions of binary form must
//    mention the above copyright notice, this
//    list of conditions and the following
//    disclaimer in a clearly visible part 
//    in associated product manual, 
//    readme, and web site of the redistributed 
//    software.
//    * Redistributions in binary form must
//    reproduce the above copyright notice,
//    this list of conditions and the
//    following disclaimer in the
//    documentation and/or other materials
//    provided with the distribution.
//    * The name of Baris Sumengen may not be
//    used to endorse or promote products
//    derived from this software without
//    specific prior written permission.
//
//THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT
//HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
//EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT
//NOT LIMITED TO, THE IMPLIED WARRANTIES OF
//MERCHANTABILITY AND FITNESS FOR A PARTICULAR
//PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
//CONTRIBUTORS 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 BUSINESS INTERRUPTION)
//HOWEVER CAUSED AND ON ANY THEORY OF
//LIABILITY, WHETHER IN CONTRACT, STRICT
//LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR
//OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
//OF THIS SOFTWARE, EVEN IF ADVISED OF THE
//POSSIBILITY OF SUCH DAMAGE.



#pragma once
#ifndef MATRIX_H
#define MATRIX_H



#include <iostream>

using std::cout;
using std::cerr;
using std::endl;
using std::ostream;
using std::right;
using std::fixed;

#include <iomanip>

using std::setw;

#include <cstdlib>
#include <ctime>
#include <math.h>

#include <string>

using std::string;

#include <sstream>

using std::ostringstream;

#include <typeinfo>


#include "./Cleaner.h"
#include "./Utility.h"
#include "./Array.h"
#include "./Vector.h"
#include "./SubMatrix.h"
#include "./RandomGen.h"


// forward declaration
template< class T > class Array;
template< class T > class Matrix;
template< class T > class SubMatrix;
template< class T > class Vector;
template< class T, int l > class List;
template< class T > class VarList;

template< class T >
class Matrix 
{

    friend class Array<T>;
    friend class Vector<T>;
    friend class SubMatrix<T>;


protected:
    T *data;
    int ndims;  // always 2
    int length; // xDim*yDim
    int xDim;
    int yDim;
    Cleaner<T> *clean;   // Reference counting Garbage collector.

    Vector<T> *columns;
    
public:
    Matrix();
    Matrix(const int rows, const int columns);
    Matrix(const int rows, const int columns, T init);
    

    Matrix(Matrix<T> &m); // copy constructor

    ~Matrix();
    
    void Clean();

    void Set(const int rows, const int columns, T *_data); 

    const T* DataPtr();
    T* Data();

    Matrix<T> Clone() const;
    
    void SwitchColumns(int i, int j);

    void SyncData2Columns();
    
    SubMatrix<T> Slice(int row1, int row2, int col1, int col2);


    const int Columns() const;
    const int Rows() const;
    const int XDim() const;
    const int YDim() const;
    const int* Dims() const;
    const List<int,2> Size() const;
    const int Length() const;
    const int NDims() const;
    void Init(const T init);
    Matrix<T>&  Rand();
    Matrix<T>&  Rand(const double max);
    static Matrix<T> Rand(const int rows, const int cols);
    static Matrix<T> Rand(const int rows, const int cols, const double max);

    void ReadFromMatrix(const Matrix<T>& m);
    void ReadFromMatrix(const Matrix<T>& m, const int rowStart, const int colStart);


    static bool IsSquare(Matrix<T>& m);
    static bool IsM2MCompatible(Matrix<T>& m1, Matrix<T>& m2);
    static Matrix<T> MMultiply(Matrix<T>& m1, Matrix<T>& m2);
    static Matrix<T> MMultiply(Matrix<T>& m1, Vector<T>& m2);
    static Matrix<T> MMultiply(Vector<T>& m1, Matrix<T>& m2);

    static Matrix<T> Transpose(Matrix<T>& m);
    Matrix<T>& Transpose();

    //static Matrix<T> Hermitian(Matrix<T>& m);
    //Matrix<T>& Hermitian();

    static bool IsCompatible(Matrix<T>& m1, Matrix<T>& m2);
    
    // Boolean Operations...
    static Matrix<int> And(Matrix<T>& m1, Matrix<T>& m2);
    static Matrix<int> Or(Matrix<T>& m1, Matrix<T>& m2);
    static Matrix<int> Lt(Matrix<T>& m1, Matrix<T>& m2);
    static Matrix<int> Gt(Matrix<T>& m1, Matrix<T>& m2);
    static Matrix<int> Le(Matrix<T>& m1, Matrix<T>& m2);
    static Matrix<int> Ge(Matrix<T>& m1, Matrix<T>& m2);
    static Matrix<int> Eq(Matrix<T>& m1, Matrix<T>& m2);
    static Matrix<int> Ne(Matrix<T>& m1, Matrix<T>& m2);

    // Boolean Operations with value type...
    static Matrix<int> And(Matrix<T>& m, T v);
    static Matrix<int> Or(Matrix<T>& m, T v);
    static Matrix<int> Lt(Matrix<T>& m, T v);
    static Matrix<int> Gt(Matrix<T>& m, T v);
    static Matrix<int> Le(Matrix<T>& m, T v);
    static Matrix<int> Ge(Matrix<T>& m, T v);
    static Matrix<int> Eq(Matrix<T>& m, T v);
    static Matrix<int> Ne(Matrix<T>& m, T v);


// Add matrix to another matrix
    static Matrix<T> Add(Matrix<T>& m1, Matrix<T>& m2);
    static Matrix<T> Subtract(Matrix<T>& m1, Matrix<T>& m2);
    static Matrix<T> Multiply(Matrix<T>& m1, Matrix<T>& m2);
    static Matrix<T> Divide(Matrix<T>& m1, Matrix<T>& m2);

// Add value to another matrix
    static Matrix<T> Add(Matrix<T>& m1, T v2);
    static Matrix<T> Subtract(Matrix<T>& m1, T v2);
    static Matrix<T> Subtract(T v2, Matrix<T>& m1);
    static Matrix<T> Multiply(Matrix<T>& m1, T v2);
    static Matrix<T> Divide(Matrix<T>& m1, T v2);
    static Matrix<T> Divide(T v2, Matrix<T>& m1);


// Add matrix to "this" matrix
    Matrix<T>& Add(Matrix<T>& m);
    Matrix<T>& Subtract(Matrix<T>& m);
    Matrix<T>& Multiply(Matrix<T>& m);
    Matrix<T>& Divide(Matrix<T>& m);

// Add value to "this" matrix
    Matrix<T>& Add(T v);
    Matrix<T>& Subtract(T v);
    Matrix<T>& Multiply(T v);
    Matrix<T>& Divide(T v);


//OPERATORS
    Matrix<T>& operator= (Matrix<T>& m);
    Matrix<T>& operator= (Array<T>& m);
    Matrix<T>& operator= (Vector<T>& m);


    Matrix<T> operator+ ();
    Matrix<T> operator- ();
    Matrix<int> operator! ();

    // Transpose
    Matrix<T> operator~ ();

    friend ostream& operator<< (ostream& output, Matrix<T>& m)
    {
        int bufferWidth = 80;
        int lm = m.Rows();
        int rowNoWidth = (int)log10((double)(lm-1))+2;
        int infoWidth = 3+rowNoWidth+2;

        int maxLength  = 1;
        for(int i=0; i<m.Length(); i++)
        {
            ostringstream oS;
            oS << m.ElemNC(i);
            int l = (int)oS.str().length();
            if(l>maxLength)
            {
                maxLength = l;
            }
        }

        int columnWidth = maxLength+3;
        int columnMax = (bufferWidth-infoWidth)/columnWidth;
        
        output << typeid(m).name() << " of size " << m.Rows() << "x" << m.Columns() << endl;
        output << "=======================" << endl;
        if(m.Columns() <= columnMax)
        {
            for(int i=0;i<m.Rows();i++)
            {
                output << "ROW" << right << setw(rowNoWidth) << i+1 << "|";
                for(int j=0;j<m.Columns();j++)
                {
                    output << right << setw(columnWidth) << m.ElemNC(i,j);
                }
                output << endl;
            }
        }
        else
        {
            int groups = m.Columns()/columnMax;
            int leftover = m.Columns() % columnMax;
            for(int k=0; k<groups; k++)
            {
                output << endl;
                output << "----------------------" << endl;
                output << "Columns " << k*columnMax+1 << " through " << k*columnMax+columnMax-1+1 << endl;
                output << "----------------------" << endl;
                for(int i=0;i<m.Rows();i++)
                {
                    output << "ROW" << right << setw(rowNoWidth) << i+1 << "|";
                    for(int j=columnMax*k;j<columnMax*(k+1);j++)
                    {
                        output << right << setw(columnWidth) << m.ElemNC(i,j);
                    }
                    output << endl;
                }
            }
            if(leftover > 0)
            {
                output << endl;
                output << "----------------------" << endl;
                output << "Columns " << groups*columnMax+1 << " through " << groups*columnMax+leftover-1+1 << endl;
                output << "----------------------" << endl;
                for(int i=0;i<m.Rows();i++)
                {
                    output << "ROW" << right << setw(rowNoWidth) << i+1 << "|";
                    for(int j=groups*columnMax;j<groups*columnMax+leftover;j++)
                    {
                        output << right << setw(columnWidth) << m.ElemNC(i,j);
                    }
                    output << endl;
                }
            }
        }
        output << "=======================" << endl;
        output << endl;
        
        return output;
    }


    T& Coor(const int x, const int y);
    T& CoorNC(const int x, const int y); // No bounds check
    T& Elem(const int _row, const int _col);
    T& ElemNC(const int _row, const int _col); // No bounds check
    T& Elem(const int i);
    T& ElemNC(const int i); // No bounds check

    T& operator() (const int i);
    T& operator() (const int i, const int j);
    Vector<T>& operator[] (const int i);
    
    Matrix<T> operator() (const int row1, const int row2, const int col1, const int col2);


    friend Matrix<T> operator& (Matrix<T>& m1, Matrix<T>& m2)
    {
        return Matrix<T>::MMultiply(m1, m2);
    }
    
    friend Matrix<T> operator& (Matrix<T>& m1, Vector<T>& m2)
    {
        return Matrix<T>::MMultiply(m1, m2);
    }

    friend Matrix<T> operator& (Vector<T>& m1, Matrix<T>& m2)
    {
        return Matrix<T>::MMultiply(m1, m2);
    }


    // Boolean operations
    friend Matrix<int> operator&& (Matrix<T>& m1, Matrix<T>& m2)
    {
        return Matrix<T>::And(m1, m2);
    }

    friend Matrix<int> operator|| (Matrix<T>& m1, Matrix<T>& m2)
    {
        return Matrix<T>::Or(m1, m2);
    }

    friend Matrix<int> operator< (Matrix<T>& m1, Matrix<T>& m2)
    {
        return Matrix<T>::Lt(m1, m2);
    }

    friend Matrix<int> operator> (Matrix<T>& m1, Matrix<T>& m2)
    {
        return Matrix<T>::Gt(m1, m2);
    }

    friend Matrix<int> operator<= (Matrix<T>& m1, Matrix<T>& m2)
    {
        return Matrix<T>::Le(m1, m2);
    }

    friend Matrix<int> operator>= (Matrix<T>& m1, Matrix<T>& m2)
    {
        return Matrix<T>::Ge(m1, m2);
    }

    friend Matrix<int> operator== (Matrix<T>& m1, Matrix<T>& m2)
    {
        return Matrix<T>::Eq(m1, m2);
    }

    friend Matrix<int> operator!= (Matrix<T>& m1, Matrix<T>& m2)
    {
        return Matrix<T>::Ne(m1, m2);
    }


    // Boolean operations with value type
    friend Matrix<int> operator&& (Matrix<T>& m, T v)
    {
        return Matrix<T>::And(m, v);
    }
    friend Matrix<int> operator&& (T v, Matrix<T>& m)
    {
        return Matrix<T>::And(m, v);
    }

    friend Matrix<int> operator|| (Matrix<T>& m, T v)
    {
        return Matrix<T>::Or(m, v);
    }
    friend Matrix<int> operator|| (T v, Matrix<T>& m)
    {
        return Matrix<T>::Or(m, v);
    }

    friend Matrix<int> operator< (Matrix<T>& m, T v)
    {
        return Matrix<T>::Lt(m, v);
    }
    friend Matrix<int> operator< (T v, Matrix<T>& m)
    {
        return Matrix<T>::Gt(m, v);
    }

    friend Matrix<int> operator> (Matrix<T>& m, T v)
    {
        return Matrix<T>::Gt(m, v);
    }
    friend Matrix<int> operator> (T v, Matrix<T>& m)
    {
        return Matrix<T>::Lt(m, v);
    }

    friend Matrix<int> operator<= (Matrix<T>& m, T v)
    {
        return Matrix<T>::Le(m, v);
    }
    friend Matrix<int> operator<= (T v, Matrix<T>& m)
    {
        return Matrix<T>::Ge(m, v);
    }

    friend Matrix<int> operator>= (Matrix<T>& m, T v)
    {
        return Matrix<T>::Ge(m, v);
    }
    friend Matrix<int> operator>= (T v, Matrix<T>& m)
    {
        return Matrix<T>::Le(m, v);
    }

    friend Matrix<int> operator== (Matrix<T>& m, T v)
    {
        return Matrix<T>::Eq(m, v);
    }
    friend Matrix<int> operator== (T v, Matrix<T>& m)
    {
        return Matrix<T>::Eq(m, v);
    }

    friend Matrix<int> operator!= (Matrix<T>& m, T v)
    {
        return Matrix<T>::Ne(m, v);
    }
    friend Matrix<int> operator!= (T v, Matrix<T>& m)
    {
        return Matrix<T>::Ne(m, v);
    }








// Add matrix to another matrix
    friend Matrix<T> operator+ (Matrix<T>& m1, Matrix<T>& m2)
    {
        return Matrix<T>::Add(m1, m2);
    }

    friend Matrix<T> operator- (Matrix<T>& m1, Matrix<T>& m2)
    {
        return Matrix<T>::Subtract(m1, m2);
    }

    friend Matrix<T> operator* (Matrix<T>& m1, Matrix<T>& m2)
    {
        return Matrix<T>::Multiply(m1, m2);
    }

    friend Matrix<T> operator/ (Matrix<T>& m1, Matrix<T>& m2)
    {
        return Matrix<T>::Divide(m1, m2);
    }


// Add value to another matrix
    friend Matrix<T> operator+ (Matrix<T>& m, T v)
    {
        return Matrix<T>::Add(m, v);
    }

    friend Matrix<T> operator+ (T v, Matrix<T>& m)
    {
        return Matrix<T>::Add(m, v);
    }

    friend Matrix<T> operator- (Matrix<T>& m, T v)
    {
        return Matrix<T>::Subtract(m, v);
    }

    friend Matrix<T> operator- (T v, Matrix<T>& m)
    {
        return Matrix<T>::Subtract(v, m);
    }

    friend Matrix<T> operator* (Matrix<T>& m, T v)
    {
        return Matrix<T>::Multiply(m, v);
    }

    friend Matrix<T> operator* (T v, Matrix<T>& m)
    {
        return Matrix<T>::Multiply(m, v);
    }

    friend Matrix<T> operator/ (Matrix<T>& m, T v)
    {
        return Matrix<T>::Divide(m, v);
    }

    friend Matrix<T> operator/ (T v, Matrix<T>& m)
    {
        return Matrix<T>::Divide(v, m);
    }




// Add inline
    Matrix<T>& operator+= (Matrix<T>& m);
    Matrix<T>& operator-= (Matrix<T>& m);
    Matrix<T>& operator*= (Matrix<T>& m);
    Matrix<T>& operator/= (Matrix<T>& m);

    Matrix<T>& operator+= (T v);
    Matrix<T>& operator-= (T v);
    Matrix<T>& operator*= (T v);
    Matrix<T>& operator/= (T v);

    
// TYPE CONVERSIONS
    Matrix(Array<T> &m);
    Matrix(Vector<T> &m);


};

#include "./Matrix.inl"




#endif

CIMPL 0.1 Code Reference. Copyright © 2004, Baris Sumengen. All rights reserved.