Матрицы

В этот разделе описаны шаблоны классов предназначенные для представления математических матриц.
Шаблон IMatrix является интерфейсным, и его классы не выделяют память для элементов матриц:

template <class Type>
class IMatrix
{
    Type * const * const data;
// Запрет оператора присваивания и конструктора копии
    void operator = ( const IMatrix & );
    IMatrix ( const IMatrix & );
public:
    const nat nRow; // к-во строк
    const nat nCol; // к-во столбцов

    IMatrix ( nat r, nat c, Type * const * const p ) : nRow(r), nCol(c), data(p) {}

    IMatrix & fill ( const Type & p )
    {
        for ( nat i = 0; i < nRow; ++i )
        {
            Type * t = data[i];
            for ( nat j = 0; j < nCol; ++j ) t[j] = p;
        }
        return *this;
    }
    IMatrix & operator *= ( const Type & p )
    {
        for ( nat i = 0; i < nRow; ++i )
        {
            Type * t = data[i];
            for ( nat j = 0; j < nCol; ++j ) t[j] *= p;
        }
        return *this;
    }
    operator       Type * const * ()       { return data; }
    operator const Type * const * () const { return data; }
};
Имеются: 1) конструктор, 2) функция-член fill, которая заполняет матрицу указанным значением, 3) оператор *=, который умножает матрицу на число, и 4) два оператора пребразования типа.

Шаблон SMatrix может создавать матрицу в стеке и применяется, когда размер матрицы известен на этапе компиляции:

template <class Type, nat R, nat C>
class SMatrix : public IMatrix<Type>
{
    Type * p[R];
    Type stor[R*C];
public:
    SMatrix () : IMatrix<Type> ( R, C, p )
    {
        p[0] = stor;
        for ( nat i = 1; i < R; ++i ) p[i] = p[i-1] + C;
    }

    SMatrix ( const SMatrix & m ) : IMatrix<Type> ( R, C, p )
    {
        nat i;
        p[0] = stor;
        for ( i = 1; i < R; ++i ) p[i] = p[i-1] + C;
        for ( i = 0; i < R*C; ++i ) stor[i] = m.stor[i];
    }

    SMatrix & operator = ( const SMatrix & m )
    {
        for ( nat i = 0; i < R*C; ++i ) stor[i] = m.stor[i];
        return *this;
    }
};
Шаблон HMatrix располагает элементы матрицы в динамической памяти ( куче ):
template <class Type, nat R, nat C>
class HMatrix : public IMatrix<Type>
{
    Type ** p;
// Запрет оператора присваивания
    HMatrix & operator = ( const HMatrix & m );
public:
    HMatrix ( nat r, nat c ) : IMatrix<Type> ( r, c, p = r > 0 ? new Type*[r] : 0 )
    {
        if ( ! nRow ) return;
        p[0] = new Type[nRow*nCol];
        for ( nat i = 1; i < nRow; ++i ) p[i] = p[i-1] + nCol;
    }

    explicit HMatrix ( const IMatrix<Type> & m ) : 
		IMatrix<Type> ( m.nRow, m.nCol, p = m.nRow > 0 ? new Type*[m.nRow] : 0 )
    {
        if ( ! nRow ) return;
        const nat n = nRow * nCol;
        Type * t = p[0] = new Type[n];
        for ( nat i = 1; i < nRow; ++i ) p[i] = p[i-1] + nCol;
        const Type * s = m[0];
        for ( nat j = 0; j < n; ++j ) t[j] = s[j];
    }

    explicit HMatrix ( const HMatrix<Type> & m ) : 
		IMatrix<Type> ( m.nRow, m.nCol, p = m.nRow > 0 ? new Type*[m.nRow] : 0 )
    {
        if ( ! nRow ) return;
        const nat n = nRow * nCol;
        Type * t = p[0] = new Type[n];
        for ( nat i = 1; i < nRow; ++i ) p[i] = p[i-1] + nCol;
        const Type * s = m[0];
        for ( nat j = 0; j < n; ++j ) t[j] = s[j];
    }

    ~HMatrix()
    {
        if ( p != 0 )
        {
            delete[] p[0];
            delete[] p;
        }
    }
};
Шаблон CmbMatrix является комбинацией двух предыдущих:
template <class Type, nat R, nat C>
class CmbMatrix : public IMatrix<Type>
{
    Type ** p;
    Type * pa[R];
    Type stor[R*C];
// Запрет оператора присваивания
    CmbMatrix & operator = ( const CmbMatrix & m );
public:
    CmbMatrix ( nat r, nat c ) : IMatrix<Type> ( r, c, p = r > R ? new Type*[r] : pa )
    {
        const nat n = nRow * nCol;
        p[0] = n > R*C ? new Type[n] : stor;
        for ( nat i = 1; i < nRow; ++i ) p[i] = p[i-1] + nCol;
    }

    explicit CmbMatrix ( const IMatrix<Type> & m ) : 
		IMatrix<Type> ( m.nRow, m.nCol, p = m.nRow > R ? new Type*[m.nRow] : pa )
    {
        const nat n = nRow * nCol;
        Type * t = p[0] = n > R*C ? new Type[n] : stor;
        for ( nat i = 1; i < nRow; ++i ) p[i] = p[i-1] + nCol;
        const Type * s = m[0];
        for ( nat j = 0; j < n; ++j ) t[j] = s[j];
    }

    explicit CmbMatrix ( const CmbMatrix<Type, R, C> & m ) : 
		IMatrix<Type> ( m.nRow, m.nCol, p = m.nRow > R ? new Type*[m.nRow] : pa )
    {
        const nat n = nRow * nCol;
        Type * t = p[0] = n > R*C ? new Type[n] : stor;
        for ( nat i = 1; i < nRow; ++i ) p[i] = p[i-1] + nCol;
        const Type * s = m[0];
        for ( nat j = 0; j < n; ++j ) t[j] = s[j];
    }

    ~CmbMatrix()
    {
        if ( p[0] != stor ) delete[] p[0];
        if ( p != pa ) delete[] p;
    }
};
Следующая группа функций-шаблонов выполняет копирование, сложение, отнимание и умножение матриц, если размеры матриц согласованы:
template <class Type>
bool copy ( const IMatrix<Type> & a, IMatrix<Type> & b ) // b = a

template <class Type>
bool plus ( const IMatrix<Type> & a, const IMatrix<Type> & b, IMatrix<Type> & c ) // c = a + b

template <class Type>
bool minus ( const IMatrix<Type> & a, const IMatrix<Type> & b, IMatrix<Type> & c ) // c = a - b

template <class Type>
bool multi ( const IMatrix<Type> & a, const IMatrix<Type> & b, IMatrix<Type> & c ) // c = a * b
Следующая группа функций определена для вещественных матриц ( исходники находятся в файле matrix.cpp ):
bool determinant ( const IMatrix<double> & m, double & d );
bool trans  ( const IMatrix<double> & a, IMatrix<double> & b ); // b = aT
bool inverse( const IMatrix<double> & a, IMatrix<double> & b ); // a * b = 1

// Матричные нормы

double norm1 ( const IMatrix<double> & ); // 1-норма
double norm2 ( const IMatrix<double> & ); // 2-норма
double normF ( const IMatrix<double> & ); // норма Фробениуса
double normU ( const IMatrix<double> & ); // бесконечная норма

// Cингулярное разложение: A = U * W * V

bool svd ( const IMatrix<double> & A, IMatrix<double> & U, IMatrix<double> & W, IMatrix<double> & V ); 
Следующая группа функций определена для комплексных матриц ( исходники находятся в файле cmatrix.cpp ):
// Транспонирование

bool trans  ( const IMatrix<Complex> & a, IMatrix<Complex> & b ); // b = aT

// Решение задачи наименьших квадратов при помощи преобразований Хаусхолдера

bool lss_h ( IMatrix<Complex> & a, const Complex * b, Complex * x ); // b[a.nRow], x[a.nCol]
Следующий шаблон создаёт матрицы размер которых можно менять в ходе выполнения программы ( метод resize ). Также у этих матриц можно менять местами строки при помощи метода swapRows.
template <class Type> class DMatrix
{
    Type * data;
    Type ** row;
    nat nRow; // к-во строк
    nat nCol; // к-во столбцов
// Запрет оператора присваивания
    DMatrix & operator = ( const DMatrix & m );
public:
    DMatrix () : nRow(0), nCol(0), row(0), data(0) {}

    DMatrix ( nat r, nat c ) : nRow(r), nCol(c)
    {
        if ( nRow > 0 )
        {
            row = new Type*[nRow];
            if ( nCol > 0 )
            {
                row[0] = data = new Type[nRow*nCol];
                for ( nat i = 1; i < nRow; ++i ) row[i] = row[i-1] + nCol;
            }
            else
            {
                data = 0;
                for ( nat i = 0; i < nRow; ++i ) row[i] = 0;
            }
        }
        else
        {
            row = 0;
            data = 0;
        }
    }

    explicit DMatrix ( const DMatrix<Type> & m ) : nRow ( m.nRow ), nCol ( m.nCol )
    {
        if ( nRow > 0 )
        {
            row = new Type*[nRow];
            if ( nCol > 0 )
            {
                nat i;
                row[0] = data = new Type[nRow*nCol];
                for ( i = 1; i < nRow; ++i ) row[i] = row[i-1] + nCol;
                for ( i = 0; i < nRow; ++i )
                {
                    Type * pi = row[i];
                    const Type * mi = m[i];
                    for ( nat j = 0; j < nCol; ++j ) pi[j] = mi[j];
                }
            }
            else
            {
                data = 0;
                for ( nat i = 0; i < nRow; ++i ) row[i] = 0;
            }
        }
        else
        {
            row = 0;
            data = 0;
        }
    }

    ~DMatrix()
    {
        delete[] data;
        delete[] row;
    }

    nat rowSize () const { return nRow; } // к-во строк
    nat colSize () const { return nCol; } // к-во столбцов

    DMatrix & resize ( nat r, nat c )
    {
        delete[] data;
        delete[] row;
        nRow = r;
        nCol = c;
        if ( r > 0 )
        {
            row = new Type*[r];
            if ( c > 0 )
            {
                row[0] = data = new Type[r*c];
                for ( nat i = 1; i < r; ++i ) row[i] = row[i-1] + c;
            }
            else
            {
                data = 0;
                for ( nat i = 0; i < r; ++i ) row[i] = 0;
            }
        }
        else
        {
            row = 0;
            data = 0;
        }
        return *this;
    }

          Type * operator [] ( nat i )       { return row[i]; }
    const Type * operator [] ( nat i ) const { return row[i]; }

    bool swapRows ( nat i1, nat i2 )
    {
        if ( i1 >= nRow || i2 >= nRow ) return false;
        _swap ( row[i1], row[i2] );
        return true;
    }

    DMatrix & fill ( const Type & p )
    {
        const nat n = nRow * nCol;
        for ( nat i = 0; i < n; ++i ) data[i] = p;
        return *this;
    }

    DMatrix & operator *= ( const Type & p )
    {
        const nat n = nRow * nCol;
        for ( nat i = 0; i < n; ++i ) data[i] *= p;
        return *this;
    }
};
Наверх