ARTS  2.4.0(git:4fb77825)
lapack Namespace Reference

Functions

void dgetrf_ (int *m, int *n, double *A, int *lda, int *ipiv, int *info)
 LU decomposition. More...
 
void zgetrf_ (int *m, int *n, std::complex< double > *A, int *lda, int *ipiv, int *info)
 
void dgetrs_ (char *trans, int *n, int *nrhs, double *A, int *lda, int *ipiv, double *b, int *ldb, int *info)
 Solve linear system of equations. More...
 
void dgetri_ (int *n, double *A, int *lda, int *ipiv, double *work, int *lwork, int *info)
 Matrix inversion. More...
 
void zgetri_ (int *n, std::complex< double > *A, int *lda, int *ipiv, std::complex< double > *work, int *lwork, int *info)
 
void ilaenv_ (int *ispec, char *name, char *opts, int *n1, int *n2, int *n3, int *n4)
 Optimal parameters for computation. More...
 
void dgesvx_ (char *fact, char *trans, int *n, int *nrhs, double *A, int *lda, double *AF, int *ldaf, int *ipiv, char *equed, double *R, double *C, double *B, int *ldb, double *X, int *ldx, double *RCOND, double *FERR, double *BERR, double *WORK, int *IWORK, int *info)
 Solve linear system. More...
 
void dgeev_ (char *jobvl, char *jobvr, int *n, double *A, int *lda, double *WR, double *WI, double *VL, int *ldvl, double *VR, int *ldvr, double *work, int *lwork, double *rwork, int *info)
 
void zgeev_ (char *jobvl, char *jobvr, int *n, std::complex< double > *A, int *lda, std::complex< double > *W, std::complex< double > *VL, int *ldvl, std::complex< double > *VR, int *ldvr, std::complex< double > *work, int *lwork, double *rwork, int *info)
 

Function Documentation

◆ dgeev_()

void lapack::dgeev_ ( char *  jobvl,
char *  jobvr,
int *  n,
double *  A,
int *  lda,
double *  WR,
double *  WI,
double *  VL,
int *  ldvl,
double *  VR,
int *  ldvr,
double *  work,
int *  lwork,
double *  rwork,
int *  info 
)

Referenced by diagonalize().

◆ dgesvx_()

void lapack::dgesvx_ ( char *  fact,
char *  trans,
int *  n,
int *  nrhs,
double *  A,
int *  lda,
double *  AF,
int *  ldaf,
int *  ipiv,
char *  equed,
double *  R,
double *  C,
double *  B,
int *  ldb,
double *  X,
int *  ldx,
double *  RCOND,
double *  FERR,
double *  BERR,
double *  WORK,
int *  IWORK,
int *  info 
)

Solve linear system.

Solves the linear system A*X = B using LAPACK's dgesv function. The use of the expert function is necessary since we need to solve the system for the inverse matrix due to the different memory layouts in arts and BLAS. Since most of the parameters are unused they are not fully documented. See LAPACK documentation.

Parameters
[in]factSpecifies whether or not the factored form of the matrix A is provided and whether the matrix A should be equilibrated.
[in]transSpecifies the form of the linear system = 'N': A * x = B = 'N': A^T * x = B = 'N': A^H * x = B
[in]nDimensionality of the system.
[in]nrhsNumber of right-hand side columns.
[in]AThe double array (column-major order) representing the matrix a A that describes the linear system.
[in]ldaThe leading dimension of A.
[in,out]AFIf a is given in factored form AF must contain the factored form of A on input. On output contains the factored form of A.
[in]ldafThe leading dimension of AF.
[in,out]ipivThe pivot vector from the LU decomposition.
[in]equedSpecifies the form of equilibration that was done.
[in]RThe low scale factors for A.
[in]CThe column scale factors for A.
[in]BThe right hand-side matrix of the system.
[in]ldbThe leading dimension of B.
[out]XThe solution matrix.
[in]ldxThe leading dimension of x.
[out]RCONDThe estimate of the reciprocal condition number of the matrix A after equilibration (if done).
[out]FERRThe estimates error bound for each solution.
[out]BERRThe componentwise relative backward error.
[out]WORK
[out]IWORK
[out]info

◆ dgetrf_()

void lapack::dgetrf_ ( int *  m,
int *  n,
double *  A,
int *  lda,
int *  ipiv,
int *  info 
)

LU decomposition.

Performs an LU decomposition of a genereal m-by-n matrix using partial pivoting with row interchanges. See LAPACK reference.

Parameters
[in]mThe number of rows of the matrix A.
[in]nThe number of column of the matrix A.
[out]AThe matrix A.
[in]ldaThe leading dimension of the matrix A.
[out]ipivInteger array to hold the pivoting indices.
[out]infoInteger indicating if operation was successful: 0 if success, otherwise failure.

Referenced by inv(), and ludcmp().

◆ dgetri_()

void lapack::dgetri_ ( int *  n,
double *  A,
int *  lda,
int *  ipiv,
double *  work,
int *  lwork,
int *  info 
)

Matrix inversion.

Inverts a given matrix using LU decompositions. See LAPACK reference.

Parameters
[in]nThe number of rows and columns of the matrix A.
[in,out]AThe matrix to be inverted.
[in]ldaThe leading dimension of A.
[in]ipivThe vector containing the pivoting indices returned from dgetrf_.
[in]workWork array for improved performance.
[in]lworkSize of the work array.
[out]infoInteger indicating if operation was successful: 0 if success, otherwise failure.

Referenced by inv().

◆ dgetrs_()

void lapack::dgetrs_ ( char *  trans,
int *  n,
int *  nrhs,
double *  A,
int *  lda,
int *  ipiv,
double *  b,
int *  ldb,
int *  info 
)

Solve linear system of equations.

Solves a linear system of equations of the form

trans(A) * x = b

using the QR decomposition obtained using dgetrf_.

Parameters
[in]transSpecifies the form of the system of equations:

trans = 'N' : trans(A) = A trans = 'T' : trans(A) = A^T trans = 'C' : trans(A) = conj( A^T )

Parameters
[in]nThe size of the system.
[in]nrhsThe number of right-hand sides.
[in]AThe matrix describing the system.
[in]ldaThe leading dimension of A.
[in]ipivThe pivot vector as returned by dgetrf_.
[in,out]bThe matrix containing the right-hand side vectors.
[in]ldbThe leading dimension of b.
[out]infoInteger indicating succes of the operation.

Referenced by lubacksub().

◆ ilaenv_()

void lapack::ilaenv_ ( int *  ispec,
char *  name,
char *  opts,
int *  n1,
int *  n2,
int *  n3,
int *  n4 
)

Optimal parameters for computation.

This function returns problem-dependent parameters for the computing environment. See LAPACK reference.

Parameters
[in,out]ispecOn input specifies which parameter to be returned and on output contains the requested parameter.
[in]nameName of the functions to be called.
[in]optsCharacter options to the function "Name"
[in]n1Problem dimension 1 of "Name".
[in]n2Problem dimension 2 of "Name".
[in]n3Problem dimension 3 of "Name".
[in]n4Problem dimension 4 of "Name".

◆ zgeev_()

void lapack::zgeev_ ( char *  jobvl,
char *  jobvr,
int *  n,
std::complex< double > *  A,
int *  lda,
std::complex< double > *  W,
std::complex< double > *  VL,
int *  ldvl,
std::complex< double > *  VR,
int *  ldvr,
std::complex< double > *  work,
int *  lwork,
double *  rwork,
int *  info 
)

Referenced by diagonalize().

◆ zgetrf_()

void lapack::zgetrf_ ( int *  m,
int *  n,
std::complex< double > *  A,
int *  lda,
int *  ipiv,
int *  info 
)

Referenced by inv(), and ComplexMatrix::inv().

◆ zgetri_()

void lapack::zgetri_ ( int *  n,
std::complex< double > *  A,
int *  lda,
int *  ipiv,
std::complex< double > *  work,
int *  lwork,
int *  info 
)

Referenced by inv(), and ComplexMatrix::inv().