To handle symmetric matrices in skyline storage format. More...

#include <SkSMatrix.h>

Inheritance diagram for SkSMatrix< T_ >:

Public Member Functions
	SkSMatrix ()
	Default constructor. More...

	SkSMatrix (size_t size, int is_diagonal=false)
	Constructor that initializes a dense symmetric matrix. More...

	SkSMatrix (Mesh &mesh, size_t dof=0, int is_diagonal=false)
	Constructor using mesh to initialize skyline structure of matrix. More...

	SkSMatrix (const Vect< size_t > &ColHt)
	Constructor that initializes skyline structure of matrix using vector of column height. More...

	SkSMatrix (const Vect< size_t > &I, const Vect< size_t > &J, int opt=1)
	Constructor for a square matrix using non zero row and column indices. More...

	SkSMatrix (const Vect< size_t > &I, const Vect< size_t > &J, const Vect< T_ > &a, int opt=1)
	Constructor for a square matrix using non zero row and column indices. More...

	SkSMatrix (const SkSMatrix< T_ > &m)
	Copy Constructor.

	~SkSMatrix ()
	Destructor.

void	setMesh (Mesh &mesh, size_t dof=0)
	Determine mesh graph and initialize matrix. More...

void	setSkyline (Mesh &mesh)
	Determine matrix structure. More...

void	setDiag ()
	Store diagonal entries in a separate internal vector.

void	set (size_t i, size_t j, const T_ &val)
	Assign a value to an entry ofthe matrix. More...

void	Axpy (T_ a, const SkSMatrix< T_ > &m)
	Add to matrix the product of a matrix by a scalar. More...

void	Axpy (T_ a, const Matrix< T_ > *m)
	Add to matrix the product of a matrix by a scalar. More...

void	MultAdd (const Vect< T_ > &x, Vect< T_ > &y) const
	Multiply matrix by vector `x` and add to `y`. More...

void	MultAdd (T_ a, const Vect< T_ > &x, Vect< T_ > &y) const
	Multiply matrix by vector `a*x` and add to `y`. More...

void	Mult (const Vect< T_ > &x, Vect< T_ > &y) const
	Multiply matrix by vector `x` and save in `y` More...

void	TMult (const Vect< T_ > &x, Vect< T_ > &y) const
	Multiply transpose of matrix by vector x and save in y. More...

void	add (size_t i, size_t j, const T_ &val)
	Add a constant to an entry of the matrix. More...

size_t	getColHeight (size_t i) const
	Return column height. More...

Vect< T_ >	getColumn (size_t j) const
	Get `j`-th column vector.

Vect< T_ >	getRow (size_t i) const
	Get `i`-th row vector.

T_ &	operator() (size_t i, size_t j)
	Operator () (Non constant version). More...

T_	operator() (size_t i, size_t j) const
	Operator () (Constant version). More...

SkSMatrix< T_ > &	operator= (const SkSMatrix< T_ > &m)
	Operator =. More...

SkSMatrix< T_ > &	operator= (const T_ &x)
	Operator =. More...

SkSMatrix< T_ > &	operator+= (const SkSMatrix< T_ > &m)
	Operator +=. More...

SkSMatrix< T_ > &	operator*= (const T_ &x)
	Operator *=. More...

int	setLDLt ()
	Factorize matrix (LDLt (Crout) factorization). More...

int	solveLDLt (const Vect< T_ > &b, Vect< T_ > &x)
	Solve a linear system using the LDLt (Crout) factorization. More...

int	solve (Vect< T_ > &b, bool fact=true)
	Solve linear system. More...

int	solve (const Vect< T_ > &b, Vect< T_ > &x, bool fact=true)
	Solve linear system. More...

T_ *	get ()
	Return C-Array. More...

void	set (size_t i, T_ x)
	Assign a value to the i-th entry of C-array containing matrix.

void	add (size_t i, const T_ &val)
	Add `val` to entry `i`.

T_	get (size_t i, size_t j) const
	Return entry `(i,j)` of matrix if this one is stored, `0` else.

Public Member Functions inherited from Matrix< T_ >
	Matrix ()
	Default constructor. More...

	Matrix (const Matrix< T_ > &m)
	Copy Constructor.

virtual	~Matrix ()
	Destructor.

virtual void	reset ()
	Set matrix to 0 and reset factorization parameter. More...

size_t	getNbRows () const
	Return number of rows.

size_t	getNbColumns () const
	Return number of columns.

void	setPenal (real_t p)
	Set Penalty Parameter (For boundary condition prescription).

void	setDiagonal ()
	Set the matrix as diagonal.

T_	getDiag (size_t k) const
	Return `k`-th diagonal entry of matrix. More...

size_t	size () const
	Return matrix dimension (Number of rows and columns).

virtual void	MultAdd (const Vect< T_ > &x, Vect< T_ > &y) const =0
	Multiply matrix by vector `x` and add to `y`

virtual void	MultAdd (T_ a, const Vect< T_ > &x, Vect< T_ > &y) const =0
	Multiply matrix by vector `a*x` and add to `y`

virtual void	Mult (const Vect< T_ > &x, Vect< T_ > &y) const =0
	Multiply matrix by vector `x` and save in `y`

virtual void	TMult (const Vect< T_ > &v, Vect< T_ > &w) const =0
	Multiply transpose of matrix by vector `x` and save in `y`

virtual void	Axpy (T_ a, const Matrix< T_ > *x)=0
	Add to matrix the product of a matrix by a scalar. More...

void	setDiagonal (Mesh &mesh)
	Initialize matrix storage in the case where only diagonal terms are stored. More...

virtual void	clear ()
	brief Set all matrix entries to zero

void	Assembly (const Element &el, T_ *a)
	Assembly of element matrix into global matrix. More...

void	Assembly (const Side &sd, T_ *a)
	Assembly of side matrix into global matrix. More...

void	Prescribe (Vect< T_ > &b, const Vect< T_ > &u, int flag=0)
	Impose by a penalty method an essential boundary condition, using the Mesh instance provided by the constructor. More...

void	Prescribe (int dof, int code, Vect< T_ > &b, const Vect< T_ > &u, int flag=0)
	Impose by a penalty method an essential boundary condition to a given degree of freedom for a given code. More...

void	Prescribe (Vect< T_ > &b, int flag=0)
	Impose by a penalty method a homegeneous (=0) essential boundary condition. More...

void	Prescribe (size_t dof, Vect< T_ > &b, const Vect< T_ > &u, int flag=0)
	Impose by a penalty method an essential boundary condition when only one DOF is treated. More...

void	PrescribeSide ()
	Impose by a penalty method an essential boundary condition when DOFs are supported by sides. More...

virtual void	add (size_t i, size_t j, const T_ &val)=0
	Add `val` to entry `(i,j)`.

virtual void	add (size_t i, const T_ &val)=0
	Add `val` to entry `i`.

virtual int	Factor ()=0
	Factorize matrix. Available only if the storage class enables it.

virtual int	solve (Vect< T_ > &b, bool fact=true)=0
	Solve the linear system. More...

virtual int	solve (const Vect< T_ > &b, Vect< T_ > &x, bool fact=true)=0
	Solve the linear system. More...

int	FactorAndSolve (Vect< T_ > &b)
	Factorize matrix and solve the linear system. More...

int	FactorAndSolve (const Vect< T_ > &b, Vect< T_ > &x)
	Factorize matrix and solve the linear system. More...

size_t	getLength () const
	Return number of stored terms in matrix.

int	isDiagonal () const
	Say if matrix is diagonal or not.

int	isFactorized () const
	Say if matrix is factorized or not. More...

virtual size_t	getColInd (size_t i) const
	Return Column index for column `i` (See the description for class SpMatrix).

virtual size_t	getRowPtr (size_t i) const
	Return Row pointer for row `i` (See the description for class SpMatrix).

virtual void	set (size_t i, size_t j, const T_ &val)=0
	Assign a value to an entry of the matrix. More...

virtual T_ &	operator() (size_t i, size_t j)=0
	Operator () (Non constant version). More...

virtual T_	operator() (size_t i, size_t j) const =0
	Operator () (Non constant version). More...

T_	operator() (size_t i) const
	Operator () with one argument (Constant version). More...

T_ &	operator() (size_t i)
	Operator () with one argument (Non Constant version). More...

T_ &	operator[] (size_t k)
	Operator [] (Non constant version). More...

T_	operator[] (size_t k) const
	Operator [] (Constant version). More...

Matrix &	operator= (Matrix< T_ > &m)
	Operator =. More...

Matrix &	operator+= (const Matrix< T_ > &m)
	Operator +=. More...

Matrix &	operator-= (const Matrix< T_ > &m)
	Operator -=. More...

Matrix &	operator= (const T_ &x)
	Operator =. More...

Matrix &	operator*= (const T_ &x)
	Operator *=. More...

Matrix &	operator+= (const T_ &x)
	Operator +=. More...

Matrix &	operator-= (const T_ &x)
	Operator -=. More...

virtual T_	get (size_t i, size_t j) const =0
	Return entry `(i,j)` of matrix if this one is stored, `0` else.

Detailed Description

template<class T_>
class OFELI::SkSMatrix< T_ >

To handle symmetric matrices in skyline storage format.

This template class allows storing and manipulating a symmetric matrix in skyline storage format.

The matrix entries are stored column by column as in the following example:

    /                       \ 
    | a0   a1    0   0   a7 |
    |      a2   a3   0   a8 |
    | ...       a4  a5   a9 |
    |               a6  a10 |
    |                   a11 |
    \                       /

Template Parameters

T_	Data type (double, float, complex<double>, ...)

Author: Rachid Touzani

Copyright: GNU Lesser Public License

Constructor & Destructor Documentation

◆ SkSMatrix() [1/6]

SkSMatrix ( )

Default constructor.

Initializes a zero-dimension matrix

◆ SkSMatrix() [2/6]

SkSMatrix	(	size_t	size,
		int	is_diagonal = `false`
	)

Constructor that initializes a dense symmetric matrix.

Normally, for a dense matrix this is not the right class.

Parameters

[in]	size	Number of matrix rows (and columns).
[in]	is_diagonal	Boolean to select if the matrix is diagonal or not [Default: false]

◆ SkSMatrix() [3/6]

SkSMatrix	(	Mesh &	mesh,
		size_t	dof = `0`,
		int	is_diagonal = `false`
	)

Constructor using mesh to initialize skyline structure of matrix.

Parameters

[in]	mesh	Mesh instance for which matrix graph is determined.
[in]	dof	Option parameter, with default value `0`. `dof=1` means that only one degree of freedom for each node (or element or side) is taken to determine matrix structure. The value `dof=0` means that matrix structure is determined using all DOFs.
[in]	is_diagonal	Boolean argument to say is the matrix is actually a diagonal matrix or not.

◆ SkSMatrix() [4/6]

SkSMatrix ( const Vect< size_t > & ColHt )

Constructor that initializes skyline structure of matrix using vector of column height.

Parameters

[in] ColHt Vect instance that contains rows lengths of matrix.

◆ SkSMatrix() [5/6]

SkSMatrix	(	const Vect< size_t > &	I,
		const Vect< size_t > &	J,
		int	opt = `1`
	)

Constructor for a square matrix using non zero row and column indices.

Parameters

[in]	I	Vector containing row indices
[in]	J	Vector containing column indices
[in]	opt	Flag indicating if vectors I and J are cleaned and ordered (opt=1) or not (opt=0). In the latter case, these vectors can contain the same contents more than once and are not necessarily ordered.

◆ SkSMatrix() [6/6]

SkSMatrix	(	const Vect< size_t > &	I,
		const Vect< size_t > &	J,
		const Vect< T_ > &	a,
		int	opt = `1`
	)

Constructor for a square matrix using non zero row and column indices.

Parameters

[in]	I	Vector containing row indices
[in]	J	Vector containing column indices
[in]	a	Vector containing matrix entries in the same order than the one given by `I` and `J`
[in]	opt	Flag indicating if vectors `I` and `J` are cleaned and ordered (`opt=1`) or not (`opt=0`). In the latter case, these vectors can contain the same contents more than once and are not necessarily ordered

Member Function Documentation

◆ add()

void add	(	size_t	i,
		size_t	j,
		const T_ &	val
	)

virtual

Add a constant to an entry of the matrix.

Parameters

[in]	i	Row index
[in]	j	Column index
[in]	val	Constant value to add to `a(i,j)`

Implements Matrix< T_ >.

◆ Axpy() [1/2]

void Axpy	(	T_	a,
		const Matrix< T_ > *	m
	)

virtual

Add to matrix the product of a matrix by a scalar.

Parameters

[in]	a	Scalar to premultiply
[in]	m	Pointer to Matrix by which `a` is multiplied. The result is added to current instance

Implements Matrix< T_ >.

◆ Axpy() [2/2]

void Axpy	(	T_	a,
		const SkSMatrix< T_ > &	m
	)

Add to matrix the product of a matrix by a scalar.

Parameters

[in]	a	Scalar to premultiply
[in]	m	Matrix by which `a` is multiplied. The result is added to current instance

◆ get()

T_ * get ( )

Return C-Array.

Skyline of matrix is stored row by row.

◆ getColHeight()

size_t getColHeight ( size_t i ) const

Return column height.

Column height at entry i is returned.

◆ Mult()

void Mult	(	const Vect< T_ > &	x,
		Vect< T_ > &	y
	)		const

virtual

Multiply matrix by vector x and save in y

Parameters

[in]	x	Vector to multiply by matrix
[out]	y	Vector that contains on output the result.

Implements Matrix< T_ >.

◆ MultAdd() [1/2]

void MultAdd	(	const Vect< T_ > &	x,
		Vect< T_ > &	y
	)		const

virtual

Multiply matrix by vector x and add to y.

Parameters

[in]	x	Vector to multiply by matrix
[in,out]	y	Vector to add to the result. `y` contains on output the result.

Implements Matrix< T_ >.

◆ MultAdd() [2/2]

void MultAdd	(	T_	a,
		const Vect< T_ > &	x,
		Vect< T_ > &	y
	)		const

virtual

Multiply matrix by vector a*x and add to y.

Parameters

[in]	a	Constant to multiply by matrix
[in]	x	Vector to multiply by matrix
[in,out]	y	Vector to add to the result. `y` contains on output the result.

Implements Matrix< T_ >.

◆ operator()() [1/2]

T_ & operator()	(	size_t	i,
		size_t	j
	)

virtual

Operator () (Non constant version).

Parameters

[in]	i	Row index
[in]	j	Column index

Warning: To modify a value of an entry of the matrix it is safer not to modify both lower and upper triangles. Otherwise, wrong values will be assigned. If not sure, use the member functions set or add.

Implements Matrix< T_ >.

◆ operator()() [2/2]

T_ operator()	(	size_t	i,
		size_t	j
	)		const

virtual

Operator () (Constant version).

Parameters

[in]	i	Row index
[in]	j	Column index

Implements Matrix< T_ >.

◆ operator*=()

SkSMatrix< T_ > & operator*= ( const T_ & x )

Operator *=.

Premultiply matrix entries by constant value x.

◆ operator+=()

SkSMatrix< T_ > & operator+= ( const SkSMatrix< T_ > & m )

Operator +=.

Add matrix m to current matrix instance.

◆ operator=() [1/2]

SkSMatrix< T_ > & operator= ( const SkSMatrix< T_ > & m )

Operator =.

Copy matrix m to current matrix instance.

◆ operator=() [2/2]

SkSMatrix< T_ > & operator= ( const T_ & x )

Operator =.

define the matrix as a diagonal one with all diagonal entries equal to x.

◆ set()

void set	(	size_t	i,
		size_t	j,
		const T_ &	val
	)

virtual

Assign a value to an entry ofthe matrix.

Parameters

[in]	i	Row index
[in]	j	Column index
[in]	val	Value to assign to `a(i,j)`

Implements Matrix< T_ >.

◆ setLDLt()

int setLDLt ( )

Factorize matrix (LDLt (Crout) factorization).

Returns

0 if factorization was normally performed
n if the n-th pivot is null

◆ setMesh()

void setMesh	(	Mesh &	mesh,
		size_t	dof = `0`
	)

Determine mesh graph and initialize matrix.

This member function is called by constructor with the same arguments

Parameters

[in]	mesh	Mesh instance for which matrix graph is determined.
[in]	dof	Option parameter, with default value `0`. `dof=1` means that only one degree of freedom for each node (or element or side) is taken to determine matrix structure. The value `dof=0` means that matrix structure is determined using all DOFs.

◆ setSkyline()

void setSkyline ( Mesh & mesh )

Determine matrix structure.

This member function calculates matrix structure using Mesh instance mesh.

◆ solve() [1/2]

int solve	(	const Vect< T_ > &	b,
		Vect< T_ > &	x,
		bool	fact = `true`
	)

virtual

Solve linear system.

The linear system having the current instance as a matrix is solved by using the LDLt decomposition. Solution is thus realized after a factorization step and a forward/backward substitution step. The factorization step is realized only if this was not already done.
Note that this function modifies the matrix contents is a factorization is performed. Naturally, if the the matrix has been modified after using this function, the user has to refactorize it using the function setLDLt. This is because the class has no non-expensive way to detect if the matrix has been modified. The function setLDLt realizes the factorization step only.

Parameters

[in]	b	Vect instance that contains right-hand side.
[out]	x	Vect instance that contains solution
[in]	fact	Set true if matrix is to be factorized (Default value), false if not

Returns

0 if solution was normally performed,
n if the n-th pivot is null.

Implements Matrix< T_ >.

◆ solve() [2/2]

int solve	(	Vect< T_ > &	b,
		bool	fact = `true`
	)

virtual

Solve linear system.

The linear system having the current instance as a matrix is solved by using the LDLt decomposition. Solution is thus realized after a factorization step and a forward/backward substitution step. The factorization step is realized only if this was not already done.
Note that this function modifies the matrix contents is a factorization is performed. Naturally, if the the matrix has been modified after using this function, the user has to refactorize it using the function setLU. This is because the class has no non-expensive way to detect if the matrix has been modified. The function setLDLt realizes the factorization step only.

Parameters

[in,out]	b	Vect instance that contains right-hand side on input and solution on output.
[in]	fact	Set true if matrix is to be factorized (Default value), false if not

Returns

0 if solution was normally performed,
n if the n-th pivot is null.

Implements Matrix< T_ >.

◆ solveLDLt()

int solveLDLt	(	const Vect< T_ > &	b,
		Vect< T_ > &	x
	)

Solve a linear system using the LDLt (Crout) factorization.

This function solves a linear system. The LDLt factorization is performed if this was not already done using the function setLU.

Parameters

[in]	b	Vect instance that contains right-hand side
[out]	x	Vect instance that contains solution

Returns

0 if solution was normally performed,
n if the n-th pivot is null

Solution is performed only is factorization has previouly been invoked.

◆ TMult()

void TMult	(	const Vect< T_ > &	x,
		Vect< T_ > &	y
	)		const

virtual

Multiply transpose of matrix by vector x and save in y.

Parameters

[in]	x	Vector to multiply by matrix
[out]	y	Vector that contains on output the result.

Implements Matrix< T_ >.

Public Member Functions

Detailed Description

Constructor & Destructor Documentation

◆ SkSMatrix() [1/6]

◆ SkSMatrix() [2/6]

◆ SkSMatrix() [3/6]

◆ SkSMatrix() [4/6]

◆ SkSMatrix() [5/6]

◆ SkSMatrix() [6/6]

Member Function Documentation

◆ add()

◆ Axpy() [1/2]

◆ Axpy() [2/2]

◆ get()

◆ getColHeight()

◆ Mult()

◆ MultAdd() [1/2]

◆ MultAdd() [2/2]

◆ operator()() [1/2]

◆ operator()() [2/2]

◆ operator*=()

◆ operator+=()

◆ operator=() [1/2]

◆ operator=() [2/2]

◆ set()

◆ setLDLt()

◆ setMesh()

◆ setSkyline()

◆ solve() [1/2]

◆ solve() [2/2]

◆ solveLDLt()

◆ TMult()