liba 0.1.15
An algorithm library based on C/C++
 
All Data Structures Files Functions Variables Typedefs Enumerator Modules Pages
Loading...
Searching...
No Matches
linear algebra functions
Collaboration diagram for linear algebra functions:

Functions

void a_real_T1 (a_uint n, a_real *A)
 transpose an n x n square matrix in-place.
 
void a_real_T2 (a_uint m, a_uint n, a_real const *__restrict A, a_real *__restrict T)
 transpose a given m x n matrix A into an n x m matrix T.
 
void a_real_mulmm (a_uint row, a_uint c_r, a_uint col, a_real const *__restrict X, a_real const *__restrict Y, a_real *__restrict Z)
 multiply two matrices X and Y, storing the result in Z.
 
void a_real_mulTm (a_uint c_r, a_uint row, a_uint col, a_real const *__restrict X, a_real const *__restrict Y, a_real *__restrict Z)
 multiply the transpose of matrix X with matrix Y, storing the result in Z.
 
void a_real_mulmT (a_uint row, a_uint col, a_uint c_r, a_real const *__restrict X, a_real const *__restrict Y, a_real *__restrict Z)
 multiply matrix X with the transpose of matrix Y, storing the result in Z.
 
void a_real_mulTT (a_uint row, a_uint c_r, a_uint col, a_real const *__restrict X, a_real const *__restrict Y, a_real *__restrict Z)
 multiply the transpose of matrix X with the transpose of matrix Y, storing the result in Z.
 
int a_real_plu (a_uint n, a_real *A, a_uint *p, int *sign)
 compute LU decomposition of a square matrix with partial pivoting.
 
void a_real_plu_P (a_uint n, a_uint const *p, a_real *P)
 construct the permutation matrix P from a permutation vector p.
 
void a_real_plu_P_ (a_uint n, a_uint const *p, a_real *P)
 
void a_real_plu_L (a_uint n, a_real const *A, a_real *L)
 extract the lower triangular matrix L from matrix A.
 
void a_real_plu_U (a_uint n, a_real const *A, a_real *U)
 extract the upper triangular matrix U from matrix A.
 
void a_real_plu_apply (a_uint n, a_uint const *p, a_real const *b, a_real *Pb)
 apply the permutation P to the vector b, producing Pb.
 
void a_real_plu_lower (a_uint n, a_real const *L, a_real *y)
 solve the lower triangular system Ly = Pb for y.
 
void a_real_plu_lower_ (a_uint n, a_real const *L, a_real *y)
 
void a_real_plu_upper (a_uint n, a_real const *U, a_real *x)
 solve the upper triangular system Ux = y for x.
 
void a_real_plu_upper_ (a_uint n, a_real const *U, a_real *x)
 
void a_real_plu_solve (a_uint n, a_real const *A, a_uint const *p, a_real const *b, a_real *x)
 solve the linear system Ax = b using LU decomposition with partial pivoting.
 
void a_real_plu_inv (a_uint n, a_real const *A, a_uint const *p, a_real *b, a_real *I)
 compute the inverse of a matrix using its LU decomposition and permutation matrix.
 
void a_real_plu_inv_ (a_uint n, a_real const *A, a_uint const *p, a_real *I)
 
a_real a_real_plu_det (a_uint n, a_real const *A, int sign)
 compute the determinant of a matrix using its LU decomposition.
 
a_real a_real_plu_lndet (a_uint n, a_real const *A)
 compute the natural logarithm of the absolute value of the determinant of a matrix using its LU decomposition.
 
int a_real_plu_sgndet (a_uint n, a_real const *A, int sign)
 compute the sign of the determinant of a matrix using its LU decomposition.
 
int a_real_llt (a_uint n, a_real *A)
 compute Cholesky decomposition of a symmetric positive-definite matrix.
 
void a_real_llt_L (a_uint n, a_real const *A, a_real *L)
 extract the lower triangular matrix L from matrix A.
 
void a_real_llt_lower (a_uint n, a_real const *L, a_real *y)
 solve the lower triangular system Ly = b for y.
 
void a_real_llt_lower_ (a_uint n, a_real const *L, a_real *y)
 
void a_real_llt_upper (a_uint n, a_real const *L, a_real *x)
 solve the upper triangular system L^T x = y for x.
 
void a_real_llt_upper_ (a_uint n, a_real const *L, a_real *x)
 
void a_real_llt_solve (a_uint n, a_real const *A, a_real *x)
 solve the linear system Ax = b using the Cholesky factorization A = LL^T.
 
void a_real_llt_inv (a_uint n, a_real const *A, a_real *b, a_real *I)
 compute the inverse of a matrix using its Cholesky factorization A = LL^T.
 
void a_real_llt_inv_ (a_uint n, a_real const *A, a_real *I)
 
a_real a_real_llt_det (a_uint n, a_real const *A)
 compute the determinant of a matrix using its Cholesky decomposition.
 
a_real a_real_llt_lndet (a_uint n, a_real const *A)
 compute the natural logarithm of the absolute value of the determinant of a matrix using its Cholesky decomposition.
 

Detailed Description

Function Documentation

◆ a_real_llt()

int a_real_llt ( a_uint n,
a_real * A )

compute Cholesky decomposition of a symmetric positive-definite matrix.

Parameters
[in]nthe order of the square matrix A (number of rows and columns).
[in,out]Aan n x n square matrix. on input, contains the matrix to decompose. on output, contains the L matrix.
Returns
0 on success, or a non-zero error code if the decomposition fails.
Return values
-1on failure, A is a singular matrix.

◆ a_real_llt_det()

a_real a_real_llt_det ( a_uint n,
a_real const * A )

compute the determinant of a matrix using its Cholesky decomposition.

Parameters
[in]nthe order of the square matrix A (number of rows and columns).
[in]Athe matrix containing L form after Cholesky decomposition.
Returns
the determinant of matrix A.

◆ a_real_llt_inv()

void a_real_llt_inv ( a_uint n,
a_real const * A,
a_real * b,
a_real * I )

compute the inverse of a matrix using its Cholesky factorization A = LL^T.

Parameters
[in]nthe order of the square matrix A (number of rows and columns).
[in]Athe matrix containing L form after Cholesky decomposition.
[in]ba pre-allocated temporary buffer of size n for intermediate computations.
[out]Ithe output matrix where the inverse of A will be stored.

◆ a_real_llt_L()

void a_real_llt_L ( a_uint n,
a_real const * A,
a_real * L )

extract the lower triangular matrix L from matrix A.

Parameters
[in]nthe order of the square matrix that was decomposed.
[in]Athe matrix containing L form after Cholesky decomposition.
[out]Lthe output matrix where the lower triangular matrix will be stored.

◆ a_real_llt_lndet()

a_real a_real_llt_lndet ( a_uint n,
a_real const * A )

compute the natural logarithm of the absolute value of the determinant of a matrix using its Cholesky decomposition.

Parameters
[in]nthe order of the square matrix A (number of rows and columns).
[in]Athe matrix containing L form after Cholesky decomposition.
Returns
the natural logarithm of the absolute value of the determinant.

◆ a_real_llt_lower()

void a_real_llt_lower ( a_uint n,
a_real const * L,
a_real * y )

solve the lower triangular system Ly = b for y.

Parameters
[in]nthe order of the square matrix L (number of rows and columns).
[in]Lthe lower triangular matrix L, stored in row-major order.
[in,out]yon input, contains the vector b. on output, contains the solution vector y.

◆ a_real_llt_solve()

void a_real_llt_solve ( a_uint n,
a_real const * A,
a_real * x )

solve the linear system Ax = b using the Cholesky factorization A = LL^T.

Parameters
[in]nthe order of the square matrix A (number of rows and columns).
[in]Athe matrix containing L form after Cholesky decomposition.
[in,out]xon input, contains the vector b. on output, contains the solution vector x.

◆ a_real_llt_upper()

void a_real_llt_upper ( a_uint n,
a_real const * L,
a_real * x )

solve the upper triangular system L^T x = y for x.

Parameters
[in]nthe order of the square matrix L (number of rows and columns).
[in]Lthe lower triangular matrix L, stored in row-major order.
[in,out]xon input, contains the vector y. on output, contains the solution vector x.

◆ a_real_mulmm()

void a_real_mulmm ( a_uint row,
a_uint c_r,
a_uint col,
a_real const *__restrict X,
a_real const *__restrict Y,
a_real *__restrict Z )

multiply two matrices X and Y, storing the result in Z.

ZZrc=XXrnYYnc=[x11x1nxr1xrn][y11y1cyn1ync]=[(x11y11++x1nyn1)(x11y1c++x1nync)(xr1y11++xrnyn1)(xr1y1c++xrnync)]

Parameters
[in]rowrows matrix Z and rows in matrix X.
[in]c_rcolumns in matrix X and rows in matrix Y.
[in]colcolumns in matrix Z and columns in matrix Y.
[in]Xthe first input matrix.
[in]Ythe second input matrix.
[out]Zthe output matrix where the result will be stored.

◆ a_real_mulmT()

void a_real_mulmT ( a_uint row,
a_uint col,
a_uint c_r,
a_real const *__restrict X,
a_real const *__restrict Y,
a_real *__restrict Z )

multiply matrix X with the transpose of matrix Y, storing the result in Z.

ZZrc=XXrnYYcnT=[x11x1nxr1xrn][y11y1nyc1ycn]T=[(x11y11++x1ny1n)(x11yc1++x1nycn)(xr1y11++xrny1n)(xr1yc1++xrnycn)]

Parameters
[in]rowrows matrix Z and rows in matrix X.
[in]colcolumns in matrix Z and rows in matrix Y.
[in]c_rcolumns in matrix X and columns in matrix Y.
[in]Xthe first input matrix.
[in]Ythe second input matrix that will be transposed during multiplication.
[out]Zthe output matrix where the result will be stored.

◆ a_real_mulTm()

void a_real_mulTm ( a_uint c_r,
a_uint row,
a_uint col,
a_real const *__restrict X,
a_real const *__restrict Y,
a_real *__restrict Z )

multiply the transpose of matrix X with matrix Y, storing the result in Z.

ZZrc=XXnrTYYnc=[x11x1rxn1xnr]T[y11y1cyn1ync]=[(x11y11++xn1yn1)(x11y1c++xn1ync)(x1ry11++xnryn1)(x1ry1c++xnrync)]=[x11y11x11y1cx1ry11x1ry1c]++[xn1yn1xn1yncxnryn1xnrync]

Parameters
[in]c_rrows in matrix X and rows in matrix Y.
[in]rowrows in matrix Z and columns in matrix X.
[in]colcolumns in matrix Z and columns in matrix Y.
[in]Xthe first input matrix that will be transposed during multiplication.
[in]Ythe second input matrix.
[out]Zthe output matrix where the result will be stored.

◆ a_real_mulTT()

void a_real_mulTT ( a_uint row,
a_uint c_r,
a_uint col,
a_real const *__restrict X,
a_real const *__restrict Y,
a_real *__restrict Z )

multiply the transpose of matrix X with the transpose of matrix Y, storing the result in Z.

ZZrc=XXnrTYYcnT=[x11x1rxn1xnr]T[y11y1nyc1ycn]T=[(x11y11++xn1y1n)(x11yc1++xn1ycn)(x1ry11++xnry1n)(x1ryc1++xnrycn)]

Parameters
[in]rowrows matrix Z and columns in matrix X.
[in]c_rrows in matrix X and columns in matrix Y.
[in]colcolumns in matrix Z and rows in matrix Y.
[in]Xthe first input matrix that will be transposed during multiplication.
[in]Ythe second input matrix that will be transposed during multiplication.
[out]Zthe output matrix where the result will be stored.

◆ a_real_plu()

int a_real_plu ( a_uint n,
a_real * A,
a_uint * p,
int * sign )

compute LU decomposition of a square matrix with partial pivoting.

This function performs an LU decomposition on the given square matrix A, where L is a lower triangular matrix, and U is an upper triangular matrix. Partial pivoting is used to improve numerical stability during the decomposition process. The result is stored in the original matrix A, with L stored below, and U stored in the diagonal and above. Additionally, it calculates a permutation matrix P that records the row exchanges made during partial pivoting, and determines the sign of the permutation (which can be used to find the determinant's sign).

Parameters
[in]nthe order of the square matrix A (number of rows and columns).
[in,out]Aan n x n square matrix. on input, contains the matrix to decompose. on output, contains the L and U matrices.
[out]pthe row permutation indices after partial pivoting.
[out]signstore the sign of the permutation (+1 or -1).
Returns
0 on success, or a non-zero error code if the decomposition fails.
Return values
-1on failure, A is a singular matrix.

◆ a_real_plu_apply()

void a_real_plu_apply ( a_uint n,
a_uint const * p,
a_real const * b,
a_real * Pb )

apply the permutation P to the vector b, producing Pb.

Parameters
[in]nthe order of the square matrix that was decomposed.
[in]pthe row permutation indices after partial pivoting.
[in]bthe input vector of size n that will be permuted.
[out]Pbthe output vector where the permuted result will be stored.

◆ a_real_plu_det()

a_real a_real_plu_det ( a_uint n,
a_real const * A,
int sign )

compute the determinant of a matrix using its LU decomposition.

Parameters
[in]nthe order of the square matrix A (number of rows and columns).
[in]Athe matrix containing L and U in a compact form after LU decomposition.
[in]signthe sign of the permutation matrix P (+1 or -1).
Returns
the determinant of matrix A.

◆ a_real_plu_inv()

void a_real_plu_inv ( a_uint n,
a_real const * A,
a_uint const * p,
a_real * b,
a_real * I )

compute the inverse of a matrix using its LU decomposition and permutation matrix.

Parameters
[in]nthe order of the square matrix A (number of rows and columns).
[in]Athe matrix containing L and U in a compact form after LU decomposition.
[in]pthe permutation indices obtained during LU decomposition.
[in]ba pre-allocated temporary buffer of size n for intermediate computations.
[out]Ithe output matrix where the inverse of A will be stored.

◆ a_real_plu_L()

void a_real_plu_L ( a_uint n,
a_real const * A,
a_real * L )

extract the lower triangular matrix L from matrix A.

Parameters
[in]nthe order of the square matrix that was decomposed.
[in]Athe matrix containing L and U in a compact form after LU decomposition.
[out]Lthe output matrix where the lower triangular matrix will be stored.

◆ a_real_plu_lndet()

a_real a_real_plu_lndet ( a_uint n,
a_real const * A )

compute the natural logarithm of the absolute value of the determinant of a matrix using its LU decomposition.

Parameters
[in]nthe order of the square matrix A (number of rows and columns).
[in]Athe matrix containing L and U in a compact form after LU decomposition.
Returns
the natural logarithm of the absolute value of the determinant.

◆ a_real_plu_lower()

void a_real_plu_lower ( a_uint n,
a_real const * L,
a_real * y )

solve the lower triangular system Ly = Pb for y.

Parameters
[in]nthe order of the square matrix L (number of rows and columns).
[in]Lthe lower triangular matrix L, stored in row-major order.
[in,out]yon input, contains the permuted vector Pb. on output, contains the solution vector y.

◆ a_real_plu_P()

void a_real_plu_P ( a_uint n,
a_uint const * p,
a_real * P )

construct the permutation matrix P from a permutation vector p.

Parameters
[in]nthe order of the square matrix that was decomposed.
[in]pthe row permutation indices after partial pivoting.
[out]Pthe output matrix where the permutation matrix will be stored.

◆ a_real_plu_sgndet()

int a_real_plu_sgndet ( a_uint n,
a_real const * A,
int sign )

compute the sign of the determinant of a matrix using its LU decomposition.

Parameters
[in]nthe order of the square matrix A (number of rows and columns).
[in]Athe matrix containing L and U in a compact form after LU decomposition.
[in]signthe sign of the permutation matrix P (+1 or -1).
Returns
the sign of the determinant: -1, 0, +1.

◆ a_real_plu_solve()

void a_real_plu_solve ( a_uint n,
a_real const * A,
a_uint const * p,
a_real const * b,
a_real * x )

solve the linear system Ax = b using LU decomposition with partial pivoting.

Parameters
[in]nthe order of the square matrix A (number of rows and columns).
[in]Athe matrix containing L and U in a compact form after LU decomposition.
[in]pthe permutation indices obtained during LU decomposition.
[in]bthe input vector b of the linear system.
[out]xthe output vector x where the solution will be stored.

◆ a_real_plu_U()

void a_real_plu_U ( a_uint n,
a_real const * A,
a_real * U )

extract the upper triangular matrix U from matrix A.

Parameters
[in]nthe order of the square matrix that was decomposed.
[in]Athe matrix containing L and U in a compact form after LU decomposition.
[out]Uthe output matrix where the upper triangular matrix will be stored.

◆ a_real_plu_upper()

void a_real_plu_upper ( a_uint n,
a_real const * U,
a_real * x )

solve the upper triangular system Ux = y for x.

Parameters
[in]nthe order of the square matrix U (number of rows and columns).
[in]Uthe upper triangular matrix U, stored in row-major order.
[in,out]xon input, contains the vector y. on output, contains the solution vector x.

◆ a_real_T1()

void a_real_T1 ( a_uint n,
a_real * A )

transpose an n x n square matrix in-place.

Parameters
[in]norder of square matrix A
[in,out]Aan n x n square matrix

◆ a_real_T2()

void a_real_T2 ( a_uint m,
a_uint n,
a_real const *__restrict A,
a_real *__restrict T )

transpose a given m x n matrix A into an n x m matrix T.

Parameters
[in]mrows in the input matrix A.
[in]ncolumns in the input matrix A.
[in]Athe input matrix A (m x n), stored in row-major order.
[out]Tthe output matrix where the transposed matrix T (n x m) will be stored.