Matrix Types

Right now, MTL4 provides three matrix types:

The type dense2D defines regular row-major and column-major matrices:

// File: dense2D.cpp

#include <iostream>
#include <boost/numeric/mtl/mtl.hpp>

int main(int argc, char* argv[])
{
    using namespace mtl;

    // a is a row-major matrix
    dense2D<double>               a(10, 10);

    // Matrices are not initialized by default
    set_to_zero(a);

    // Assign a value to a matrix element
    a(2, 3)= 7.0;

    // You can also use a more C-like notation
    a[2][4]= 3.0;

    std::cout << "a is \n" << a << "\n";
    
    // b is a column-major matrix
    dense2D<float, matrix::parameters<tag::col_major> > b(10, 10);

    // Assign a three times the identity to b
    b= 3;
    std::cout << "b is \n" << b << "\n";

    return 0;
}

If no matrix parameters are defined, dense matrices are by default row-major. There are more matrix parameters besides the orientation. As they are not yet fully supported we refrain from discussing them.

Matrix elements can be accessed by a(i, j) or in the more familiar form a[i][j]. The second form is internally transformed into the first one at compile time so that the run-time performance is not affected (unless the compiler does not inline completely which we never observed so far). Also, the compile time is not conceivably increased by this transformation.

Please notice that overwriting single matrix elements is only defined for dense matrix types. For a generic way to modify matrices see Matrix Insertion.

Assigning a scalar value to a matrix stores a multiple of the identity matrix, i.e. the scalar is assigned to all diagonal elements and all off-diagonal elements are 0. If the matrix is not square this assignment throws an exception. This operation is generic (i.e. applicable to all matrix types including sparse).

Just in case you wonder why the scalar value is only assigned to the diagonal elements of the matrix not to all entries, this becomes quite clear when you think of a matrix as a linear operator (from one vector space to another one). For instance, consider the multiplication of vector x with the scalar alpha:

    y= alpha * x;
where y is a vector too. This operation is equivalent to assigning alpha to the matrix A and multiplying x with A:
    A= alpha;
    y= A * x;
In other words, the matrix A has the same impact on x as the scalar alpha itself.

Assigning the scalar value to the diagonal requires of course that the matrix is square. In the special case that the scalar value is 0 (more precisely the multiplicative identity element of the matrix's value_type) the matrix can be non-square. This is consistent with the linear operator characteristic: applying the zero operator on some vector results in the zero vector with the dimension of the operators image. From a more pragmatic prospective

    A= 0; 
can be used to clear any matrix, square or rectangular, sparse and dense.

Dense matrices with a recursively designed memory layout can be defined with the type morton_dense:

// File: morton_dense.cpp

#include <iostream>
#include <boost/numeric/mtl/mtl.hpp>

int main(int argc, char* argv[])
{
    using namespace mtl;

    // Z-order matrix
    morton_dense<double, recursion::morton_z_mask>  a(10, 10);

    set_to_zero(a);
    a(2, 3)= 7.0;
    a[2][4]= 3.0;
    std::cout << "a is \n" << a << "\n";
    
    // b is a N-order matrix with column-major 4x4 blocks
    morton_dense<float, recursion::doppled_4_col_mask> b(10, 10);

    // Assign a three times the identity to b
    b= 3;
    std::cout << "b is \n" << b << "\n";

    return 0;
}

A detailed description will be added soon.

Sparse matrices are defined with the type compressed2D:

// File: compressed2D.cpp

#include <iostream>
#include <boost/numeric/mtl/mtl.hpp>

int main(int argc, char* argv[])
{
    using namespace mtl;

    // CRS matrix
    compressed2D<double>   a(12, 12);

    // Laplace operator discretized on a 3x4 grid
    matrix::laplacian_setup(a, 3, 4);
    std::cout << "a is \n" << a;
    
    // Element access is allowed for reading
    std::cout << "a[3][2] is " << a[3][2] << "\n\n";
    
    // CCS matrix
    compressed2D<float, matrix::parameters<tag::col_major> > b(10, 10);
    // Assign a three times the identity to b
    b= 3;
    std::cout << "b is \n" << b << "\n";

    return 0;
}

Matrix a is stored as compressed row storage (CRS). Its assigned values correspond to a discretized Laplace operator. To change or insert single elements of a compressed matrix is not supported. Especially for very large matrices, this would result in an unbearable performance burden.

However, it is allowed to assign a scalar value to the entire matrix given it is square as in the example. Matrix b is stored in compressed column storage (CCS).

Which orientation is favorable dependents on the performed operations and might require some experimentation. All operations are provided in the same way for both formats

How to fill sparse matrices is shown in the following chapter.

Return to Vector Types                                Table of Content                                Proceed to Vector Insertion






Matrix Types -- MTL 4 -- Peter Gottschling and Andrew Lumsdaine -- Generated on 19 May 2009 by Doxygen 1.5.5 -- Copyright 2007 by the Trustees of Indiana University.