GammaLib  2.1.0.dev
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
GMatrixSymmetric Class Reference

Symmetric matrix class interface definition. More...

#include <GMatrixSymmetric.hpp>

Inheritance diagram for GMatrixSymmetric:
GMatrixBase GBase

Public Member Functions

 GMatrixSymmetric (void)
 Void constructor. More...
 
 GMatrixSymmetric (const int &rows, const int &columns)
 Matrix constructor. More...
 
 GMatrixSymmetric (const GMatrix &matrix)
 GMatrix to GMatrixSymmetric storage class convertor. More...
 
 GMatrixSymmetric (const GMatrixSparse &matrix)
 GMatrixSparse to GMatrixSymmetric storage class convertor. More...
 
 GMatrixSymmetric (const GMatrixSymmetric &matrix)
 Copy constructor. More...
 
virtual ~GMatrixSymmetric (void)
 Destructor. More...
 
virtual double & operator() (const int &row, const int &column)
 Return reference to matrix element. More...
 
virtual const double & operator() (const int &row, const int &column) const
 Return reference to matrix element (const version) More...
 
virtual GVector operator* (const GVector &v) const
 Vector multiplication. More...
 
virtual GMatrixSymmetricoperator= (const GMatrixSymmetric &matrix)
 Matrix assignment operator. More...
 
virtual GMatrixSymmetricoperator= (const double &value)
 Value assignment operator. More...
 
virtual GMatrixSymmetric operator+ (const GMatrixSymmetric &matrix) const
 Binary matrix addition. More...
 
virtual GMatrixSymmetric operator+ (const double &scalar) const
 Binary matrix scalar addition. More...
 
virtual GMatrixSymmetric operator- (const GMatrixSymmetric &matrix) const
 Binary matrix subtraction. More...
 
virtual GMatrixSymmetric operator- (const double &scalar) const
 Binary matrix scalar subtraction. More...
 
virtual GMatrix operator* (const GMatrixSymmetric &matrix) const
 Binary matrix multiplication operator. More...
 
virtual GMatrixSymmetric operator- (void) const
 Negate matrix elements. More...
 
virtual GMatrixSymmetricoperator+= (const GMatrixSymmetric &matrix)
 Unary matrix addition operator. More...
 
virtual GMatrixSymmetricoperator+= (const double &scalar)
 Unary matrix scalar addition operator. More...
 
virtual GMatrixSymmetricoperator-= (const GMatrixSymmetric &matrix)
 Unary matrix subtraction operator. More...
 
virtual GMatrixSymmetricoperator-= (const double &scalar)
 Unary matrix scalar subtraction operator. More...
 
virtual GMatrixSymmetricoperator*= (const double &scaler)
 Scale matrix elements. More...
 
virtual GMatrixSymmetricoperator/= (const double &scalar)
 Divide matrix elements. More...
 
virtual void clear (void)
 Clear matrix. More...
 
virtual GMatrixSymmetricclone (void) const
 Clone matrix. More...
 
virtual std::string classname (void) const
 Return class name. More...
 
virtual double & at (const int &row, const int &column)
 Return reference to matrix element. More...
 
virtual const double & at (const int &row, const int &column) const
 Return reference to matrix element (const version) More...
 
virtual GVector row (const int &row) const
 Extract row as vector from matrix. More...
 
virtual void row (const int &row, const GVector &vector)
 Set row in matrix. More...
 
virtual GVector column (const int &column) const
 Extract column as vector from matrix. More...
 
virtual void column (const int &column, const GVector &vector)
 Set matrix column from vector. More...
 
virtual void add_to_row (const int &row, const GVector &vector)
 Add row to matrix elements. More...
 
virtual void add_to_column (const int &column, const GVector &vector)
 Add vector column into matrix. More...
 
virtual double fill (void) const
 Determine fill of matrix. More...
 
virtual double min (void) const
 Return minimum matrix element. More...
 
virtual double max (void) const
 Return maximum matrix element. More...
 
virtual double sum (void) const
 Sum matrix elements. More...
 
virtual std::string print (const GChatter &chatter=NORMAL) const
 Print matrix. More...
 
GMatrixSymmetric transpose (void) const
 Return transposed matrix. More...
 
GMatrixSymmetric invert (void) const
 Return inverted matrix. More...
 
GVector solve (const GVector &vector) const
 Solves linear matrix equation. More...
 
GMatrixSymmetric abs (void) const
 Return absolute of matrix. More...
 
GMatrix extract_lower_triangle (void) const
 Extract lower triangle of matrix. More...
 
GMatrix extract_upper_triangle (void) const
 Extract upper triangle of matrix. More...
 
GMatrixSymmetric cholesky_decompose (const bool &compress=true) const
 Return Cholesky decomposition. More...
 
GVector cholesky_solver (const GVector &vector, const bool &compress=true) const
 Cholesky solver. More...
 
GMatrixSymmetric cholesky_invert (const bool &compress=true) const
 Invert matrix using a Cholesky decomposition. More...
 
- Public Member Functions inherited from GMatrixBase
 GMatrixBase (void)
 Void constructor. More...
 
 GMatrixBase (const GMatrixBase &matrix)
 Copy constructor. More...
 
virtual ~GMatrixBase (void)
 Destructor. More...
 
virtual GMatrixBaseoperator= (const GMatrixBase &matrix)
 Assignment operator. More...
 
virtual bool operator== (const GMatrixBase &matrix) const
 Equalty operator. More...
 
virtual bool operator!= (const GMatrixBase &matrix) const
 Non-equality operator. More...
 
bool is_empty (void) const
 Check if matrix is empty. More...
 
const int & size (void) const
 Return number of matrix elements. More...
 
const int & columns (void) const
 Return number of matrix columns. More...
 
const int & rows (void) const
 Return number of matrix rows. More...
 
- Public Member Functions inherited from GBase
virtual ~GBase (void)
 Destructor. More...
 

Private Member Functions

void init_members (void)
 Initialise class mambers. More...
 
void copy_members (const GMatrixSymmetric &matrix)
 Copy class members. More...
 
void free_members (void)
 Delete class members. More...
 
void alloc_members (const int &rows, const int &columns)
 Allocates matrix memory. More...
 
void set_inx (void)
 Set index selection. More...
 

Private Attributes

int m_num_inx
 Number of indices in array. More...
 
int * m_inx
 Index array of non-zero rows/columns. More...
 

Additional Inherited Members

- Protected Member Functions inherited from GMatrixBase
void init_members (void)
 Initialise class members. More...
 
void copy_members (const GMatrixBase &matrix)
 Copy class members. More...
 
void free_members (void)
 Delete class members. More...
 
void select_non_zero (void)
 Setup compressed matrix factorisation. More...
 
void scale_elements (const double &scalar)
 Scale all matrix elements with a scalar. More...
 
void set_all_elements (const double &value)
 Set all elements to a specific value. More...
 
double get_min_element (void) const
 Return minimum matrix element. More...
 
double get_max_element (void) const
 Returns maximum matrix element. More...
 
double get_element_sum (void) const
 Returns sum over all matrix elements. More...
 
std::string print_elements (const GChatter &chatter=NORMAL, const int &num=15) const
 Print all matrix elements. More...
 
std::string print_row_compression (const GChatter &chatter=NORMAL) const
 Print row compression scheme if it exists. More...
 
std::string print_col_compression (const GChatter &chatter=NORMAL) const
 Print column compression scheme if it exists. More...
 
- Protected Attributes inherited from GMatrixBase
int m_rows
 Number of rows. More...
 
int m_cols
 Number of columns. More...
 
int m_elements
 Number of elements stored in matrix. More...
 
int m_alloc
 Size of allocated matrix memory. More...
 
int m_num_rowsel
 Number of selected rows (for compressed decomposition) More...
 
int m_num_colsel
 Number of selected columns (for compressed decomposition) More...
 
int * m_colstart
 Column start indices (m_cols+1) More...
 
int * m_rowsel
 Row selection (for compressed decomposition) More...
 
int * m_colsel
 Column selection (for compressed decomposition) More...
 
double * m_data
 Matrix data. More...
 

Detailed Description

Symmetric matrix class interface definition.

This class implements a symmetric matrix. The class stores only a triangle physically, imposing thus strict matrix symmetry.

For a description of the common matrix methods, please refer to the GMatrixBase class.

Matrix allocation is done using the constructors

GMatrixSymmetric symmetricmatrix(rows, columns);
GMatrixSymmetric symmetricmatrix(matrix);
GMatrixSymmetric symmetricmatrix(sparsematrix);
GMatrixSymmetric symmetricmatrix(symmetricmatrix);

where rows and columns specify the number of rows and columns of the matrix. Storage conversion constructors exist that allow allocating a symmetric matrix by copying from a sparse matrix of type GMatrixSparse and a general matrix of type GMatrix. Exceptions will be thrown in case that the matrix from which the object should be allocated is not symmetric.

Methods are available to extract the lower or the upper triangle of the matrix:

matrix.extract_lower_triangle();
matrix.extract_upper_triangle();

Definition at line 71 of file GMatrixSymmetric.hpp.

Constructor & Destructor Documentation

GMatrixSymmetric::GMatrixSymmetric ( void  )

Void constructor.

Constructs an empty matrix. An empty matrix has zero rows and columns.

Definition at line 72 of file GMatrixSymmetric.cpp.

References init_members().

Referenced by clone().

GMatrixSymmetric::GMatrixSymmetric ( const int &  rows,
const int &  columns 
)
explicit

Matrix constructor.

Parameters
[in]rowsNumber of rows [>=0].
[in]columnsNumber of columns [>=0].
Exceptions
GException::invalid_argumentNumber of rows or columns is negative.

Constructs a matrix with the specified number of rows and columns.

Definition at line 92 of file GMatrixSymmetric.cpp.

References alloc_members(), G_CONSTRUCTOR, init_members(), and gammalib::str().

GMatrixSymmetric::GMatrixSymmetric ( const GMatrix matrix)

GMatrix to GMatrixSymmetric storage class convertor.

Parameters
[in]matrixGeneral matrix (GMatrix).
Exceptions
GException::invalid_argumentMatrix is not symmetric.

Converts a general matrix into the symmetric storage class. If the input matrix is not symmetric, an exception is thrown.

Definition at line 133 of file GMatrixSymmetric.cpp.

References alloc_members(), GMatrixBase::columns(), G_MATRIX, init_members(), GMatrixBase::m_cols, GMatrixBase::m_rows, row(), GMatrixBase::rows(), and gammalib::str().

GMatrixSymmetric::GMatrixSymmetric ( const GMatrixSparse matrix)

GMatrixSparse to GMatrixSymmetric storage class convertor.

Parameters
[in]matrixSparse matrix (GMatrixSparse).
Exceptions
GException::invalid_argumentSparse matrix is not symmetric.

Converts a sparse matrix into the symmetric storage class. If the input matrix is not symmetric, an exception is thrown.

Definition at line 178 of file GMatrixSymmetric.cpp.

References alloc_members(), GMatrixBase::columns(), G_MATRIX, init_members(), GMatrixBase::m_cols, GMatrixBase::m_rows, row(), GMatrixBase::rows(), and gammalib::str().

GMatrixSymmetric::GMatrixSymmetric ( const GMatrixSymmetric matrix)

Copy constructor.

Parameters
[in]matrixMatrix.

Definition at line 217 of file GMatrixSymmetric.cpp.

References copy_members(), and init_members().

GMatrixSymmetric::~GMatrixSymmetric ( void  )
virtual

Destructor.

Definition at line 234 of file GMatrixSymmetric.cpp.

References free_members().

Member Function Documentation

GMatrixSymmetric GMatrixSymmetric::abs ( void  ) const

Return absolute of matrix.

Returns
Absolute of matrix

Returns matrix where all elements of the matrix have been replaced by their absolute values.

Definition at line 894 of file GMatrixSymmetric.cpp.

References abs(), GMatrixBase::m_cols, GMatrixBase::m_data, GMatrixBase::m_elements, and GMatrixBase::m_rows.

void GMatrixSymmetric::add_to_column ( const int &  column,
const GVector vector 
)
virtual

Add vector column into matrix.

Parameters
[in]columnMatrix column [0,...,columns()[.
[in]vectorVector.
Exceptions
GException::out_of_rangeInvalid matrix column specified.
GException::matrix_vector_mismatchVector size does not match number of matrix rows.

Adds the content of a vector to a matrix column.

Implements GMatrixBase.

Definition at line 804 of file GMatrixSymmetric.cpp.

References column(), G_ADD_TO_COLUMN, GMatrixBase::m_cols, GMatrixBase::m_rows, row(), GVector::size(), and gammalib::str().

void GMatrixSymmetric::add_to_row ( const int &  row,
const GVector vector 
)
virtual

Add row to matrix elements.

Parameters
[in]rowMatrix row [0,...,rows()[.
[in]vectorVector of matrix row elements (columns() elements).
Exceptions
GException::out_of_rangeInvalid matrix row specified.
GException::matrix_vector_mismatchVector does not match the matrix dimensions.
Todo:
To be implemented.

Implements GMatrixBase.

Definition at line 766 of file GMatrixSymmetric.cpp.

References G_ADD_TO_ROW, GMatrixBase::m_cols, GMatrixBase::m_rows, GVector::size(), and gammalib::str().

void GMatrixSymmetric::alloc_members ( const int &  rows,
const int &  columns 
)
private

Allocates matrix memory.

Parameters
[in]rowsNumber of rows.
[in]columnsNumber of columns.
Exceptions
GException::invalid_argumentMatrix is not symmetric.

This method is the main constructor code that allocates and initialises memory for matrix elements. The method assumes that no memory has been allocated for the matrix elements, the column start index array and the index array. The method allocates the memory for matrix elements, the column start indices and the index array, sets all matrix elements to 0.0, and sets the column start indices. The content of the index array is undefined.

Todo:
Verify if the index array m_inx should be initialized.

Definition at line 1492 of file GMatrixSymmetric.cpp.

References GMatrixBase::columns(), G_ALLOC_MEMBERS, GMatrixBase::m_alloc, GMatrixBase::m_cols, GMatrixBase::m_colstart, GMatrixBase::m_data, GMatrixBase::m_elements, m_inx, GMatrixBase::m_rows, GMatrixBase::rows(), and gammalib::str().

Referenced by GMatrixSymmetric().

double & GMatrixSymmetric::at ( const int &  row,
const int &  column 
)
virtual

Return reference to matrix element.

Parameters
[in]rowMatrix row [0,...,rows()[.
[in]columnMatrix column [0,...,columns()[.
Returns
Reference to matrix element.
Exceptions
GException::out_of_rangeMatrix row or column out of range.

Implements GMatrixBase.

Definition at line 572 of file GMatrixSymmetric.cpp.

References column(), G_AT, GMatrixBase::m_cols, GMatrixBase::m_colstart, GMatrixBase::m_data, GMatrixBase::m_rows, and row().

const double & GMatrixSymmetric::at ( const int &  row,
const int &  column 
) const
virtual

Return reference to matrix element (const version)

Parameters
[in]rowMatrix row [0,...,rows()[.
[in]columnMatrix column [0,...,columns()[.
Returns
Reference to matrix element.
Exceptions
GException::out_of_rangeMatrix row or column out of range.

Implements GMatrixBase.

Definition at line 601 of file GMatrixSymmetric.cpp.

References column(), G_AT, GMatrixBase::m_cols, GMatrixBase::m_colstart, GMatrixBase::m_data, GMatrixBase::m_rows, and row().

GMatrixSymmetric GMatrixSymmetric::cholesky_decompose ( const bool &  compress = true) const

Return Cholesky decomposition.

Parameters
[in]compressUse zero-row/column compression (defaults to true).
Returns
Cholesky decomposition of matrix
Exceptions
GException::runtime_errorMatrix is not positive definite.

Returns the Cholesky decomposition of a sparse matrix. The decomposition is stored within a GMatrixSymmetric object.

The method is inspired by the algorithm found in Numerical Recipes. The decomposition, which is a matrix occupying only the lower triange, is stored in the elements of the symmetric matrix. To visualise the matrix one has to use 'lower_triangle()' to extract the relevant part. Case A operates on a full matrix, Case B operates on a (logically) compressed matrix where zero rows/columns have been removed.

Definition at line 1037 of file GMatrixSymmetric.cpp.

References G_CHOL_DECOMP, GMatrixBase::m_cols, GMatrixBase::m_colstart, GMatrixBase::m_data, m_inx, m_num_inx, GMatrixBase::m_rows, row(), set_inx(), sqrt(), gammalib::str(), and sum().

Referenced by cholesky_invert(), and solve().

GMatrixSymmetric GMatrixSymmetric::cholesky_invert ( const bool &  compress = true) const

Invert matrix using a Cholesky decomposition.

Parameters
[in]compressUse zero-row/column compression (defaults to true).
Returns
Inverted matrix.
Exceptions
GException::runtime_errorAll matrix elements are zero.

Inverts the matrix using a Cholesky decomposition.

The method distinguish two cases. Case A operates on a full matrix while Case B operates on a (logically) compressed matrix where all zero rows/columns are skipped.

Definition at line 1250 of file GMatrixSymmetric.cpp.

References cholesky_decompose(), G_CHOL_INVERT, GMatrixBase::m_cols, GMatrixBase::m_colstart, GMatrixBase::m_data, m_inx, m_num_inx, GMatrixBase::m_rows, row(), and sum().

Referenced by invert().

GVector GMatrixSymmetric::cholesky_solver ( const GVector vector,
const bool &  compress = true 
) const

Cholesky solver.

Parameters
[in]vectorVector for which should be solved.
[in]compressUse zero-row/column compression (default: true).
Exceptions
GException::invalid_argumentMatrix and vector do not match.
GException::runtime_errorAll matrix elements are zero.

Solves the linear equation A*x=b using a Cholesky decomposition of A. This function is to be applied on a decomposition GMatrixSymmetric matrix that is produced by 'cholesky_decompose'. Case A operates on a full matrix, Case B on a zero rows/columns (logically) compressed matrix.

Definition at line 1152 of file GMatrixSymmetric.cpp.

References G_CHOL_SOLVE, GMatrixBase::m_colstart, GMatrixBase::m_data, m_inx, m_num_inx, GMatrixBase::m_rows, row(), GVector::size(), gammalib::str(), and sum().

Referenced by solve().

std::string GMatrixSymmetric::classname ( void  ) const
inlinevirtual

Return class name.

Returns
String containing the class name ("GMatrixSymmetric").

Implements GMatrixBase.

Definition at line 152 of file GMatrixSymmetric.hpp.

void GMatrixSymmetric::clear ( void  )
virtual

Clear matrix.

Implements GMatrixBase.

Definition at line 537 of file GMatrixSymmetric.cpp.

References free_members(), and init_members().

GMatrixSymmetric * GMatrixSymmetric::clone ( void  ) const
virtual

Clone matrix.

Returns
Pointer to deep copy of matrix.

Implements GMatrixBase.

Definition at line 555 of file GMatrixSymmetric.cpp.

References GMatrixSymmetric().

GVector GMatrixSymmetric::column ( const int &  column) const
virtual

Extract column as vector from matrix.

Parameters
[in]columnMatrix column [0,...,columns()[.
Returns
Vector of matrix column elements (rows() elements).
Exceptions
GException::out_of_rangeInvalid matrix column specified.

This method extracts a matrix column into a vector.

Implements GMatrixBase.

Definition at line 686 of file GMatrixSymmetric.cpp.

References G_EXTRACT_COLUMN, GMatrixBase::m_cols, GMatrixBase::m_rows, and row().

Referenced by add_to_column(), at(), column(), GMatrixSparse::GMatrixSparse(), and operator()().

void GMatrixSymmetric::column ( const int &  column,
const GVector vector 
)
virtual

Set matrix column from vector.

Parameters
[in]columnMatrix column [0,...,columns()[.
[in]vectorVector.
Exceptions
GException::out_of_rangeInvalid matrix column specified.
GException::invalid_argumentVector size does not match number of matrix rows.

Inserts the content of a vector into a matrix column. Any previous content in the matrix column will be overwritted.

Implements GMatrixBase.

Definition at line 723 of file GMatrixSymmetric.cpp.

References column(), G_SET_COLUMN, GMatrixBase::m_cols, GMatrixBase::m_rows, row(), GVector::size(), and gammalib::str().

void GMatrixSymmetric::copy_members ( const GMatrixSymmetric matrix)
private

Copy class members.

Parameters
[in]matrixMatrix.

Definition at line 1442 of file GMatrixSymmetric.cpp.

References GMatrixBase::m_cols, m_inx, and m_num_inx.

Referenced by GMatrixSymmetric(), and operator=().

GMatrix GMatrixSymmetric::extract_lower_triangle ( void  ) const

Extract lower triangle of matrix.

This method extracts the lower triangle of a matrix into another matrix. (including the diagonal elements). All remaining matrix elements will be zero.

Definition at line 977 of file GMatrixSymmetric.cpp.

References GMatrixBase::m_cols, GMatrixBase::m_colstart, GMatrixBase::m_data, GMatrixBase::m_rows, and row().

GMatrix GMatrixSymmetric::extract_upper_triangle ( void  ) const

Extract upper triangle of matrix.

This method extracts the upper triangle of a matrix into another matrix. (including the diagonal elements). All remaining matrix elements will be zero.

Definition at line 1001 of file GMatrixSymmetric.cpp.

References GMatrixBase::m_cols, GMatrixBase::m_colstart, GMatrixBase::m_data, GMatrixBase::m_rows, and row().

double GMatrixSymmetric::fill ( void  ) const
virtual

Determine fill of matrix.

Returns
Matrix fill (between 0 and 1).

The fill of a matrix is defined as the number of non-zero elements devided by the total number of matrix elements.

Implements GMatrixBase.

Definition at line 919 of file GMatrixSymmetric.cpp.

References GMatrixBase::m_cols, GMatrixBase::m_data, GMatrixBase::m_elements, GMatrixBase::m_rows, and row().

void GMatrixSymmetric::free_members ( void  )
private

Delete class members.

Definition at line 1463 of file GMatrixSymmetric.cpp.

References m_inx.

Referenced by clear(), operator=(), and ~GMatrixSymmetric().

void GMatrixSymmetric::init_members ( void  )
private

Initialise class mambers.

Definition at line 1426 of file GMatrixSymmetric.cpp.

References m_inx, and m_num_inx.

Referenced by clear(), GMatrixSymmetric(), and operator=().

GMatrixSymmetric GMatrixSymmetric::invert ( void  ) const

Return inverted matrix.

Returns
Inverted matrix.

Returns inverse of matrix. Inversion is done for the moment using Cholesky decomposition. This does not work on any kind of matrix.

Todo:
Specify in documentation for which kind of matrix the method works.

Definition at line 845 of file GMatrixSymmetric.cpp.

References cholesky_invert(), GMatrixBase::m_cols, and GMatrixBase::m_rows.

double GMatrixSymmetric::max ( void  ) const
inlinevirtual

Return maximum matrix element.

Returns
Maximum matrix element

Returns the largest element in the matrix.

Implements GMatrixBase.

Definition at line 196 of file GMatrixSymmetric.hpp.

References GMatrixBase::get_max_element().

double GMatrixSymmetric::min ( void  ) const
inlinevirtual

Return minimum matrix element.

Returns
Minimum matrix element

Returns the smallest element in the matrix.

Implements GMatrixBase.

Definition at line 182 of file GMatrixSymmetric.hpp.

References GMatrixBase::get_min_element().

double & GMatrixSymmetric::operator() ( const int &  row,
const int &  column 
)
virtual

Return reference to matrix element.

Parameters
[in]rowMatrix row [0,...,rows()-1].
[in]columnMatrix column [0,...,columns()-1].
Returns
Reference to matrix element.

Implements GMatrixBase.

Definition at line 312 of file GMatrixSymmetric.cpp.

References column(), GMatrixBase::m_colstart, GMatrixBase::m_data, and row().

const double & GMatrixSymmetric::operator() ( const int &  row,
const int &  column 
) const
virtual

Return reference to matrix element (const version)

Parameters
[in]rowMatrix row [0,...,rows()-1].
[in]columnMatrix column [0,...,columns()-1].
Returns
Reference to matrix element.

Implements GMatrixBase.

Definition at line 330 of file GMatrixSymmetric.cpp.

References column(), GMatrixBase::m_colstart, GMatrixBase::m_data, and row().

GVector GMatrixSymmetric::operator* ( const GVector vector) const
virtual

Vector multiplication.

Parameters
[in]vectorVector.
Exceptions
GException::invalid_argumentVector size differs from number of columns in matrix.

This method performs a vector multiplication of a matrix. The vector multiplication will produce a vector. The multiplication can only be performed when the number of matrix columns is identical to the size of the vector. If the size of the vector differs an exception will be thrown.

Implements GMatrixBase.

Definition at line 356 of file GMatrixSymmetric.cpp.

References G_OP_MUL_VEC, GMatrixBase::m_cols, GMatrixBase::m_rows, row(), GVector::size(), gammalib::str(), and sum().

GMatrix GMatrixSymmetric::operator* ( const GMatrixSymmetric matrix) const
virtual

Binary matrix multiplication operator.

Parameters
[in]matrixMatrix to be multiplied.
Exceptions
GException::invalid_argumentNumber of rows in matrix is different from number of matrix columns.

This method performs a matrix multiplication. Since the product of two symmetric matrices is not necessarily symmetric, this method returns a GMatrix object.

The operation can only succeed when the dimensions of both matrices are compatible.

Definition at line 1616 of file GMatrixSymmetric.cpp.

References G_OP_MAT_MUL, GMatrixBase::m_cols, GMatrixBase::m_rows, row(), gammalib::str(), and sum().

GMatrixSymmetric & GMatrixSymmetric::operator*= ( const double &  scalar)
inlinevirtual

Scale matrix elements.

Parameters
[in]scalarScale factor.
Returns
Matrix with elements multiplied by scalar.

Returns a matrix where all elements have been multiplied by the specified scalar value.

Definition at line 283 of file GMatrixSymmetric.hpp.

References GMatrixBase::scale_elements().

GMatrixSymmetric GMatrixSymmetric::operator+ ( const GMatrixSymmetric matrix) const
inlinevirtual

Binary matrix addition.

Parameters
[in]matrixMatrix.
Returns
Result of matrix addition.

Returns the sum of two matrices. The method makes use of the unary addition operator.

Definition at line 212 of file GMatrixSymmetric.hpp.

GMatrixSymmetric GMatrixSymmetric::operator+ ( const double &  scalar) const
inlinevirtual

Binary matrix scalar addition.

Parameters
[in]scalarScalar.
Returns
Matrix with scalar added.

Returns a matrix where a scalar has been added to each matrix element.

Definition at line 229 of file GMatrixSymmetric.hpp.

GMatrixSymmetric & GMatrixSymmetric::operator+= ( const GMatrixSymmetric matrix)
virtual

Unary matrix addition operator.

Parameters
[in]matrixMatrix.
Exceptions
GException::invalid_argumentIncompatible number of matrix rows or columns.

This method performs a matrix addition. The operation can only succeed when the dimensions of both matrices are identical.

Definition at line 417 of file GMatrixSymmetric.cpp.

References G_OP_ADD, GMatrixBase::m_cols, GMatrixBase::m_data, GMatrixBase::m_elements, GMatrixBase::m_rows, and gammalib::str().

GMatrixSymmetric & GMatrixSymmetric::operator+= ( const double &  scalar)
virtual

Unary matrix scalar addition operator.

Parameters
[in]scalarScalar.

Adds a scalar to each matrix element.

Definition at line 454 of file GMatrixSymmetric.cpp.

References GMatrixBase::m_data, and GMatrixBase::m_elements.

GMatrixSymmetric GMatrixSymmetric::operator- ( const GMatrixSymmetric matrix) const
inlinevirtual

Binary matrix subtraction.

Parameters
[in]matrixMatrix.
Returns
Result of matrix subtraction.

Returns the difference between two matrices. The method makes use of the unary subtraction operator.

Definition at line 247 of file GMatrixSymmetric.hpp.

GMatrixSymmetric GMatrixSymmetric::operator- ( const double &  scalar) const
inlinevirtual

Binary matrix scalar subtraction.

Parameters
[in]scalarScalar.
Returns
Matrix with scalar subtracted.

Returns a matrix where a scalar has been subtracted from each matrix element.

Definition at line 265 of file GMatrixSymmetric.hpp.

GMatrixSymmetric GMatrixSymmetric::operator- ( void  ) const
virtual

Negate matrix elements.

Returns
Matrix with negated elements.

Returns a matrix where each element has been replaced by its negative element.

Definition at line 390 of file GMatrixSymmetric.cpp.

References GMatrixBase::m_data, and GMatrixBase::m_elements.

GMatrixSymmetric & GMatrixSymmetric::operator-= ( const GMatrixSymmetric matrix)
virtual

Unary matrix subtraction operator.

Parameters
[in]matrixMatrix.
Exceptions
GException::invalid_argumentIncompatible number of matrix rows or columns.

This method performs a matrix addition. The operation can only succeed when the dimensions of both matrices are identical.

Definition at line 478 of file GMatrixSymmetric.cpp.

References G_OP_SUB, GMatrixBase::m_cols, GMatrixBase::m_data, GMatrixBase::m_elements, GMatrixBase::m_rows, and gammalib::str().

GMatrixSymmetric & GMatrixSymmetric::operator-= ( const double &  scalar)
virtual

Unary matrix scalar subtraction operator.

Parameters
[in]scalarScalar.

Subtracts a scalar from each matrix element.

Definition at line 515 of file GMatrixSymmetric.cpp.

References GMatrixBase::m_data, and GMatrixBase::m_elements.

GMatrixSymmetric & GMatrixSymmetric::operator/= ( const double &  scalar)
inlinevirtual

Divide matrix elements.

Parameters
[in]scalarScalar.
Returns
Matrix with elements divided by scalar.

Returns a matrix where all elements have been divided by the specified scalar value.

Definition at line 300 of file GMatrixSymmetric.hpp.

GMatrixSymmetric & GMatrixSymmetric::operator= ( const GMatrixSymmetric matrix)
virtual

Matrix assignment operator.

Parameters
[in]matrixMatrix.
Returns
Matrix.

Assigns the content of another matrix to the actual matrix instance.

Definition at line 258 of file GMatrixSymmetric.cpp.

References copy_members(), free_members(), init_members(), and GMatrixBase::operator=().

GMatrixSymmetric & GMatrixSymmetric::operator= ( const double &  value)
virtual

Value assignment operator.

Parameters
[in]valueValue.
Returns
Matrix.

Assigns the specified value to all elements of the matrix.

Definition at line 292 of file GMatrixSymmetric.cpp.

References GMatrixBase::m_data, and GMatrixBase::m_elements.

std::string GMatrixSymmetric::print ( const GChatter chatter = NORMAL) const
virtual
GVector GMatrixSymmetric::row ( const int &  row) const
virtual

Extract row as vector from matrix.

Parameters
[in]rowMatrix row [0,...,rows()[.
Returns
Vector of matrix row elements (columns() elements).
Exceptions
GException::out_of_rangeInvalid matrix row specified.

This method extracts a matrix row into a vector.

Implements GMatrixBase.

Definition at line 631 of file GMatrixSymmetric.cpp.

References G_EXTRACT_ROW, GMatrixBase::m_cols, and GMatrixBase::m_rows.

Referenced by add_to_column(), at(), cholesky_decompose(), cholesky_invert(), cholesky_solver(), column(), extract_lower_triangle(), extract_upper_triangle(), fill(), GMatrixSymmetric(), operator()(), operator*(), set_inx(), and sum().

void GMatrixSymmetric::row ( const int &  row,
const GVector vector 
)
virtual

Set row in matrix.

Parameters
[in]rowMatrix row [0,...,rows()[.
[in]vectorVector of matrix row elements (columns() elements).
Todo:
To be implemented.

Implements GMatrixBase.

Definition at line 661 of file GMatrixSymmetric.cpp.

References G_SET_ROW, and GMatrixBase::m_rows.

void GMatrixSymmetric::set_inx ( void  )
private

Set index selection.

Determines the non-zero rows/columns in matrix and set up index array that points to these rows/columns. This array is used for compressed matrix calculations (Case B in the above routines).

Definition at line 1550 of file GMatrixSymmetric.cpp.

References GMatrixBase::m_cols, GMatrixBase::m_colstart, GMatrixBase::m_data, m_inx, m_num_inx, GMatrixBase::m_rows, and row().

Referenced by cholesky_decompose().

GVector GMatrixSymmetric::solve ( const GVector vector) const

Solves linear matrix equation.

Parameters
[in]vectorSolution vector.

Solves the linear equation

\[M \times {\tt solution} = {\tt vector} \]

where \(M\) is the matrix, \({\tt vector}\) is the result, and \({\tt solution}\) is the solution. Solving is done using Cholesky decomposition. This does not work on any kind of matrix.

Todo:
Specify in documentation for which kind of matrix the method works.

Definition at line 873 of file GMatrixSymmetric.cpp.

References cholesky_decompose(), and cholesky_solver().

double GMatrixSymmetric::sum ( void  ) const
virtual

Sum matrix elements.

Returns
Sum of all matrix elements.

Implements GMatrixBase.

Definition at line 944 of file GMatrixSymmetric.cpp.

References GMatrixBase::m_cols, GMatrixBase::m_colstart, GMatrixBase::m_data, GMatrixBase::m_rows, and row().

Referenced by cholesky_decompose(), cholesky_invert(), cholesky_solver(), and operator*().

GMatrixSymmetric GMatrixSymmetric::transpose ( void  ) const
inline

Return transposed matrix.

Returns
Transposed matrix.

Returns transposed matrix of the matrix. As the transposed matrix of a symmetric matrix is identical to the original matrix, this method simply returns the actual matrix.

Definition at line 168 of file GMatrixSymmetric.hpp.

Member Data Documentation

int* GMatrixSymmetric::m_inx
private

Index array of non-zero rows/columns.

Definition at line 142 of file GMatrixSymmetric.hpp.

Referenced by alloc_members(), cholesky_decompose(), cholesky_invert(), cholesky_solver(), copy_members(), free_members(), init_members(), and set_inx().

int GMatrixSymmetric::m_num_inx
private

Number of indices in array.

Definition at line 141 of file GMatrixSymmetric.hpp.

Referenced by cholesky_decompose(), cholesky_invert(), cholesky_solver(), copy_members(), init_members(), and set_inx().


The documentation for this class was generated from the following files: