46 #ifndef _PEXSI_LAPACK_HPP_
47 #define _PEXSI_LAPACK_HPP_
59 typedef std::complex<float> scomplex;
60 typedef std::complex<double> dcomplex;
67 void Potrf(
char uplo, Int n,
const float* A, Int lda );
68 void Potrf(
char uplo, Int n,
const double* A, Int lda );
69 void Potrf(
char uplo, Int n,
const scomplex* A, Int lda );
70 void Potrf(
char uplo, Int n,
const dcomplex* A, Int lda );
77 void Getrf( Int m, Int n,
float* A, Int lda, Int* p );
78 void Getrf( Int m, Int n,
double* A, Int lda, Int* p );
79 void Getrf( Int m, Int n, scomplex* A, Int lda, Int* p );
80 void Getrf( Int m, Int n, dcomplex* A, Int lda, Int* p );
88 ( Int itype,
char uplo,
89 Int n,
float* A, Int lda,
const float* B, Int ldb );
91 ( Int itype,
char uplo,
92 Int n,
double* A, Int lda,
const double* B, Int ldb );
94 ( Int itype,
char uplo,
95 Int n, scomplex* A, Int lda,
const scomplex* B, Int ldb );
97 ( Int itype,
char uplo,
98 Int n, dcomplex* A, Int lda,
const dcomplex* B, Int ldb );
106 (
char jobz,
char uplo, Int n,
double* A, Int lda,
double* eigs );
114 (
int itype,
char jobz,
char uplo, Int n,
double* A, Int lda,
115 double* B, Int ldb,
double* eigs );
124 (
char uplo,
char diag, Int n,
const float* A, Int lda );
126 (
char uplo,
char diag, Int n,
const double* A, Int lda );
128 (
char uplo,
char diag, Int n,
const scomplex* A, Int lda );
130 (
char uplo,
char diag, Int n,
const dcomplex* A, Int lda );
137 void DivideAndConquerSVD
138 ( Int m, Int n,
float* A, Int lda,
139 float* s,
float* U, Int ldu,
float* VTrans, Int ldvt );
140 void DivideAndConquerSVD
141 ( Int m, Int n,
double* A, Int lda,
142 double* s,
double* U, Int ldu,
double* VTrans, Int ldvt );
143 void DivideAndConquerSVD
144 ( Int m, Int n, scomplex* A, Int lda,
145 float* s, scomplex* U, Int ldu, scomplex* VAdj, Int ldva );
146 void DivideAndConquerSVD
147 ( Int m, Int n, dcomplex* A, Int lda,
148 double* s, dcomplex* U, Int ldu, dcomplex* VAdj, Int ldva );
155 ( Int m, Int n,
float* A, Int lda,
156 float* s,
float* U, Int ldu,
float* VTrans, Int ldvt );
158 ( Int m, Int n,
double* A, Int lda,
159 double* s,
double* U, Int ldu,
double* VTrans, Int ldvt );
161 ( Int m, Int n, scomplex* A, Int lda,
162 float* s, scomplex* U, Int ldu, scomplex* VAdj, Int ldva );
164 ( Int m, Int n, dcomplex* A, Int lda,
165 double* s, dcomplex* U, Int ldu, dcomplex* VAdj, Int ldva );
172 void SingularValues( Int m, Int n,
float* A, Int lda,
float* s );
173 void SingularValues( Int m, Int n,
double* A, Int lda,
double* s );
174 void SingularValues( Int m, Int n, scomplex* A, Int lda,
float* s );
175 void SingularValues( Int m, Int n, dcomplex* A, Int lda,
double* s );
182 (
char uplo, Int n, Int numColsVTrans, Int numRowsU,
183 float* d,
float* e,
float* VTrans, Int ldVTrans,
float* U, Int ldU );
185 (
char uplo, Int n, Int numColsVTrans, Int numRowsU,
186 double* d,
double* e,
double* VTrans, Int ldVTrans,
double* U, Int ldU );
188 (
char uplo, Int n, Int numColsVAdj, Int numRowsU,
189 float* d,
float* e, scomplex* VAdj, Int ldVAdj, scomplex* U, Int ldU );
191 (
char uplo, Int n, Int numColsVAdj, Int numRowsU,
192 double* d,
double* e, dcomplex* VAdj, Int ldVAdj, dcomplex* U, Int ldU );
197 void SVDLeastSquare( Int m, Int n, Int nrhs,
float * A, Int lda,
198 float * B, Int ldb,
float * S,
float rcond,
200 void SVDLeastSquare( Int m, Int n, Int nrhs,
double * A, Int lda,
201 double * B, Int ldb,
double * S,
double rcond,
203 void SVDLeastSquare( Int m, Int n, Int nrhs, scomplex * A, Int lda,
204 scomplex * B, Int ldb,
float * S,
float rcond,
206 void SVDLeastSquare( Int m, Int n, Int nrhs, dcomplex * A, Int lda,
207 dcomplex * B, Int ldb,
double * S,
double rcond,
214 void Lacpy(
char uplo, Int m, Int n,
const double* A, Int lda,
215 double* B, Int ldb );
217 void Lacpy(
char uplo, Int m, Int n,
const dcomplex* A, Int lda,
218 dcomplex* B, Int ldb );
225 void Getri ( Int n,
double* A, Int lda,
const Int* ipiv );
227 void Getri ( Int n, dcomplex* A, Int lda,
const Int* ipiv );
234 double Lange (
char norm, Int m, Int n,
float * A, Int lda,
float* work);
235 double Lange (
char norm, Int m, Int n,
double * A, Int lda,
double* work);
236 double Lange (
char norm, Int m, Int n, scomplex * A, Int lda, scomplex* work);
237 double Lange (
char norm, Int m, Int n, dcomplex * A, Int lda, dcomplex* work);
242 #endif //_PEXSI_LAPACK_HPP_