template class mrpt::math::CMatrixFixed

A compile-time fixed-size numeric matrix container.

It uses a RowMajor element memory layout.

For a complete introduction to Matrices and vectors in MRPT, see: https://www.mrpt.org/Matrices_vectors_arrays_and_Linear_Algebra_MRPT_and_Eigen_classes

See also:

CMatrixDynamic (for dynamic-size matrices)

#include <mrpt/math/CMatrixFixed.h>

template <typename T, std::size_t ROWS, std::size_t COLS>
class CMatrixFixed: public mrpt::math::MatrixBase
{
public:
    // typedefs

    typedef T value_type;

    // construction

    CMatrixFixed();
    CMatrixFixed(TConstructorFlags_Matrices);

    template <size_t N>
    CMatrixFixed(const T(&) vals [N]);

    CMatrixFixed(const T* data);

    template <class Derived>
    CMatrixFixed(const Eigen::MatrixBase<Derived>& m);

    template <typename _Lhs, typename _Rhs, int Option>
    CMatrixFixed(const Eigen::Product<_Lhs, _Rhs, Option>& p);

    template <typename Op, typename Lhs, typename Rhs>
    CMatrixFixed(const Eigen::CwiseBinaryOp<Op, Lhs, Rhs>& p);

    template <typename VectorType, int Size>
    CMatrixFixed(const Eigen::VectorBlock<VectorType, Size>& m);

    CMatrixFixed(const size_type rows, const size_type cols);

    //
methods

    template <class Derived>
    CMatrixFixed& operator = (const Eigen::MatrixBase<Derived>& m);

    template <typename VectorType, int Size>
    CMatrixFixed& operator = (const Eigen::VectorBlock<VectorType, Size>& m);

    template <typename U>
    CMatrixFixed& operator = (const CMatrixDynamic<U>& m);

    void loadFromRawPointer(const T* data);
    void setSize(size_t row, size_t col, ] bool zeroNewElements = false);
    void resize(const matrix_size_t& siz, ] bool zeroNewElements = false);
    constexpr size_type rows() const;
    constexpr size_type cols() const;
    constexpr matrix_size_t size() const;

    template <
        typename EIGEN_MATRIX = eigen_t,
        typename EIGEN_MAP = Eigen::Map<EIGEN_MATRIX, MRPT_MAX_STATIC_ALIGN_BYTES, Eigen::InnerStride<1>>
        >
    EIGEN_MAP asEigen();

    template <
        typename EIGEN_MATRIX = eigen_t,
        typename EIGEN_MAP = Eigen::Map<const EIGEN_MATRIX, MRPT_MAX_STATIC_ALIGN_BYTES,            Eigen::InnerStride<1>>
        >
    EIGEN_MAP asEigen() const;

    const T* data() const;
    T* data();
    T& operator () (int row, int col);
    T& operator () (int i);
    T& operator [] (int i);
    CMatrixFixed<T, ROWS, 1> llt_solve(const CMatrixFixed<T, ROWS, 1>& b) const;
    CMatrixFixed<T, ROWS, 1> lu_solve(const CMatrixFixed<T, ROWS, 1>& b) const;
    void sum_At(const CMatrixFixed<Scalar, ROWS, COLS>& A);
};

// direct descendants

template <class T>
class CQuaternion;

Inherited Members

public:
    //
methods

    void setConstant(const Scalar value);
    void setConstant(size_t nrows, size_t ncols, const Scalar value);
    void setConstant(size_t nrows, const Scalar value);
    void assign(const std::size_t N, const Scalar value);
    void setZero();
    void setZero(size_t nrows, size_t ncols);
    void setZero(size_t nrows);
    static Derived Constant(const Scalar value);
    static Derived Constant(size_t nrows, size_t ncols, const Scalar value);
    static Derived Zero();
    static Derived Zero(size_t nrows, size_t ncols);
    auto block(int start_row, int start_col, int BLOCK_ROWS, int BLOCK_COLS);
    auto block(int start_row, int start_col, int BLOCK_ROWS, int BLOCK_COLS) const;
    auto transpose();
    auto transpose() const;
    auto array();
    auto array() const;
    auto operator - () const;

    template <typename S2, class D2>
    auto operator + (const MatrixVectorBase<S2, D2>& m2) const;

    template <typename S2, class D2>
    void operator += (const MatrixVectorBase<S2, D2>& m2);

    template <typename S2, class D2>
    auto operator - (const MatrixVectorBase<S2, D2>& m2) const;

    template <typename S2, class D2>
    void operator -= (const MatrixVectorBase<S2, D2>& m2);

    template <typename S2, class D2>
    auto operator * (const MatrixVectorBase<S2, D2>& m2) const;

    auto operator * (const Scalar s) const;

    template <int N>
    CMatrixFixed<Scalar, N, 1> tail() const;

    template <int N>
    CMatrixFixed<Scalar, N, 1> head() const;

    Scalar& coeffRef(int r, int c);
    const Scalar& coeff(int r, int c) const;
    Scalar minCoeff(std::size_t& outIndexOfMin) const;
    Scalar minCoeff(std::size_t& rowIdx, std::size_t& colIdx) const;
    Scalar maxCoeff(std::size_t& outIndexOfMax) const;
    Scalar maxCoeff(std::size_t& rowIdx, std::size_t& colIdx) const;
    void operator += (Scalar s);
    void operator -= (Scalar s);
    void operator *= (Scalar s);
    CMatrixDynamic<Scalar> operator * (const CMatrixDynamic<Scalar>& v);
    Derived operator + (const Derived& m2) const;
    void operator += (const Derived& m2);
    Derived operator - (const Derived& m2) const;
    void operator -= (const Derived& m2);
    Derived operator * (const Derived& m2) const;
    Scalar dot(const MatrixVectorBase<Scalar, Derived>& v) const;

    template <typename OTHERMATVEC>
    bool operator == (const OTHERMATVEC& o) const;

    template <typename OTHERMATVEC>
    bool operator != (const OTHERMATVEC& o) const;

    Derived& mvbDerived();
    const Derived& mvbDerived() const;
    auto col(int colIdx);
    auto col(int colIdx) const;
    auto row(int rowIdx);
    auto row(int rowIdx) const;

    template <typename VECTOR_LIKE>
    void extractRow(int rowIdx, VECTOR_LIKE& v) const;

    template <typename VECTOR_LIKE>
    VECTOR_LIKE extractRow(int rowIdx) const;

    template <typename VECTOR_LIKE>
    void extractColumn(int colIdx, VECTOR_LIKE& v) const;

    template <typename VECTOR_LIKE>
    VECTOR_LIKE extractColumn(int colIdx) const;

    template <int BLOCK_ROWS, int BLOCK_COLS>
    CMatrixFixed<Scalar, BLOCK_ROWS, BLOCK_COLS> extractMatrix(
        const int start_row = 0,
        const int start_col = 0
        ) const;

    CMatrixDynamic<Scalar> extractMatrix(
        const int BLOCK_ROWS,
        const int BLOCK_COLS,
        const int start_row,
        const int start_col
        ) const;

    Derived& mbDerived();
    const Derived& mbDerived() const;
    void setIdentity();
    void setIdentity(const std::size_t N);
    static Derived Identity();
    static Derived Identity(const std::size_t N);

Typedefs

typedef T value_type

The type of the matrix elements.

Construction

CMatrixFixed()

Default constructor, initializes all elements to zero.

CMatrixFixed(TConstructorFlags_Matrices)

Constructor which leaves the matrix uninitialized.

Example of usage: CMatrixFixed<double,3,2> M(mrpt::math::UNINITIALIZED_MATRIX);

template <size_t N>
CMatrixFixed(const T(&) vals [N])

Initializes from a C array with RowMajor values.

CMatrixFixed(const T* data)

Initializes from a plain buffer with RowMajor values.

template <class Derived>
CMatrixFixed(const Eigen::MatrixBase<Derived>& m)

Convert from Eigen matrix.

template <typename _Lhs, typename _Rhs, int Option>
CMatrixFixed(const Eigen::Product<_Lhs, _Rhs, Option>& p)

Convert from Eigen product.

template <typename Op, typename Lhs, typename Rhs>
CMatrixFixed(const Eigen::CwiseBinaryOp<Op, Lhs, Rhs>& p)

Convert from Eigen binary op.

template <typename VectorType, int Size>
CMatrixFixed(const Eigen::VectorBlock<VectorType, Size>& m)

Convert from Eigen block.

CMatrixFixed(const size_type rows, const size_type cols)

Convenient ctor from size: in this class, it throws if size does not match compile-time size.

It is provided for the sake of offering a uniform API with CMatrixDynamic.

Methods

template <class Derived>
CMatrixFixed& operator = (const Eigen::MatrixBase<Derived>& m)

Assignment from an Eigen matrix.

template <typename VectorType, int Size>
CMatrixFixed& operator = (const Eigen::VectorBlock<VectorType, Size>& m)

Assignment from an Eigen vector block.

template <typename U>
CMatrixFixed& operator = (const CMatrixDynamic<U>& m)

Assignment from a Dynamic matrix.

void loadFromRawPointer(const T* data)

Initializes from a plain buffer with RowMajor values.

Unsafe, prefer loadFromArray() wherever possible, to ensure buffer length checks.

void setSize(size_t row, size_t col, ] bool zeroNewElements = false)

Throws if size does not match with the fixed matrix size.

void resize(const matrix_size_t& siz, ] bool zeroNewElements = false)

Throws if size does not match with the fixed matrix size.

constexpr size_type rows() const

Number of rows in the matrix.

See also:

rows()

constexpr size_type cols() const

Number of columns in the matrix.

See also:

rows()

constexpr matrix_size_t size() const

Get a 2-vector with [NROWS NCOLS] (as in MATLAB command size(x))

template <
    typename EIGEN_MATRIX = eigen_t,
    typename EIGEN_MAP = Eigen::Map<EIGEN_MATRIX, MRPT_MAX_STATIC_ALIGN_BYTES, Eigen::InnerStride<1>>
    >
EIGEN_MAP asEigen()

Get as an Eigen-compatible Eigen::Map object.

template <
    typename EIGEN_MATRIX = eigen_t,
    typename EIGEN_MAP = Eigen::Map<const EIGEN_MATRIX, MRPT_MAX_STATIC_ALIGN_BYTES,            Eigen::InnerStride<1>>
    >
EIGEN_MAP asEigen() const

This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts.

const T* data() const

Return raw pointer to row-major data buffer.

All matrix cells can be assumed to be stored contiguously in memory, i.e. row stride = column count.

T* data()

This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts.

T& operator () (int row, int col)

Access (row,col), without out-of-bounds check (except in Debug builds)

T& operator () (int i)

Access the i-th element, Row-Major order, without out-of-bounds check (except in Debug builds)

T& operator [] (int i)

Access the [i-th] element (for 1xN or Nx1 matrices)

CMatrixFixed<T, ROWS, 1> llt_solve(const CMatrixFixed<T, ROWS, 1>& b) const

Solves the linear system Ax=b, returns x, with A this symmetric matrix.

See also:

lu_solve()

CMatrixFixed<T, ROWS, 1> lu_solve(const CMatrixFixed<T, ROWS, 1>& b) const

Solves the linear system Ax=b, returns x, with A this asymmetric matrix.

See also:

llt_solve()

void sum_At(const CMatrixFixed<Scalar, ROWS, COLS>& A)

this += A T