42 #ifndef _TEUCHOS_LAPACK_HPP_
43 #define _TEUCHOS_LAPACK_HPP_
89 struct UndefinedLAPACKRoutine
92 static inline T notDefined() {
return T::LAPACK_routine_not_defined_for_this_type(); }
95 std::complex<double> convert_Fortran_complex_to_CXX_complex(_Complex
double val);
97 std::complex<float> convert_Fortran_complex_to_CXX_complex(_Complex
float val);
99 template<
typename OrdinalType,
typename ScalarType>
123 void PTTRF(
const OrdinalType& n, MagnitudeType* d, ScalarType* e, OrdinalType* info)
const;
126 void PTTRS(
const OrdinalType& n,
const OrdinalType& nrhs,
const MagnitudeType* d,
const ScalarType* e, ScalarType* B,
const OrdinalType& ldb, OrdinalType* info)
const;
129 void POTRF(
const char& UPLO,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, OrdinalType* info)
const;
132 void POTRS(
const char& UPLO,
const OrdinalType& n,
const OrdinalType& nrhs,
const ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, OrdinalType* info)
const;
135 void POTRI(
const char& UPLO,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, OrdinalType* info)
const;
138 void POCON(
const char& UPLO,
const OrdinalType& n,
const ScalarType* A,
const OrdinalType& lda,
const ScalarType& anorm, ScalarType* rcond, ScalarType* WORK, OrdinalType* IWORK, OrdinalType* info)
const;
141 void POSV(
const char& UPLO,
const OrdinalType& n,
const OrdinalType& nrhs, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, OrdinalType* info)
const;
144 void POEQU(
const OrdinalType& n,
const ScalarType* A,
const OrdinalType& lda, MagnitudeType* S, MagnitudeType* scond, MagnitudeType* amax, OrdinalType* info)
const;
147 void PORFS(
const char& UPLO,
const OrdinalType& n,
const OrdinalType& nrhs,
const ScalarType* A,
const OrdinalType& lda,
const ScalarType* AF,
const OrdinalType& ldaf,
const ScalarType* B,
const OrdinalType& ldb, ScalarType* X,
const OrdinalType& ldx, ScalarType* FERR, ScalarType* BERR, ScalarType* WORK, OrdinalType* IWORK, OrdinalType* info)
const;
150 void POSVX(
const char& FACT,
const char& UPLO,
const OrdinalType& n,
const OrdinalType& nrhs, ScalarType* A,
const OrdinalType& lda, ScalarType* AF,
const OrdinalType& ldaf,
char* EQUED, ScalarType* S, ScalarType* B,
const OrdinalType& ldb, ScalarType* X,
const OrdinalType& ldx, ScalarType* rcond, ScalarType* FERR, ScalarType* BERR, ScalarType* WORK, OrdinalType* IWORK, OrdinalType* info)
const;
157 void GELS(
const char&
TRANS,
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& nrhs, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const;
192 void GELSS(
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& nrhs, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, MagnitudeType* S,
const MagnitudeType rcond, OrdinalType* rank, ScalarType* WORK,
const OrdinalType& lwork, MagnitudeType* RWORK, OrdinalType* info)
const;
195 void GELSS(
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& nrhs, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, ScalarType* S,
const ScalarType& rcond, OrdinalType* rank, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const;
198 void GGLSE(
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& p, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, ScalarType* C, ScalarType* D, ScalarType* X, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const;
201 void GEQRF (
const OrdinalType& m,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* TAU, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const;
204 void GEQR2 (
const OrdinalType& m,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* TAU, ScalarType* WORK, OrdinalType*
const info)
const;
207 void GETRF(
const OrdinalType& m,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, OrdinalType* IPIV, OrdinalType* info)
const;
210 void GETRS(
const char&
TRANS,
const OrdinalType& n,
const OrdinalType& nrhs,
const ScalarType* A,
const OrdinalType& lda,
const OrdinalType* IPIV, ScalarType* B,
const OrdinalType& ldb, OrdinalType* info)
const;
213 void LASCL(
const char& TYPE,
const OrdinalType& kl,
const OrdinalType& ku,
const MagnitudeType cfrom,
const MagnitudeType cto,
const OrdinalType& m,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, OrdinalType* info)
const;
216 void GEQP3(
const OrdinalType& m,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, OrdinalType* jpvt, ScalarType* TAU, ScalarType* WORK,
const OrdinalType& lwork, MagnitudeType* RWORK, OrdinalType* info )
const;
219 void LASWP (
const OrdinalType& N, ScalarType* A,
const OrdinalType& LDA,
const OrdinalType& K1,
const OrdinalType& K2,
const OrdinalType* IPIV,
const OrdinalType& INCX)
const;
222 void GBTRF(
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& kl,
const OrdinalType& ku, ScalarType* A,
const OrdinalType& lda, OrdinalType* IPIV, OrdinalType* info)
const;
225 void GBTRS(
const char&
TRANS,
const OrdinalType& n,
const OrdinalType& kl,
const OrdinalType& ku,
const OrdinalType& nrhs,
const ScalarType* A,
const OrdinalType& lda,
const OrdinalType* IPIV, ScalarType* B,
const OrdinalType& ldb, OrdinalType* info)
const;
228 void GTTRF(
const OrdinalType& n, ScalarType* dl, ScalarType* d, ScalarType* du, ScalarType* du2, OrdinalType* IPIV, OrdinalType* info)
const;
231 void GTTRS(
const char&
TRANS,
const OrdinalType& n,
const OrdinalType& nrhs,
const ScalarType* dl,
const ScalarType* d,
const ScalarType* du,
const ScalarType* du2,
const OrdinalType* IPIV, ScalarType* B,
const OrdinalType& ldb, OrdinalType* info)
const;
234 void GETRI(
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda,
const OrdinalType* IPIV, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const;
240 void LATRS (
const char& UPLO,
const char&
TRANS,
const char& DIAG,
const char& NORMIN,
const OrdinalType& N,
const ScalarType* A,
const OrdinalType& LDA, ScalarType* X, MagnitudeType* SCALE, MagnitudeType* CNORM, OrdinalType* INFO)
const;
243 void GECON(
const char& NORM,
const OrdinalType& n,
const ScalarType* A,
const OrdinalType& lda,
const ScalarType& anorm, ScalarType* rcond, ScalarType* WORK, OrdinalType* IWORK, OrdinalType* info)
const;
246 void GBCON(
const char& NORM,
const OrdinalType& n,
const OrdinalType& kl,
const OrdinalType& ku,
const ScalarType* A,
const OrdinalType& lda,
const OrdinalType* IPIV,
const ScalarType& anorm, ScalarType* rcond, ScalarType* WORK, OrdinalType* IWORK, OrdinalType* info)
const;
252 void GESV(
const OrdinalType& n,
const OrdinalType& nrhs, ScalarType* A,
const OrdinalType& lda, OrdinalType* IPIV, ScalarType* B,
const OrdinalType& ldb, OrdinalType* info)
const;
255 void GEEQU(
const OrdinalType& m,
const OrdinalType& n,
const ScalarType* A,
const OrdinalType& lda, ScalarType* R, ScalarType* C, ScalarType* rowcond, ScalarType* colcond, ScalarType* amax, OrdinalType* info)
const;
258 void GERFS(
const char&
TRANS,
const OrdinalType& n,
const OrdinalType& nrhs,
const ScalarType* A,
const OrdinalType& lda,
const ScalarType* AF,
const OrdinalType& ldaf,
const OrdinalType* IPIV,
const ScalarType* B,
const OrdinalType& ldb, ScalarType* X,
const OrdinalType& ldx, ScalarType* FERR, ScalarType* BERR, ScalarType* WORK, OrdinalType* IWORK, OrdinalType* info)
const;
261 void GBEQU(
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& kl,
const OrdinalType& ku,
const ScalarType* A,
const OrdinalType& lda, MagnitudeType* R, MagnitudeType* C, MagnitudeType* rowcond, MagnitudeType* colcond, MagnitudeType* amax, OrdinalType* info)
const;
264 void GBRFS(
const char&
TRANS,
const OrdinalType& n,
const OrdinalType& kl,
const OrdinalType& ku,
const OrdinalType& nrhs,
const ScalarType* A,
const OrdinalType& lda,
const ScalarType* AF,
const OrdinalType& ldaf,
const OrdinalType* IPIV,
const ScalarType* B,
const OrdinalType& ldb, ScalarType* X,
const OrdinalType& ldx, ScalarType* FERR, ScalarType* BERR, ScalarType* WORK, OrdinalType* IWORK, OrdinalType* info)
const;
268 void GESVX(
const char& FACT,
const char&
TRANS,
const OrdinalType& n,
const OrdinalType& nrhs, ScalarType* A,
const OrdinalType& lda, ScalarType* AF,
const OrdinalType& ldaf, OrdinalType* IPIV,
char* EQUED, ScalarType* R, ScalarType* C, ScalarType* B,
const OrdinalType& ldb, ScalarType* X,
const OrdinalType& ldx, ScalarType* rcond, ScalarType* FERR, ScalarType* BERR, ScalarType* WORK, OrdinalType* IWORK, OrdinalType* info)
const;
273 void SYTRD(
const char& UPLO,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* D, ScalarType* E, ScalarType* TAU, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const;
276 void GEHRD(
const OrdinalType& n,
const OrdinalType& ilo,
const OrdinalType& ihi, ScalarType* A,
const OrdinalType& lda, ScalarType* TAU, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const;
279 void TRTRS(
const char& UPLO,
const char&
TRANS,
const char& DIAG,
const OrdinalType& n,
const OrdinalType& nrhs,
const ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, OrdinalType* info)
const;
282 void TRTRI(
const char& UPLO,
const char& DIAG,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, OrdinalType* info)
const;
290 void SPEV(
const char& JOBZ,
const char& UPLO,
const OrdinalType& n, ScalarType* AP, ScalarType* W, ScalarType* Z,
const OrdinalType& ldz, ScalarType* WORK, OrdinalType* info)
const;
295 void SYEV(
const char& JOBZ,
const char& UPLO,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* W, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const;
300 void SYGV(
const OrdinalType& itype,
const char& JOBZ,
const char& UPLO,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, ScalarType* W, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const;
305 void HEEV(
const char& JOBZ,
const char& UPLO,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, MagnitudeType* W, ScalarType* WORK,
const OrdinalType& lwork, MagnitudeType* RWORK, OrdinalType* info)
const;
310 void HEGV(
const OrdinalType& itype,
const char& JOBZ,
const char& UPLO,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, MagnitudeType* W, ScalarType* WORK,
const OrdinalType& lwork, MagnitudeType *RWORK, OrdinalType* info)
const;
313 void STEQR(
const char& COMPZ,
const OrdinalType& n, MagnitudeType* D, MagnitudeType* E, ScalarType* Z,
const OrdinalType& ldz, MagnitudeType* WORK, OrdinalType* info)
const;
316 void PTEQR(
const char& COMPZ,
const OrdinalType& n, MagnitudeType* D, MagnitudeType* E, ScalarType* Z,
const OrdinalType& ldz, MagnitudeType* WORK, OrdinalType* info)
const;
321 void HSEQR(
const char& JOB,
const char& COMPZ,
const OrdinalType& n,
const OrdinalType& ilo,
const OrdinalType& ihi, ScalarType* H,
const OrdinalType& ldh, ScalarType* WR, ScalarType* WI, ScalarType* Z,
const OrdinalType& ldz, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const;
327 void GEES(
const char& JOBVS,
const char& SORT, OrdinalType& (*ptr2func)(ScalarType*, ScalarType*),
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, OrdinalType* sdim, ScalarType* WR, ScalarType* WI, ScalarType* VS,
const OrdinalType& ldvs, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* BWORK, OrdinalType* info)
const;
332 void GEES(
const char& JOBVS,
const char& SORT, OrdinalType& (*ptr2func)(ScalarType*),
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, OrdinalType* sdim, ScalarType* W, ScalarType* VS,
const OrdinalType& ldvs, ScalarType* WORK,
const OrdinalType& lwork, MagnitudeType* RWORK, OrdinalType* BWORK, OrdinalType* info)
const;
337 void GEES(
const char& JOBVS,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, OrdinalType* sdim, MagnitudeType* WR, MagnitudeType* WI, ScalarType* VS,
const OrdinalType& ldvs, ScalarType* WORK,
const OrdinalType& lwork, MagnitudeType* RWORK, OrdinalType* BWORK, OrdinalType* info)
const;
344 void GEEV(
const char& JOBVL,
const char& JOBVR,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, MagnitudeType* WR, MagnitudeType* WI, ScalarType* VL,
const OrdinalType& ldvl, ScalarType* VR,
const OrdinalType& ldvr, ScalarType* WORK,
const OrdinalType& lwork, MagnitudeType* RWORK, OrdinalType* info)
const;
350 void GEEVX(
const char& BALANC,
const char& JOBVL,
const char& JOBVR,
const char& SENSE,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* WR, ScalarType* WI, ScalarType* VL,
const OrdinalType& ldvl, ScalarType* VR,
const OrdinalType& ldvr, OrdinalType* ilo, OrdinalType* ihi, MagnitudeType* SCALE, MagnitudeType* abnrm, MagnitudeType* RCONDE, MagnitudeType* RCONDV, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* IWORK, OrdinalType* info)
const;
356 void GGEVX(
const char& BALANC,
const char& JOBVL,
const char& JOBVR,
const char& SENSE,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, MagnitudeType* ALPHAR, MagnitudeType* ALPHAI, ScalarType* BETA, ScalarType* VL,
const OrdinalType& ldvl, ScalarType* VR,
const OrdinalType& ldvr, OrdinalType* ilo, OrdinalType* ihi, MagnitudeType* lscale, MagnitudeType* rscale, MagnitudeType* abnrm, MagnitudeType* bbnrm, MagnitudeType* RCONDE, MagnitudeType* RCONDV, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* IWORK, OrdinalType* BWORK, OrdinalType* info)
const;
361 void GGEV(
const char& JOBVL,
const char& JOBVR,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, MagnitudeType *ALPHAR, MagnitudeType *ALPHAI, ScalarType* BETA, ScalarType* VL,
const OrdinalType& ldvl, ScalarType* VR,
const OrdinalType& ldvr, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const;
367 void TRSEN(
const char& JOB,
const char& COMPQ,
const OrdinalType* SELECT,
const OrdinalType& n, ScalarType* T,
const OrdinalType& ldt, ScalarType* Q,
const OrdinalType& ldq, MagnitudeType *WR, MagnitudeType *WI, OrdinalType* M, ScalarType* S, MagnitudeType *SEP, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* IWORK,
const OrdinalType& liwork, OrdinalType* info )
const;
372 void TGSEN(
const OrdinalType& ijob,
const OrdinalType& wantq,
const OrdinalType& wantz,
const OrdinalType* SELECT,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, MagnitudeType *ALPHAR, MagnitudeType *ALPHAI, MagnitudeType *BETA, ScalarType* Q,
const OrdinalType& ldq, ScalarType* Z,
const OrdinalType& ldz, OrdinalType* M, MagnitudeType *PL, MagnitudeType *PR, MagnitudeType *DIF, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* IWORK,
const OrdinalType& liwork, OrdinalType* info )
const;
377 void GGES(
const char& JOBVL,
const char& JOBVR,
const char& SORT, OrdinalType& (*ptr2func)(ScalarType* , ScalarType* , ScalarType* ),
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, OrdinalType* sdim, MagnitudeType *ALPHAR, MagnitudeType *ALPHAI, MagnitudeType *BETA, ScalarType* VL,
const OrdinalType& ldvl, ScalarType* VR,
const OrdinalType& ldvr, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* BWORK, OrdinalType* info )
const;
384 void GESVD(
const char& JOBU,
const char& JOBVT,
const OrdinalType& m,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, MagnitudeType* S, ScalarType* U,
const OrdinalType& ldu, ScalarType* V,
const OrdinalType& ldv, ScalarType* WORK,
const OrdinalType& lwork, MagnitudeType* RWORK, OrdinalType* info)
const;
401 void ORMQR(
const char& SIDE,
const char&
TRANS,
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& k,
const ScalarType* A,
const OrdinalType& lda,
const ScalarType* TAU, ScalarType* C,
const OrdinalType& ldc, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const;
408 void ORM2R(
const char& SIDE,
const char&
TRANS,
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& k,
const ScalarType* A,
const OrdinalType& lda,
const ScalarType* TAU, ScalarType* C,
const OrdinalType& ldc, ScalarType* WORK, OrdinalType*
const info)
const;
418 void UNMQR(
const char& SIDE,
const char&
TRANS,
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& k,
const ScalarType* A,
const OrdinalType& lda,
const ScalarType* TAU, ScalarType* C,
const OrdinalType& ldc, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const;
424 void UNM2R(
const char& SIDE,
const char&
TRANS,
const OrdinalType& M,
const OrdinalType& N,
const OrdinalType& K,
const ScalarType* A,
const OrdinalType& LDA,
const ScalarType* TAU, ScalarType* C,
const OrdinalType& LDC, ScalarType* WORK, OrdinalType*
const INFO)
const;
435 void ORGQR(
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& k, ScalarType* A,
const OrdinalType& lda,
const ScalarType* TAU, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const;
445 void UNGQR(
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& k, ScalarType* A,
const OrdinalType& lda,
const ScalarType* TAU, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const;
450 void ORGHR(
const OrdinalType& n,
const OrdinalType& ilo,
const OrdinalType& ihi, ScalarType* A,
const OrdinalType& lda,
const ScalarType* TAU, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const;
455 void ORMHR(
const char& SIDE,
const char&
TRANS,
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& ilo,
const OrdinalType& ihi,
const ScalarType* A,
const OrdinalType& lda,
const ScalarType* TAU, ScalarType* C,
const OrdinalType& ldc, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const;
463 void TREVC(
const char& SIDE,
const char& HOWMNY, OrdinalType* select,
const OrdinalType& n,
const ScalarType* T,
const OrdinalType& ldt, ScalarType* VL,
const OrdinalType& ldvl, ScalarType* VR,
const OrdinalType& ldvr,
const OrdinalType& mm, OrdinalType* m, ScalarType* WORK, OrdinalType* info)
const;
468 void TREVC(
const char& SIDE,
const OrdinalType& n,
const ScalarType* T,
const OrdinalType& ldt, ScalarType* VL,
const OrdinalType& ldvl, ScalarType* VR,
const OrdinalType& ldvr,
const OrdinalType& mm, OrdinalType* m, ScalarType* WORK, MagnitudeType* RWORK, OrdinalType* info)
const;
473 void TREXC(
const char& COMPQ,
const OrdinalType& n, ScalarType* T,
const OrdinalType& ldt, ScalarType* Q,
const OrdinalType& ldq, OrdinalType* ifst, OrdinalType* ilst, ScalarType* WORK, OrdinalType* info)
const;
478 void TGEVC(
const char& SIDE,
const char& HOWMNY,
const OrdinalType* SELECT,
const OrdinalType& n,
const ScalarType* S,
const OrdinalType& lds,
const ScalarType* P,
const OrdinalType& ldp, ScalarType* VL,
const OrdinalType& ldvl, ScalarType* VR,
const OrdinalType& ldvr,
const OrdinalType& mm, OrdinalType* M, ScalarType* WORK, OrdinalType* info)
const;
487 void LARTG(
const ScalarType& f,
const ScalarType& g, MagnitudeType* c, ScalarType* s, ScalarType* r )
const;
490 void LARFG(
const OrdinalType& n, ScalarType* alpha, ScalarType* x,
const OrdinalType& incx, ScalarType* tau )
const;
499 void GEBAL(
const char& JOBZ,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, OrdinalType* ilo, OrdinalType* ihi, MagnitudeType* scale, OrdinalType* info)
const;
502 void GEBAK(
const char& JOBZ,
const char& SIDE,
const OrdinalType& n,
const OrdinalType& ilo,
const OrdinalType& ihi,
const MagnitudeType* scale ,
const OrdinalType& m, ScalarType* V,
const OrdinalType& ldv, OrdinalType* info)
const;
508 ScalarType
LARND(
const OrdinalType& idist, OrdinalType* seed )
const;
512 void LARNV(
const OrdinalType& idist, OrdinalType* seed,
const OrdinalType& n, ScalarType* v )
const;
520 ScalarType
LAMCH(
const char& CMACH)
const;
526 OrdinalType
ILAENV(
const OrdinalType& ispec,
const std::string& NAME,
const std::string& OPTS,
const OrdinalType& N1 = -1,
const OrdinalType& N2 = -1,
const OrdinalType& N3 = -1,
const OrdinalType& N4 = -1 )
const;
534 ScalarType
LAPY2(
const ScalarType& x,
const ScalarType& y)
const;
543 template<
typename OrdinalType,
typename ScalarType>
546 UndefinedLAPACKRoutine<ScalarType>::notDefined();
549 template<
typename OrdinalType,
typename ScalarType>
550 void LAPACK<OrdinalType, ScalarType>::PTTRS(
const OrdinalType& n,
const OrdinalType& nrhs,
const MagnitudeType* d,
const ScalarType* e, ScalarType* B,
const OrdinalType& ldb, OrdinalType* info)
const
552 UndefinedLAPACKRoutine<ScalarType>::notDefined();
555 template<
typename OrdinalType,
typename ScalarType>
558 UndefinedLAPACKRoutine<ScalarType>::notDefined();
561 template<
typename OrdinalType,
typename ScalarType>
562 void LAPACK<OrdinalType, ScalarType>::POTRS(
const char& UPLO,
const OrdinalType& n,
const OrdinalType& nrhs,
const ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, OrdinalType* info)
const
564 UndefinedLAPACKRoutine<ScalarType>::notDefined();
567 template<
typename OrdinalType,
typename ScalarType>
570 UndefinedLAPACKRoutine<ScalarType>::notDefined();
573 template<
typename OrdinalType,
typename ScalarType>
574 void LAPACK<OrdinalType, ScalarType>::POCON(
const char& UPLO,
const OrdinalType& n,
const ScalarType* A,
const OrdinalType& lda,
const ScalarType& anorm, ScalarType* rcond, ScalarType* WORK, OrdinalType* IWORK, OrdinalType* info)
const
576 UndefinedLAPACKRoutine<ScalarType>::notDefined();
579 template<
typename OrdinalType,
typename ScalarType>
580 void LAPACK<OrdinalType, ScalarType>::POSV(
const char& UPLO,
const OrdinalType& n,
const OrdinalType& nrhs, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, OrdinalType* info)
const
582 UndefinedLAPACKRoutine<ScalarType>::notDefined();
585 template<
typename OrdinalType,
typename ScalarType>
586 void LAPACK<OrdinalType, ScalarType>::POEQU(
const OrdinalType& n,
const ScalarType* A,
const OrdinalType& lda, MagnitudeType* S, MagnitudeType* scond, MagnitudeType* amax, OrdinalType* info)
const
592 }
else if (lda < TEUCHOS_MAX(1, n)) {
613 MagnitudeType smin = S[0];
615 for (OrdinalType i=0; i<n; ++i) {
617 smin = TEUCHOS_MIN( smin, S[i] );
618 *amax = TEUCHOS_MAX( *amax, S[i] );
623 for (OrdinalType i=0; i<n; ++i) {
629 for (OrdinalType i=0; i<n; ++i) {
637 template<
typename OrdinalType,
typename ScalarType>
638 void LAPACK<OrdinalType, ScalarType>::PORFS(
const char& UPLO,
const OrdinalType& n,
const OrdinalType& nrhs,
const ScalarType* A,
const OrdinalType& lda,
const ScalarType* AF,
const OrdinalType& ldaf,
const ScalarType* B,
const OrdinalType& ldb, ScalarType* X,
const OrdinalType& ldx, ScalarType* FERR, ScalarType* BERR, ScalarType* WORK, OrdinalType* IWORK, OrdinalType* info)
const
640 UndefinedLAPACKRoutine<ScalarType>::notDefined();
643 template<
typename OrdinalType,
typename ScalarType>
644 void LAPACK<OrdinalType, ScalarType>::POSVX(
const char& FACT,
const char& UPLO,
const OrdinalType& n,
const OrdinalType& nrhs, ScalarType* A,
const OrdinalType& lda, ScalarType* AF,
const OrdinalType& ldaf,
char* EQUED, ScalarType* S, ScalarType* B,
const OrdinalType& ldb, ScalarType* X,
const OrdinalType& ldx, ScalarType* rcond, ScalarType* FERR, ScalarType* BERR, ScalarType* WORK, OrdinalType* IWORK, OrdinalType* info)
const
646 UndefinedLAPACKRoutine<ScalarType>::notDefined();
649 template<
typename OrdinalType,
typename ScalarType>
650 void LAPACK<OrdinalType,ScalarType>::GELS(
const char&
TRANS,
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& nrhs, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const
652 UndefinedLAPACKRoutine<ScalarType>::notDefined();
655 template<
typename OrdinalType,
typename ScalarType>
656 void LAPACK<OrdinalType, ScalarType>::GELSS(
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& nrhs, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, MagnitudeType* S,
const MagnitudeType rcond, OrdinalType* rank, ScalarType* WORK,
const OrdinalType& lwork, MagnitudeType* RWORK, OrdinalType* info)
const
658 UndefinedLAPACKRoutine<ScalarType>::notDefined();
661 template<
typename OrdinalType,
typename ScalarType>
662 void LAPACK<OrdinalType,ScalarType>::GELSS(
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& nrhs, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, ScalarType* S,
const ScalarType& rcond, OrdinalType* rank, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const
664 UndefinedLAPACKRoutine<ScalarType>::notDefined();
667 template<
typename OrdinalType,
typename ScalarType>
668 void LAPACK<OrdinalType,ScalarType>::GGLSE(
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& p, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, ScalarType* C, ScalarType* D, ScalarType* X, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const
670 UndefinedLAPACKRoutine<ScalarType>::notDefined();
673 template<
typename OrdinalType,
typename ScalarType>
674 void LAPACK<OrdinalType,ScalarType>::GEQRF(
const OrdinalType& m,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* TAU, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const
676 UndefinedLAPACKRoutine<ScalarType>::notDefined();
679 template<
typename OrdinalType,
typename ScalarType>
680 void LAPACK<OrdinalType,ScalarType>::GEQR2 (
const OrdinalType& m,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* TAU, ScalarType* WORK, OrdinalType*
const info)
const
682 UndefinedLAPACKRoutine<ScalarType>::notDefined();
685 template<
typename OrdinalType,
typename ScalarType>
688 UndefinedLAPACKRoutine<ScalarType>::notDefined();
691 template<
typename OrdinalType,
typename ScalarType>
692 void LAPACK<OrdinalType,ScalarType>::GETRS(
const char&
TRANS,
const OrdinalType& n,
const OrdinalType& nrhs,
const ScalarType* A,
const OrdinalType& lda,
const OrdinalType* IPIV, ScalarType* B,
const OrdinalType& ldb, OrdinalType* info)
const
694 UndefinedLAPACKRoutine<ScalarType>::notDefined();
697 template<
typename OrdinalType,
typename ScalarType>
698 void LAPACK<OrdinalType,ScalarType>::LASCL(
const char& TYPE,
const OrdinalType& kl,
const OrdinalType& ku,
const MagnitudeType cfrom,
const MagnitudeType cto,
const OrdinalType& m,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, OrdinalType* info)
const
714 MagnitudeType cfromc = cfrom;
715 MagnitudeType ctoc = cto;
716 MagnitudeType cfrom1;
721 cfrom1 = cfromc*smlnum;
722 if (cfrom1 == cfromc) {
728 cto1 = ctoc / bignum;
748 for (j=0; j<n; j++) {
750 for (i=0; i<m; i++) { *ptr = mul * (*ptr); ptr++; }
756 template<
typename OrdinalType,
typename ScalarType>
757 void LAPACK<OrdinalType,ScalarType>::GEQP3 (
const OrdinalType& m,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, OrdinalType* jpvt, ScalarType* TAU, ScalarType* WORK,
const OrdinalType& lwork, MagnitudeType* RWORK, OrdinalType* info)
const
759 UndefinedLAPACKRoutine<ScalarType>::notDefined();
762 template<
typename OrdinalType,
typename ScalarType>
763 void LAPACK<OrdinalType, ScalarType>::LASWP (
const OrdinalType& N, ScalarType* A,
const OrdinalType& LDA,
const OrdinalType& K1,
const OrdinalType& K2,
const OrdinalType* IPIV,
const OrdinalType& INCX)
const
765 UndefinedLAPACKRoutine<ScalarType>::notDefined();
768 template<
typename OrdinalType,
typename ScalarType>
769 void LAPACK<OrdinalType,ScalarType>::GBTRF(
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& kl,
const OrdinalType& ku, ScalarType* A,
const OrdinalType& lda, OrdinalType* IPIV, OrdinalType* info)
const
771 UndefinedLAPACKRoutine<ScalarType>::notDefined();
774 template<
typename OrdinalType,
typename ScalarType>
775 void LAPACK<OrdinalType,ScalarType>::GBTRS(
const char&
TRANS,
const OrdinalType& n,
const OrdinalType& kl,
const OrdinalType& ku,
const OrdinalType& nrhs,
const ScalarType* A,
const OrdinalType& lda,
const OrdinalType* IPIV, ScalarType* B,
const OrdinalType& ldb, OrdinalType* info)
const
777 UndefinedLAPACKRoutine<ScalarType>::notDefined();
780 template<
typename OrdinalType,
typename ScalarType>
783 UndefinedLAPACKRoutine<ScalarType>::notDefined();
786 template<
typename OrdinalType,
typename ScalarType>
787 void LAPACK<OrdinalType,ScalarType>::GTTRS(
const char&
TRANS,
const OrdinalType& n,
const OrdinalType& nrhs,
const ScalarType* dl,
const ScalarType* d,
const ScalarType* du,
const ScalarType* du2,
const OrdinalType* IPIV, ScalarType* B,
const OrdinalType& ldb, OrdinalType* info)
const
789 UndefinedLAPACKRoutine<ScalarType>::notDefined();
792 template<
typename OrdinalType,
typename ScalarType>
793 void LAPACK<OrdinalType,ScalarType>::GETRI(
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda,
const OrdinalType* IPIV, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const
795 UndefinedLAPACKRoutine<ScalarType>::notDefined();
798 template<
typename OrdinalType,
typename ScalarType>
799 void LAPACK<OrdinalType,ScalarType>::LATRS (
const char& UPLO,
const char&
TRANS,
const char& DIAG,
const char& NORMIN,
const OrdinalType& N,
const ScalarType* A,
const OrdinalType& LDA, ScalarType* X, MagnitudeType* SCALE, MagnitudeType* CNORM, OrdinalType* INFO)
const
801 UndefinedLAPACKRoutine<ScalarType>::notDefined();
804 template<
typename OrdinalType,
typename ScalarType>
805 void LAPACK<OrdinalType,ScalarType>::GECON(
const char& NORM,
const OrdinalType& n,
const ScalarType* A,
const OrdinalType& lda,
const ScalarType& anorm, ScalarType* rcond, ScalarType* WORK, OrdinalType* IWORK, OrdinalType* info)
const
807 UndefinedLAPACKRoutine<ScalarType>::notDefined();
810 template<
typename OrdinalType,
typename ScalarType>
811 void LAPACK<OrdinalType,ScalarType>::GBCON(
const char& NORM,
const OrdinalType& n,
const OrdinalType& kl,
const OrdinalType& ku,
const ScalarType* A,
const OrdinalType& lda,
const OrdinalType* IPIV,
const ScalarType& anorm, ScalarType* rcond, ScalarType* WORK, OrdinalType* IWORK, OrdinalType* info)
const
813 UndefinedLAPACKRoutine<ScalarType>::notDefined();
816 template<
typename OrdinalType,
typename ScalarType>
819 UndefinedLAPACKRoutine<ScalarType>::notDefined();
822 template<
typename OrdinalType,
typename ScalarType>
823 void LAPACK<OrdinalType,ScalarType>::GESV(
const OrdinalType& n,
const OrdinalType& nrhs, ScalarType* A,
const OrdinalType& lda, OrdinalType* IPIV, ScalarType* B,
const OrdinalType& ldb, OrdinalType* info)
const
825 UndefinedLAPACKRoutine<ScalarType>::notDefined();
828 template<
typename OrdinalType,
typename ScalarType>
829 void LAPACK<OrdinalType,ScalarType>::GEEQU(
const OrdinalType& m,
const OrdinalType& n,
const ScalarType* A,
const OrdinalType& lda, ScalarType* R, ScalarType* C, ScalarType* rowcond, ScalarType* colcond, ScalarType* amax, OrdinalType* info)
const
838 }
else if (lda < TEUCHOS_MAX(1, m)) {
851 if (m == 0 || n == 0) {
863 for (OrdinalType i=0; i<m; i++) {
868 for (OrdinalType j=0; j<n; j++) {
869 for (OrdinalType i=0; i<m; i++) {
875 MagnitudeType rcmin = bignum;
876 MagnitudeType rcmax = mZero;
877 for (OrdinalType i=0; i<m; i++) {
878 rcmax = TEUCHOS_MAX( rcmax, R[i] );
879 rcmin = TEUCHOS_MIN( rcmin, R[i] );
883 if (rcmin == mZero) {
885 for (OrdinalType i=0; i<m; i++) {
891 for (OrdinalType i=0; i<m; i++) {
892 R[i] = mOne / TEUCHOS_MIN( TEUCHOS_MAX( R[i], smlnum ), bignum );
895 *rowcond = TEUCHOS_MAX( rcmin, smlnum ) / TEUCHOS_MIN( rcmax, bignum );
899 for (OrdinalType j=0; j<n; j++) {
904 for (OrdinalType j=0; j<n; j++) {
905 for (OrdinalType i=0; i<m; i++) {
913 for (OrdinalType j=0; j<n; j++) {
914 rcmax = TEUCHOS_MAX( rcmax, C[j] );
915 rcmin = TEUCHOS_MIN( rcmin, C[j] );
918 if (rcmin == mZero) {
920 for (OrdinalType j=0; j<n; j++) {
926 for (OrdinalType j=0; j<n; j++) {
927 C[j] = mOne / TEUCHOS_MIN( TEUCHOS_MAX( C[j], smlnum ), bignum );
930 *colcond = TEUCHOS_MAX( rcmin, smlnum ) / TEUCHOS_MIN( rcmax, bignum );
934 template<
typename OrdinalType,
typename ScalarType>
935 void LAPACK<OrdinalType,ScalarType>::GERFS(
const char&
TRANS,
const OrdinalType& n,
const OrdinalType& nrhs,
const ScalarType* A,
const OrdinalType& lda,
const ScalarType* AF,
const OrdinalType& ldaf,
const OrdinalType* IPIV,
const ScalarType* B,
const OrdinalType& ldb, ScalarType* X,
const OrdinalType& ldx, ScalarType* FERR, ScalarType* BERR, ScalarType* WORK, OrdinalType* IWORK, OrdinalType* info)
const
937 UndefinedLAPACKRoutine<ScalarType>::notDefined();
940 template<
typename OrdinalType,
typename ScalarType>
941 void LAPACK<OrdinalType,ScalarType>::GBEQU(
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& kl,
const OrdinalType& ku,
const ScalarType* A,
const OrdinalType& lda, MagnitudeType* R, MagnitudeType* C, MagnitudeType* rowcond, MagnitudeType* colcond, MagnitudeType* amax, OrdinalType* info)
const
954 }
else if (lda < kl+ku+1) {
967 if (m == 0 || n == 0) {
979 for (OrdinalType i=0; i<m; i++) {
984 for (OrdinalType j=0; j<n; j++) {
985 for (OrdinalType i=TEUCHOS_MAX(j-ku,0); i<TEUCHOS_MIN(j+kl,m-1); i++) {
991 MagnitudeType rcmin = bignum;
992 MagnitudeType rcmax = mZero;
993 for (OrdinalType i=0; i<m; i++) {
994 rcmax = TEUCHOS_MAX( rcmax, R[i] );
995 rcmin = TEUCHOS_MIN( rcmin, R[i] );
999 if (rcmin == mZero) {
1001 for (OrdinalType i=0; i<m; i++) {
1007 for (OrdinalType i=0; i<m; i++) {
1008 R[i] = mOne / TEUCHOS_MIN( TEUCHOS_MAX( R[i], smlnum ), bignum );
1011 *rowcond = TEUCHOS_MAX( rcmin, smlnum ) / TEUCHOS_MIN( rcmax, bignum );
1015 for (OrdinalType j=0; j<n; j++) {
1020 for (OrdinalType j=0; j<n; j++) {
1021 for (OrdinalType i=TEUCHOS_MAX(j-ku,0); i<TEUCHOS_MIN(j+kl,m-1); i++) {
1029 for (OrdinalType j=0; j<n; j++) {
1030 rcmax = TEUCHOS_MAX( rcmax, C[j] );
1031 rcmin = TEUCHOS_MIN( rcmin, C[j] );
1034 if (rcmin == mZero) {
1036 for (OrdinalType j=0; j<n; j++) {
1042 for (OrdinalType j=0; j<n; j++) {
1043 C[j] = mOne / TEUCHOS_MIN( TEUCHOS_MAX( C[j], smlnum ), bignum );
1046 *colcond = TEUCHOS_MAX( rcmin, smlnum ) / TEUCHOS_MIN( rcmax, bignum );
1050 template<
typename OrdinalType,
typename ScalarType>
1051 void LAPACK<OrdinalType,ScalarType>::GBRFS(
const char&
TRANS,
const OrdinalType& n,
const OrdinalType& kl,
const OrdinalType& ku,
const OrdinalType& nrhs,
const ScalarType* A,
const OrdinalType& lda,
const ScalarType* AF,
const OrdinalType& ldaf,
const OrdinalType* IPIV,
const ScalarType* B,
const OrdinalType& ldb, ScalarType* X,
const OrdinalType& ldx, ScalarType* FERR, ScalarType* BERR, ScalarType* WORK, OrdinalType* IWORK, OrdinalType* info)
const
1053 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1056 template<
typename OrdinalType,
typename ScalarType>
1057 void LAPACK<OrdinalType,ScalarType>::GESVX(
const char& FACT,
const char&
TRANS,
const OrdinalType& n,
const OrdinalType& nrhs, ScalarType* A,
const OrdinalType& lda, ScalarType* AF,
const OrdinalType& ldaf, OrdinalType* IPIV,
char* EQUED, ScalarType* R, ScalarType* C, ScalarType* B,
const OrdinalType& ldb, ScalarType* X,
const OrdinalType& ldx, ScalarType* rcond, ScalarType* FERR, ScalarType* BERR, ScalarType* WORK, OrdinalType* IWORK, OrdinalType* info)
const
1059 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1062 template<
typename OrdinalType,
typename ScalarType>
1063 void LAPACK<OrdinalType,ScalarType>::SYTRD(
const char& UPLO,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* D, ScalarType* E, ScalarType* TAU, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const
1065 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1068 template<
typename OrdinalType,
typename ScalarType>
1069 void LAPACK<OrdinalType,ScalarType>::GEHRD(
const OrdinalType& n,
const OrdinalType& ilo,
const OrdinalType& ihi, ScalarType* A,
const OrdinalType& lda, ScalarType* TAU, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const
1071 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1074 template<
typename OrdinalType,
typename ScalarType>
1075 void LAPACK<OrdinalType,ScalarType>::TRTRS(
const char& UPLO,
const char&
TRANS,
const char& DIAG,
const OrdinalType& n,
const OrdinalType& nrhs,
const ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, OrdinalType* info)
const
1077 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1080 template<
typename OrdinalType,
typename ScalarType>
1083 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1086 template<
typename OrdinalType,
typename ScalarType>
1087 void LAPACK<OrdinalType,ScalarType>::SPEV(
const char& JOBZ,
const char& UPLO,
const OrdinalType& n, ScalarType* AP, ScalarType* W, ScalarType* Z,
const OrdinalType& ldz, ScalarType* WORK, OrdinalType* info)
const
1089 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1092 template<
typename OrdinalType,
typename ScalarType>
1093 void LAPACK<OrdinalType,ScalarType>::SYEV(
const char& JOBZ,
const char& UPLO,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* W, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const
1095 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1098 template<
typename OrdinalType,
typename ScalarType>
1099 void LAPACK<OrdinalType,ScalarType>::SYGV(
const OrdinalType& itype,
const char& JOBZ,
const char& UPLO,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, ScalarType* W, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const
1101 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1104 template<
typename OrdinalType,
typename ScalarType>
1105 void LAPACK<OrdinalType,ScalarType>::HEEV(
const char& JOBZ,
const char& UPLO,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, MagnitudeType* W, ScalarType* WORK,
const OrdinalType& lwork, MagnitudeType* RWORK, OrdinalType* info)
const
1107 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1110 template<
typename OrdinalType,
typename ScalarType>
1111 void LAPACK<OrdinalType,ScalarType>::HEGV(
const OrdinalType& itype,
const char& JOBZ,
const char& UPLO,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, MagnitudeType* W, ScalarType* WORK,
const OrdinalType& lwork, MagnitudeType* RWORK, OrdinalType* info)
const
1113 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1116 template<
typename OrdinalType,
typename ScalarType>
1117 void LAPACK<OrdinalType,ScalarType>::STEQR(
const char& COMPZ,
const OrdinalType& n, MagnitudeType* D, MagnitudeType* E, ScalarType* Z,
const OrdinalType& ldz, MagnitudeType* WORK, OrdinalType* info)
const
1119 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1122 template<
typename OrdinalType,
typename ScalarType>
1123 void LAPACK<OrdinalType,ScalarType>::PTEQR(
const char& COMPZ,
const OrdinalType& n, MagnitudeType* D, MagnitudeType* E, ScalarType* Z,
const OrdinalType& ldz, MagnitudeType* WORK, OrdinalType* info)
const
1125 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1128 template<
typename OrdinalType,
typename ScalarType>
1129 void LAPACK<OrdinalType, ScalarType>::HSEQR(
const char& JOB,
const char& COMPZ,
const OrdinalType& n,
const OrdinalType& ilo,
const OrdinalType& ihi, ScalarType* H,
const OrdinalType& ldh, ScalarType* WR, ScalarType* WI, ScalarType* Z,
const OrdinalType& ldz, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const
1131 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1134 template<
typename OrdinalType,
typename ScalarType>
1135 void LAPACK<OrdinalType, ScalarType>::GEES(
const char& JOBVS,
const char& SORT, OrdinalType& (*ptr2func)(ScalarType*, ScalarType*),
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, OrdinalType* sdim, ScalarType* WR, ScalarType* WI, ScalarType* VS,
const OrdinalType& ldvs, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* BWORK, OrdinalType* info)
const
1137 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1140 template<
typename OrdinalType,
typename ScalarType>
1141 void LAPACK<OrdinalType, ScalarType>::GEES(
const char& JOBVS,
const char& SORT, OrdinalType& (*ptr2func)(ScalarType*),
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, OrdinalType* sdim, ScalarType* W, ScalarType* VS,
const OrdinalType& ldvs, ScalarType* WORK,
const OrdinalType& lwork, MagnitudeType *RWORK, OrdinalType* BWORK, OrdinalType* info)
const
1143 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1146 template<
typename OrdinalType,
typename ScalarType>
1147 void LAPACK<OrdinalType, ScalarType>::GEES(
const char& JOBVS,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, OrdinalType* sdim, MagnitudeType* WR, MagnitudeType* WI, ScalarType* VS,
const OrdinalType& ldvs, ScalarType* WORK,
const OrdinalType& lwork, MagnitudeType *RWORK, OrdinalType* BWORK, OrdinalType* info)
const
1149 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1152 template<
typename OrdinalType,
typename ScalarType>
1153 void LAPACK<OrdinalType, ScalarType>::GEEV(
const char& JOBVL,
const char& JOBVR,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, MagnitudeType* WR, MagnitudeType* WI, ScalarType* VL,
const OrdinalType& ldvl, ScalarType* VR,
const OrdinalType& ldvr, ScalarType* WORK,
const OrdinalType& lwork, MagnitudeType* rwork, OrdinalType* info)
const
1155 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1158 template<
typename OrdinalType,
typename ScalarType>
1159 void LAPACK<OrdinalType, ScalarType>::GEEVX(
const char& BALANC,
const char& JOBVL,
const char& JOBVR,
const char& SENSE,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* WR, ScalarType* WI, ScalarType* VL,
const OrdinalType& ldvl, ScalarType* VR,
const OrdinalType& ldvr, OrdinalType* ilo, OrdinalType* ihi, MagnitudeType* SCALE, MagnitudeType* abnrm, MagnitudeType* RCONDE, MagnitudeType* RCONDV, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* IWORK, OrdinalType* info)
const
1161 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1164 template<
typename OrdinalType,
typename ScalarType>
1165 void LAPACK<OrdinalType, ScalarType>::GESVD(
const char& JOBU,
const char& JOBVT,
const OrdinalType& m,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, MagnitudeType* S, ScalarType* U,
const OrdinalType& ldu, ScalarType* V,
const OrdinalType& ldv, ScalarType* WORK,
const OrdinalType& lwork, MagnitudeType* RWORK, OrdinalType* info)
const
1167 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1170 template<
typename OrdinalType,
typename ScalarType>
1171 void LAPACK<OrdinalType, ScalarType>::GGEVX(
const char& BALANC,
const char& JOBVL,
const char& JOBVR,
const char& SENSE,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, MagnitudeType* ALPHAR, MagnitudeType* ALPHAI, ScalarType* BETA, ScalarType* VL,
const OrdinalType& ldvl, ScalarType* VR,
const OrdinalType& ldvr, OrdinalType* ilo, OrdinalType* ihi, MagnitudeType* lscale, MagnitudeType* rscale, MagnitudeType* abnrm, MagnitudeType* bbnrm, MagnitudeType* RCONDE, MagnitudeType* RCONDV, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* IWORK, OrdinalType* BWORK, OrdinalType* info)
const
1173 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1176 template<
typename OrdinalType,
typename ScalarType>
1177 void LAPACK<OrdinalType, ScalarType>::GGEV(
const char& JOBVL,
const char& JOBVR,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, MagnitudeType *ALPHAR, MagnitudeType *ALPHAI, ScalarType* BETA, ScalarType* VL,
const OrdinalType& ldvl, ScalarType* VR,
const OrdinalType& ldvr, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const
1179 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1183 template<
typename OrdinalType,
typename ScalarType>
1184 void LAPACK<OrdinalType,ScalarType>::TRSEN(
const char& JOB,
const char& COMPQ,
const OrdinalType* SELECT,
const OrdinalType& n, ScalarType* T,
const OrdinalType& ldt, ScalarType* Q,
const OrdinalType& ldq, MagnitudeType *WR, MagnitudeType *WI, OrdinalType* M, ScalarType* S, MagnitudeType *SEP, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* IWORK,
const OrdinalType& liwork, OrdinalType* info )
const
1186 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1190 template<
typename OrdinalType,
typename ScalarType>
1191 void LAPACK<OrdinalType,ScalarType>::TGSEN(
const OrdinalType& ijob,
const OrdinalType& wantq,
const OrdinalType& wantz,
const OrdinalType* SELECT,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, MagnitudeType *ALPHAR, MagnitudeType *ALPHAI, MagnitudeType *BETA, ScalarType* Q,
const OrdinalType& ldq, ScalarType* Z,
const OrdinalType& ldz, OrdinalType* M, MagnitudeType *PL, MagnitudeType *PR, MagnitudeType *DIF, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* IWORK,
const OrdinalType& liwork, OrdinalType* info )
const
1193 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1197 template<
typename OrdinalType,
typename ScalarType>
1198 void LAPACK<OrdinalType, ScalarType>::GGES(
const char& JOBVL,
const char& JOBVR,
const char& SORT, OrdinalType& (*ptr2func)(ScalarType*, ScalarType*, ScalarType*),
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, ScalarType* B,
const OrdinalType& ldb, OrdinalType* sdim, MagnitudeType* ALPHAR, MagnitudeType* ALPHAI, MagnitudeType* BETA, ScalarType* VL,
const OrdinalType& ldvl, ScalarType* VR,
const OrdinalType& ldvr, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* BWORK, OrdinalType* info )
const
1200 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1203 template<
typename OrdinalType,
typename ScalarType>
1204 void LAPACK<OrdinalType, ScalarType>::ORMQR(
const char& SIDE,
const char&
TRANS,
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& k,
const ScalarType* A,
const OrdinalType& lda,
const ScalarType* TAU, ScalarType* C,
const OrdinalType& ldc, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const
1206 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1209 template<
typename OrdinalType,
typename ScalarType>
1210 void LAPACK<OrdinalType, ScalarType>::ORM2R(
const char& SIDE,
const char&
TRANS,
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& k,
const ScalarType* A,
const OrdinalType& lda,
const ScalarType* TAU, ScalarType* C,
const OrdinalType& ldc, ScalarType* WORK, OrdinalType*
const info)
const
1212 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1215 template<
typename OrdinalType,
typename ScalarType>
1216 void LAPACK<OrdinalType, ScalarType>::UNMQR(
const char& SIDE,
const char&
TRANS,
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& k,
const ScalarType* A,
const OrdinalType& lda,
const ScalarType* TAU, ScalarType* C,
const OrdinalType& ldc, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const
1218 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1221 template<
typename OrdinalType,
typename ScalarType>
1222 void LAPACK<OrdinalType, ScalarType>::UNM2R(
const char& SIDE,
const char&
TRANS,
const OrdinalType& M,
const OrdinalType& N,
const OrdinalType& K,
const ScalarType* A,
const OrdinalType& LDA,
const ScalarType* TAU, ScalarType* C,
const OrdinalType& LDC, ScalarType* WORK, OrdinalType*
const INFO)
const
1224 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1227 template<
typename OrdinalType,
typename ScalarType>
1228 void LAPACK<OrdinalType, ScalarType>::ORGQR(
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& k, ScalarType* A,
const OrdinalType& lda,
const ScalarType* TAU, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const
1230 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1233 template<
typename OrdinalType,
typename ScalarType>
1234 void LAPACK<OrdinalType, ScalarType>::UNGQR(
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& k, ScalarType* A,
const OrdinalType& lda,
const ScalarType* TAU, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const
1236 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1239 template<
typename OrdinalType,
typename ScalarType>
1240 void LAPACK<OrdinalType, ScalarType>::ORGHR(
const OrdinalType& n,
const OrdinalType& ilo,
const OrdinalType& ihi, ScalarType* A,
const OrdinalType& lda,
const ScalarType* TAU, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const
1242 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1245 template<
typename OrdinalType,
typename ScalarType>
1246 void LAPACK<OrdinalType, ScalarType>::ORMHR(
const char& SIDE,
const char&
TRANS,
const OrdinalType& m,
const OrdinalType& n,
const OrdinalType& ilo,
const OrdinalType& ihi,
const ScalarType* A,
const OrdinalType& lda,
const ScalarType* TAU, ScalarType* C,
const OrdinalType& ldc, ScalarType* WORK,
const OrdinalType& lwork, OrdinalType* info)
const
1248 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1251 template<
typename OrdinalType,
typename ScalarType>
1252 void LAPACK<OrdinalType, ScalarType>::TREVC(
const char& SIDE,
const char& HOWMNY, OrdinalType* select,
const OrdinalType& n,
const ScalarType* T,
const OrdinalType& ldt, ScalarType* VL,
const OrdinalType& ldvl, ScalarType* VR,
const OrdinalType& ldvr,
const OrdinalType& mm, OrdinalType* m, ScalarType* WORK, OrdinalType* info)
const
1254 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1257 template<
typename OrdinalType,
typename ScalarType>
1258 void LAPACK<OrdinalType, ScalarType>::TREVC(
const char& SIDE,
const OrdinalType& n,
const ScalarType* T,
const OrdinalType& ldt, ScalarType* VL,
const OrdinalType& ldvl, ScalarType* VR,
const OrdinalType& ldvr,
const OrdinalType& mm, OrdinalType* m, ScalarType* WORK, MagnitudeType* RWORK, OrdinalType* info)
const
1260 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1263 template<
typename OrdinalType,
typename ScalarType>
1264 void LAPACK<OrdinalType, ScalarType>::TREXC(
const char& COMPQ,
const OrdinalType& n, ScalarType* T,
const OrdinalType& ldt, ScalarType* Q,
const OrdinalType& ldq, OrdinalType* ifst, OrdinalType* ilst, ScalarType* WORK, OrdinalType* info)
const
1266 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1270 template<
typename OrdinalType,
typename ScalarType>
1271 void LAPACK<OrdinalType, ScalarType>::TGEVC(
const char& SIDE,
const char& HOWMNY,
const OrdinalType* SELECT,
const OrdinalType& n,
const ScalarType* S,
const OrdinalType& lds,
const ScalarType* P,
const OrdinalType& ldp, ScalarType* VL,
const OrdinalType& ldvl, ScalarType* VR,
const OrdinalType& ldvr,
const OrdinalType& mm, OrdinalType* M, ScalarType* WORK, OrdinalType* info)
const
1273 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1277 template<
typename OrdinalType,
typename ScalarType>
1280 return UndefinedLAPACKRoutine<ScalarType>::notDefined();
1283 template<
typename OrdinalType,
typename ScalarType>
1284 OrdinalType
LAPACK<OrdinalType, ScalarType>::ILAENV(
const OrdinalType& ispec,
const std::string& NAME,
const std::string& OPTS,
const OrdinalType& N1,
const OrdinalType& N2,
const OrdinalType& N3,
const OrdinalType& N4 )
const
1286 return UndefinedLAPACKRoutine<OrdinalType>::notDefined();
1289 template<
typename OrdinalType,
typename ScalarType>
1292 return UndefinedLAPACKRoutine<ScalarType>::notDefined();
1295 template<
typename OrdinalType,
typename ScalarType>
1298 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1301 template<
typename OrdinalType,
typename ScalarType>
1304 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1307 template<
typename OrdinalType,
typename ScalarType>
1308 void LAPACK<OrdinalType, ScalarType>::GEBAL(
const char& JOBZ,
const OrdinalType& n, ScalarType* A,
const OrdinalType& lda, OrdinalType* ilo, OrdinalType* ihi, MagnitudeType* scale, OrdinalType* info )
const
1310 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1314 template<
typename OrdinalType,
typename ScalarType>
1315 void LAPACK<OrdinalType, ScalarType>::GEBAK(
const char& JOBZ,
const char& SIDE,
const OrdinalType& n,
const OrdinalType& ilo,
const OrdinalType& ihi,
const MagnitudeType* scale,
const OrdinalType& m, ScalarType* V,
const OrdinalType& ldv, OrdinalType* info )
const
1317 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1320 template<
typename OrdinalType,
typename ScalarType>
1323 return UndefinedLAPACKRoutine<ScalarType>::notDefined();
1326 template<
typename OrdinalType,
typename ScalarType>
1329 UndefinedLAPACKRoutine<ScalarType>::notDefined();
1334 #ifndef DOXYGEN_SHOULD_SKIP_THIS
1339 class TEUCHOSNUMERICS_LIB_DLL_EXPORT LAPACK<int, float>
1342 inline LAPACK(
void) {}
1343 inline LAPACK(
const LAPACK<int, float>& ) {}
1344 inline virtual ~LAPACK(
void) {}
1347 void PTTRF(
const int& n,
float* d,
float* e,
int* info)
const;
1348 void PTTRS(
const int& n,
const int& nrhs,
const float* d,
const float* e,
float* B,
const int& ldb,
int* info)
const;
1349 void POTRF(
const char& UPLO,
const int& n,
float* A,
const int& lda,
int* info)
const;
1350 void POTRS(
const char& UPLO,
const int& n,
const int& nrhs,
const float* A,
const int& lda,
float* B,
const int& ldb,
int* info)
const;
1351 void POTRI(
const char& UPLO,
const int& n,
float* A,
const int& lda,
int* info)
const;
1352 void POCON(
const char& UPLO,
const int& n,
const float* A,
const int& lda,
const float& anorm,
float* rcond,
float* WORK,
int* IWORK,
int* info)
const;
1353 void POSV(
const char& UPLO,
const int& n,
const int& nrhs,
float* A,
const int& lda,
float* B,
const int& ldb,
int* info)
const;
1354 void POEQU(
const int& n,
const float* A,
const int& lda,
float* S,
float* scond,
float* amax,
int* info)
const;
1355 void PORFS(
const char& UPLO,
const int& n,
const int& nrhs,
const float* A,
const int& lda,
const float* AF,
const int& ldaf,
const float* B,
const int& ldb,
float* X,
const int& ldx,
float* FERR,
float* BERR,
float* WORK,
int* IWORK,
int* info)
const;
1357 void POSVX(
const char& FACT,
const char& UPLO,
const int& n,
const int& nrhs,
float* A,
const int& lda,
float* AF,
const int& ldaf,
char* EQUED,
float* S,
float* B,
const int& ldb,
float* X,
const int& ldx,
float* rcond,
float* FERR,
float* BERR,
float* WORK,
int* IWORK,
int* info)
const;
1360 void GELS(
const char&
TRANS,
const int& m,
const int& n,
const int& nrhs,
float* A,
const int& lda,
float* B,
const int& ldb,
float* WORK,
const int& lwork,
int* info)
const;
1361 void GELSS(
const int& m,
const int& n,
const int& nrhs,
float* A,
const int& lda,
float* B,
const int& ldb,
float* S,
const float& rcond,
int* rank,
float* WORK,
const int& lwork,
float* RWORK,
int* info)
const;
1362 void GELSS(
const int& m,
const int& n,
const int& nrhs,
float* A,
const int& lda,
float* B,
const int& ldb,
float* S,
const float& rcond,
int* rank,
float* WORK,
const int& lwork,
int* info)
const;
1363 void GGLSE(
const int& m,
const int& n,
const int& p,
float* A,
const int& lda,
float* B,
const int& ldb,
float* C,
float* D,
float* X,
float* WORK,
const int& lwork,
int* info)
const;
1364 void GEQRF(
const int& m,
const int& n,
float* A,
const int& lda,
float* TAU,
float* WORK,
const int& lwork,
int* info)
const;
1365 void GEQR2(
const int& m,
const int& n,
float* A,
const int& lda,
float* TAU,
float* WORK,
int*
const info)
const;
1367 void GETRF(
const int& m,
const int& n,
float* A,
const int& lda,
int* IPIV,
int* info)
const;
1368 void GETRS(
const char&
TRANS,
const int& n,
const int& nrhs,
const float* A,
const int& lda,
const int* IPIV,
float* B,
const int& ldb,
int* info)
const;
1369 void LASCL(
const char& TYPE,
const int& kl,
const int& ku,
const float& cfrom,
const float& cto,
const int& m,
const int& n,
float* A,
const int& lda,
int* info)
const;
1371 void GEQP3 (
const int& m,
const int& n,
float* A,
const int& lda,
int* jpvt,
float* TAU,
float* WORK,
const int& lwork,
float* RWORK,
int* info)
const;
1372 void LASWP (
const int& N,
float* A,
const int& LDA,
const int& K1,
const int& K2,
const int* IPIV,
const int& INCX)
const;
1374 void GBTRF(
const int& m,
const int& n,
const int& kl,
const int& ku,
float* A,
const int& lda,
int* IPIV,
int* info)
const;
1375 void GBTRS(
const char&
TRANS,
const int& n,
const int& kl,
const int& ku,
const int& nrhs,
const float* A,
const int& lda,
const int* IPIV,
float* B,
const int& ldb,
int* info)
const;
1376 void GTTRF(
const int& n,
float* dl,
float* d,
float* du,
float* du2,
int* IPIV,
int* info)
const;
1377 void GTTRS(
const char&
TRANS,
const int& n,
const int& nrhs,
const float* dl,
const float* d,
const float* du,
const float* du2,
const int* IPIV,
float* B,
const int& ldb,
int* info)
const;
1380 void GETRI(
const int& n,
float* A,
const int& lda,
const int* IPIV,
float* WORK,
const int& lwork,
int* info)
const;
1381 void LATRS (
const char& UPLO,
const char&
TRANS,
const char& DIAG,
const char& NORMIN,
const int& N,
const float* A,
const int& LDA,
float* X,
float* SCALE,
float* CNORM,
int* INFO)
const;
1382 void GECON(
const char& NORM,
const int& n,
const float* A,
const int& lda,
const float& anorm,
float* rcond,
float* WORK,
int* IWORK,
int* info)
const;
1383 void GBCON(
const char& NORM,
const int& n,
const int& kl,
const int& ku,
const float* A,
const int& lda,
const int* IPIV,
const float& anorm,
float* rcond,
float* WORK,
int* IWORK,
int* info)
const;
1384 float LANGB(
const char& NORM,
const int& n,
const int& kl,
const int& ku,
const float* A,
const int& lda,
float* WORK)
const;
1385 void GESV(
const int& n,
const int& nrhs,
float* A,
const int& lda,
int* IPIV,
float* B,
const int& ldb,
int* info)
const;
1386 void GEEQU(
const int& m,
const int& n,
const float* A,
const int& lda,
float* R,
float* C,
float* rowcond,
float* colcond,
float* amax,
int* info)
const;
1387 void GERFS(
const char&
TRANS,
const int& n,
const int& nrhs,
const float* A,
const int& lda,
const float* AF,
const int& ldaf,
const int* IPIV,
const float* B,
const int& ldb,
float* X,
const int& ldx,
float* FERR,
float* BERR,
float* WORK,
int* IWORK,
int* info)
const;
1388 void GBEQU(
const int& m,
const int& n,
const int& kl,
const int& ku,
const float* A,
const int& lda,
float* R,
float* C,
float* rowcond,
float* colcond,
float* amax,
int* info)
const;
1389 void GBRFS(
const char&
TRANS,
const int& n,
const int& kl,
const int& ku,
const int& nrhs,
const float* A,
const int& lda,
const float* AF,
const int& ldaf,
const int* IPIV,
const float* B,
const int& ldb,
float* X,
const int& ldx,
float* FERR,
float* BERR,
float* WORK,
int* IWORK,
int* info)
const;
1391 void GESVX(
const char& FACT,
const char&
TRANS,
const int& n,
const int& nrhs,
float* A,
const int& lda,
float* AF,
const int& ldaf,
int* IPIV,
char* EQUED,
float* R,
float* C,
float* B,
const int& ldb,
float* X,
const int& ldx,
float* rcond,
float* FERR,
float* BERR,
float* WORK,
int* IWORK,
int* info)
const;
1393 void SYTRD(
const char& UPLO,
const int& n,
float* A,
const int& lda,
float* D,
float* E,
float* TAU,
float* WORK,
const int& lwork,
int* info)
const;
1394 void GEHRD(
const int& n,
const int& ilo,
const int& ihi,
float* A,
const int& lda,
float* TAU,
float* WORK,
const int& lwork,
int* info)
const;
1395 void TRTRS(
const char& UPLO,
const char&
TRANS,
const char& DIAG,
const int& n,
const int& nrhs,
const float* A,
const int& lda,
float* B,
const int& ldb,
int* info)
const;
1396 void TRTRI(
const char& UPLO,
const char& DIAG,
const int& n,
float* A,
const int& lda,
int* info)
const;
1399 void STEQR(
const char& COMPZ,
const int& n,
float* D,
float* E,
float* Z,
const int& ldz,
float* WORK,
int* info)
const;
1400 void PTEQR(
const char& COMPZ,
const int& n,
float* D,
float* E,
float* Z,
const int& ldz,
float* WORK,
int* info)
const;
1401 void SPEV(
const char& JOBZ,
const char& UPLO,
const int& n,
float* AP,
float* W,
float* Z,
const int& ldz,
float* WORK,
int* info)
const;
1402 void SYEV(
const char& JOBZ,
const char& UPLO,
const int& n,
float* A,
const int& lda,
float* W,
float* WORK,
const int& lwork,
int* info)
const;
1403 void SYGV(
const int& itype,
const char& JOBZ,
const char& UPLO,
const int& n,
float* A,
const int& lda,
float* B,
const int& ldb,
float* W,
float* WORK,
const int& lwork,
int* info)
const;
1404 void HEEV(
const char& JOBZ,
const char& UPLO,
const int& n,
float* A,
const int& lda,
float* W,
float* WORK,
const int& lwork,
float* RWORK,
int* info)
const;
1405 void HEGV(
const int& itype,
const char& JOBZ,
const char& UPLO,
const int& n,
float* A,
const int& lda,
float* B,
const int& ldb,
float* W,
float* WORK,
const int& lwork,
float* RWORK,
int* info)
const;
1408 void HSEQR(
const char& JOB,
const char& COMPZ,
const int& n,
const int& ilo,
const int& ihi,
float* H,
const int& ldh,
float* WR,
float* WI,
float* Z,
const int& ldz,
float* WORK,
const int& lwork,
int* info)
const;
1409 void GEES(
const char& JOBVS,
const char& SORT,
int (*ptr2func)(
float*,
float*),
const int& n,
float* A,
const int& lda,
int* sdim,
float* WR,
float* WI,
float* VS,
const int& ldvs,
float* WORK,
const int& lwork,
int* BWORK,
int* info)
const;
1410 void GEES(
const char& JOBVS,
const int& n,
float* A,
const int& lda,
int* sdim,
float* WR,
float* WI,
float* VS,
const int& ldvs,
float* WORK,
const int& lwork,
float* RWORK,
int* BWORK,
int* info)
const;
1412 void GEEV(
const char& JOBVL,
const char& JOBVR,
const int& n,
float* A,
const int& lda,
float* WR,
float* WI,
float* VL,
const int& ldvl,
float* VR,
const int& ldvr,
float* WORK,
const int& lwork,
int* info)
const;
1413 void GEEV(
const char& JOBVL,
const char& JOBVR,
const int& n,
float* A,
const int& lda,
float* WR,
float* WI,
float* VL,
const int& ldvl,
float* VR,
const int& ldvr,
float* WORK,
const int& lwork,
float* rwork,
int* info)
const;
1415 void GEEVX(
const char& BALANC,
const char& JOBVL,
const char& JOBVR,
const char& SENSE,
const int& n,
float* A,
const int& lda,
float* WR,
float* WI,
float* VL,
const int& ldvl,
float* VR,
const int& ldvr,
int* ilo,
int* ihi,
float* SCALE,
float* abnrm,
float* RCONDE,
float* RCONDV,
float* WORK,
const int& lwork,
int* IWORK,
int* info)
const;
1416 void GGEVX(
const char& BALANC,
const char& JOBVL,
const char& JOBVR,
const char& SENSE,
const int& n,
float* A,
const int& lda,
float* B,
const int& ldb,
float* ALPHAR,
float* ALPHAI,
float* BETA,
float* VL,
const int& ldvl,
float* VR,
const int& ldvr,
int* ilo,
int* ihi,
float* lscale,
float* rscale,
float* abnrm,
float* bbnrm,
float* RCONDE,
float* RCONDV,
float* WORK,
const int& lwork,
int* IWORK,
int* BWORK,
int* info)
const;
1417 void GGEVX(
const char& BALANC,
const char& JOBVL,
const char& JOBVR,
const char& SENSE,
const int& n,
float* A,
const int& lda,
float* B,
const int& ldb,
float* ALPHAR,
float* ALPHAI,
float* BETA,
float* VL,
const int& ldvl,
float* VR,
const int& ldvr,
int* ilo,
int* ihi,
float* lscale,
float* rscale,
float* abnrm,
float* bbnrm,
float* RCONDE,
float* RCONDV,
float* WORK,
const int& lwork,
float* rwork,
int* IWORK,
int* BWORK,
int* info)
const;
1418 void GGEV(
const char& JOBVL,
const char& JOBVR,
const int& n,
float* A,
const int& lda,
float* B,
const int& ldb,
float* ALPHAR,
float* ALPHAI,
float* BETA,
float* VL,
const int& ldvl,
float* VR,
const int& ldvr,
float* WORK,
const int& lwork,
int* info)
const;
1419 void TRSEN(
const char& JOB,
const char& COMPQ,
const int* SELECT,
const int& n,
float* T,
const int& ldt,
float* Q,
const int& ldq,
float* WR,
float* WI,
int* M,
float* S,
float* SEP,
float* WORK,
const int& lwork,
int* IWORK,
const int& liwork,
int* info )
const;
1420 void TGSEN(
const int& ijob,
const int& wantq,
const int& wantz,
const int* SELECT,
const int& n,
float* A,
const int& lda,
float* B,
const int& ldb,
float* ALPHAR,
float* ALPHAI,
float* BETA,
float* Q,
const int& ldq,
float* Z,
const int& ldz,
int* M,
float* PL,
float* PR,
float* DIF,
float* WORK,
const int& lwork,
int* IWORK,
const int& liwork,
int* info )
const;
1421 void GGES(
const char& JOBVL,
const char& JOBVR,
const char& SORT,
int (*ptr2func)(
float*,
float*,
float*),
const int& n,
float* A,
const int& lda,
float* B,
const int& ldb,
int* sdim,
float* ALPHAR,
float* ALPHAI,
float* BETA,
float* VL,
const int& ldvl,
float* VR,
const int& ldvr,
float* WORK,
const int& lwork,
int* bwork,
int* info )
const;
1424 void GESVD(
const char& JOBU,
const char& JOBVT,
const int& m,
const int& n,
float* A,
const int& lda,
float* S,
float* U,
const int& ldu,
float* V,
const int& ldv,
float* WORK,
const int& lwork,
float* RWORK,
int* info)
const;
1427 void ORMQR(
const char& SIDE,
const char&
TRANS,
const int& m,
const int& n,
const int& k,
const float* A,
const int& lda,
const float* TAU,
float* C,
const int& ldc,
float* WORK,
const int& lwork,
int* info)
const;
1428 void ORM2R(
const char& SIDE,
const char&
TRANS,
const int& m,
const int& n,
const int& k,
const float* A,
const int& lda,
const float* TAU,
float* C,
const int& ldc,
float* WORK,
int*
const info)
const;
1429 void UNMQR(
const char& SIDE,
const char&
TRANS,
const int& m,
const int& n,
const int& k,
const float* A,
const int& lda,
const float* TAU,
float* C,
const int& ldc,
float* WORK,
const int& lwork,
int* info)
const;
1430 void UNM2R(
const char& SIDE,
const char&
TRANS,
const int& M,
const int& N,
const int& K,
const float* A,
const int& LDA,
const float* TAU,
float* C,
const int& LDC,
float* WORK,
int*
const INFO)
const;
1431 void ORGQR(
const int& m,
const int& n,
const int& k,
float* A,
const int& lda,
const float* TAU,
float* WORK,
const int& lwork,
int* info)
const;
1432 void UNGQR(
const int& m,
const int& n,
const int& k,
float* A,
const int& lda,
const float* TAU,
float* WORK,
const int& lwork,
int* info)
const;
1433 void ORGHR(
const int& n,
const int& ilo,
const int& ihi,
float* A,
const int& lda,
const float* TAU,
float* WORK,
const int& lwork,
int* info)
const;
1434 void ORMHR(
const char& SIDE,
const char&
TRANS,
const int& m,
const int& n,
const int& ilo,
const int& ihi,
const float* A,
const int& lda,
const float* TAU,
float* C,
const int& ldc,
float* WORK,
const int& lwork,
int* info)
const;
1437 void TREVC(
const char& SIDE,
const char& HOWMNY,
int* select,
const int& n,
const float* T,
const int& ldt,
float* VL,
const int& ldvl,
float* VR,
const int& ldvr,
const int& mm,
int* m,
float* WORK,
int* info)
const;
1438 void TREVC(
const char& SIDE,
const int& n,
const float* T,
const int& ldt,
float* VL,
const int& ldvl,
float* VR,
const int& ldvr,
const int& mm,
int* m,
float* WORK,
float* RWORK,
int* info)
const;
1440 void TREXC(
const char& COMPQ,
const int& n,
float* T,
const int& ldt,
float* Q,
const int& ldq,
int* ifst,
int* ilst,
float* WORK,
int* info)
const;
1442 void TGEVC(
const char& SIDE,
const char& HOWMNY,
const int* SELECT,
const int& n,
const float* S,
const int& lds,
const float* P,
const int& ldp,
float* VL,
const int& ldvl,
float* VR,
const int& ldvr,
const int& mm,
int* M,
float* WORK,
int* info)
const;
1445 void LARTG(
const float& f,
const float& g,
float* c,
float* s,
float* r )
const;
1446 void LARFG(
const int& n,
float* alpha,
float* x,
const int& incx,
float* tau )
const;
1450 void GEBAL(
const char& JOBZ,
const int& n,
float* A,
const int& lda,
int* ilo,
int* ihi,
float* scale,
int* info)
const;
1452 void GEBAK(
const char& JOBZ,
const char& SIDE,
const int& n,
const int& ilo,
const int& ihi,
const float* scale,
const int& m,
float* V,
const int& ldv,
int* info)
const;
1455 float LARND(
const int& idist,
int* seed )
const;
1456 void LARNV(
const int& idist,
int* seed,
const int& n,
float* v )
const;
1459 float LAMCH(
const char& CMACH)
const;
1460 int ILAENV(
const int& ispec,
const std::string& NAME,
const std::string& OPTS,
const int& N1 = -1,
const int& N2 = -1,
const int& N3 = -1,
const int& N4 = -1 )
const;
1463 float LAPY2(
const float& x,
const float& y)
const;
1472 class TEUCHOSNUMERICS_LIB_DLL_EXPORT LAPACK<int, double>
1475 inline LAPACK(
void) {}
1476 inline LAPACK(
const LAPACK<int, double>& ) {}
1477 inline virtual ~LAPACK(
void) {}
1480 void PTTRF(
const int& n,
double* d,
double* e,
int* info)
const;
1481 void PTTRS(
const int& n,
const int& nrhs,
const double* d,
const double* e,
double* B,
const int& ldb,
int* info)
const;
1482 void POTRF(
const char& UPLO,
const int& n,
double* A,
const int& lda,
int* info)
const;
1483 void POTRS(
const char& UPLO,
const int& n,
const int& nrhs,
const double* A,
const int& lda,
double* B,
const int& ldb,
int* info)
const;
1484 void POTRI(
const char& UPLO,
const int& n,
double* A,
const int& lda,
int* info)
const;
1485 void POCON(
const char& UPLO,
const int& n,
const double* A,
const int& lda,
const double& anorm,
double* rcond,
double* WORK,
int* IWORK,
int* info)
const;
1486 void POSV(
const char& UPLO,
const int& n,
const int& nrhs,
double* A,
const int& lda,
double* B,
const int& ldb,
int* info)
const;
1487 void POEQU(
const int& n,
const double* A,
const int& lda,
double* S,
double* scond,
double* amax,
int* info)
const;
1488 void PORFS(
const char& UPLO,
const int& n,
const int& nrhs,
const double* A,
const int& lda,
const double* AF,
const int& ldaf,
const double* B,
const int& ldb,
double* X,
const int& ldx,
double* FERR,
double* BERR,
double* WORK,
int* IWORK,
int* info)
const;
1490 void POSVX(
const char& FACT,
const char& UPLO,
const int& n,
const int& nrhs,
double* A,
const int& lda,
double* AF,
const int& ldaf,
char* EQUED,
double* S,
double* B,
const int& ldb,
double* X,
const int& ldx,
double* rcond,
double* FERR,
double* BERR,
double* WORK,
int* IWORK,
int* info)
const;
1493 void GELS(
const char&
TRANS,
const int& m,
const int& n,
const int& nrhs,
double* A,
const int& lda,
double* B,
const int& ldb,
double* WORK,
const int& lwork,
int* info)
const;
1494 void GELSS(
const int& m,
const int& n,
const int& nrhs,
double* A,
const int& lda,
double* B,
const int& ldb,
double* S,
const double& rcond,
int* rank,
double* WORK,
const int& lwork,
double* RWORK,
int* info)
const;
1495 void GELSS(
const int& m,
const int& n,
const int& nrhs,
double* A,
const int& lda,
double* B,
const int& ldb,
double* S,
const double& rcond,
int* rank,
double* WORK,
const int& lwork,
int* info)
const;
1496 void GGLSE(
const int& m,
const int& n,
const int& p,
double* A,
const int& lda,
double* B,
const int& ldb,
double* C,
double* D,
double* X,
double* WORK,
const int& lwork,
int* info)
const;
1497 void GEQRF(
const int& m,
const int& n,
double* A,
const int& lda,
double* TAU,
double* WORK,
const int& lwork,
int* info)
const;
1498 void GEQR2(
const int& m,
const int& n,
double* A,
const int& lda,
double* TAU,
double* WORK,
int*
const info)
const;
1499 void GETRF(
const int& m,
const int& n,
double* A,
const int& lda,
int* IPIV,
int* info)
const;
1500 void GETRS(
const char&
TRANS,
const int& n,
const int& nrhs,
const double* A,
const int& lda,
const int* IPIV,
double* B,
const int& ldb,
int* info)
const;
1501 void LASCL(
const char& TYPE,
const int& kl,
const int& ku,
const double& cfrom,
const double& cto,
const int& m,
const int& n,
double* A,
const int& lda,
int* info)
const;
1503 void GEQP3 (
const int& m,
const int& n,
double* A,
const int& lda,
int* jpvt,
double* TAU,
double* WORK,
const int& lwork,
double* RWORK,
int* info)
const;
1504 void LASWP (
const int& N,
double* A,
const int& LDA,
const int& K1,
const int& K2,
const int* IPIV,
const int& INCX)
const;
1506 void GBTRF(
const int& m,
const int& n,
const int& kl,
const int& ku,
double* A,
const int& lda,
int* IPIV,
int* info)
const;
1507 void GBTRS(
const char&
TRANS,
const int& n,
const int& kl,
const int& ku,
const int& nrhs,
const double* A,
const int& lda,
const int* IPIV,
double* B,
const int& ldb,
int* info)
const;
1508 void GTTRF(
const int& n,
double* dl,
double* d,
double* du,
double* du2,
int* IPIV,
int* info)
const;
1509 void GTTRS(
const char&
TRANS,
const int& n,
const int& nrhs,
const double* dl,
const double* d,
const double* du,
const double* du2,
const int* IPIV,
double* B,
const int& ldb,
int* info)
const;
1510 void GETRI(
const int& n,
double* A,
const int& lda,
const int* IPIV,
double* WORK,
const int& lwork,
int* info)
const;
1511 void LATRS (
const char& UPLO,
const char&
TRANS,
const char& DIAG,
const char& NORMIN,
const int& N,
const double* A,
const int& LDA,
double* X,
double* SCALE,
double* CNORM,
int* INFO)
const;
1512 void GECON(
const char& NORM,
const int& n,
const double* A,
const int& lda,
const double& anorm,
double* rcond,
double* WORK,
int* IWORK,
int* info)
const;
1513 void GBCON(
const char& NORM,
const int& n,
const int& kl,
const int& ku,
const double* A,
const int& lda,
const int* IPIV,
const double& anorm,
double* rcond,
double* WORK,
int* IWORK,
int* info)
const;
1514 double LANGB(
const char& NORM,
const int& n,
const int& kl,
const int& ku,
const double* A,
const int& lda,
double* WORK)
const;
1515 void GESV(
const int& n,
const int& nrhs,
double* A,
const int& lda,
int* IPIV,
double* B,
const int& ldb,
int* info)
const;
1516 void GEEQU(
const int& m,
const int& n,
const double* A,
const int& lda,
double* R,
double* C,
double* rowcond,
double* colcond,
double* amax,
int* info)
const;
1517 void GERFS(
const char&
TRANS,
const int& n,
const int& nrhs,
const double* A,
const int& lda,
const double* AF,
const int& ldaf,
const int* IPIV,
const double* B,
const int& ldb,
double* X,
const int& ldx,
double* FERR,
double* BERR,
double* WORK,
int* IWORK,
int* info)
const;
1518 void GBEQU(
const int& m,
const int& n,
const int& kl,
const int& ku,
const double* A,
const int& lda,
double* R,
double* C,
double* rowcond,
double* colcond,
double* amax,
int* info)
const;
1519 void GBRFS(
const char&
TRANS,
const int& n,
const int& kl,
const int& ku,
const int& nrhs,
const double* A,
const int& lda,
const double* AF,
const int& ldaf,
const int* IPIV,
const double* B,
const int& ldb,
double* X,
const int& ldx,
double* FERR,
double* BERR,
double* WORK,
int* IWORK,
int* info)
const;
1521 void GESVX(
const char& FACT,
const char&
TRANS,
const int& n,
const int& nrhs,
double* A,
const int& lda,
double* AF,
const int& ldaf,
int* IPIV,
char* EQUED,
double* R,
double* C,
double* B,
const int& ldb,
double* X,
const int& ldx,
double* rcond,
double* FERR,
double* BERR,
double* WORK,
int* IWORK,
int* info)
const;
1523 void SYTRD(
const char& UPLO,
const int& n,
double* A,
const int& lda,
double* D,
double* E,
double* TAU,
double* WORK,
const int& lwork,
int* info)
const;
1524 void GEHRD(
const int& n,
const int& ilo,
const int& ihi,
double* A,
const int& lda,
double* TAU,
double* WORK,
const int& lwork,
int* info)
const;
1525 void TRTRS(
const char& UPLO,
const char&
TRANS,
const char& DIAG,
const int& n,
const int& nrhs,
const double* A,
const int& lda,
double* B,
const int& ldb,
int* info)
const;
1526 void TRTRI(
const char& UPLO,
const char& DIAG,
const int& n,
double* A,
const int& lda,
int* info)
const;
1529 void STEQR(
const char& COMPZ,
const int& n,
double* D,
double* E,
double* Z,
const int& ldz,
double* WORK,
int* info)
const;
1530 void PTEQR(
const char& COMPZ,
const int& n,
double* D,
double* E,
double* Z,
const int& ldz,
double* WORK,
int* info)
const;
1531 void SPEV(
const char& JOBZ,
const char& UPLO,
const int& n,
double* AP,
double* W,
double* Z,
const int& ldz,
double* WORK,
int* info)
const;
1532 void SYEV(
const char& JOBZ,
const char& UPLO,
const int& n,
double* A,
const int& lda,
double* W,
double* WORK,
const int& lwork,
int* info)
const;
1533 void SYGV(
const int& itype,
const char& JOBZ,
const char& UPLO,
const int& n,
double* A,
const int& lda,
double* B,
const int& ldb,
double* W,
double* WORK,
const int& lwork,
int* info)
const;
1534 void HEEV(
const char& JOBZ,
const char& UPLO,
const int& n,
double* A,
const int& lda,
double* W,
double* WORK,
const int& lwork,
double* RWORK,
int* info)
const;
1535 void HEGV(
const int& itype,
const char& JOBZ,
const char& UPLO,
const int& n,
double* A,
const int& lda,
double* B,
const int& ldb,
double* W,
double* WORK,
const int& lwork,
double* RWORK,
int* info)
const;
1538 void HSEQR(
const char& JOB,
const char& COMPZ,
const int& n,
const int& ilo,
const int& ihi,
double* H,
const int& ldh,
double* WR,
double* WI,
double* Z,
const int& ldz,
double* WORK,
const int& lwork,
int* info)
const;
1539 void GEES(
const char& JOBVS,
const char& SORT,
int (*ptr2func)(
double*,
double*),
const int& n,
double* A,
const int& lda,
int* sdim,
double* WR,
double* WI,
double* VS,
const int& ldvs,
double* WORK,
const int& lwork,
int* BWORK,
int* info)
const;
1540 void GEES(
const char& JOBVS,
const int& n,
double* A,
const int& lda,
int* sdim,
double* WR,
double* WI,
double* VS,
const int& ldvs,
double* WORK,
const int& lwork,
double* RWORK,
int* BWORK,
int* info)
const;
1542 void GEEV(
const char& JOBVL,
const char& JOBVR,
const int& n,
double* A,
const int& lda,
double* WR,
double* WI,
double* VL,
const int& ldvl,
double* VR,
const int& ldvr,
double* WORK,
const int& lwork,
int* info)
const;
1543 void GEEV(
const char& JOBVL,
const char& JOBVR,
const int& n,
double* A,
const int& lda,
double* WR,
double* WI,
double* VL,
const int& ldvl,
double* VR,
const int& ldvr,
double* WORK,
const int& lwork,
double* RWORK,
int* info)
const;
1545 void GEEVX(
const char& BALANC,
const char& JOBVL,
const char& JOBVR,
const char& SENSE,
const int& n,
double* A,
const int& lda,
double* WR,
double* WI,
double* VL,
const int& ldvl,
double* VR,
const int& ldvr,
int* ilo,
int* ihi,
double* SCALE,
double* abnrm,
double* RCONDE,
double* RCONDV,
double* WORK,
const int& lwork,
int* IWORK,
int* info)
const;
1546 void GGEVX(
const char& BALANC,
const char& JOBVL,
const char& JOBVR,
const char& SENSE,
const int& n,
double* A,
const int& lda,
double* B,
const int& ldb,
double* ALPHAR,
double* ALPHAI,
double* BETA,
double* VL,
const int& ldvl,
double* VR,
const int& ldvr,
int* ilo,
int* ihi,
double* lscale,
double* rscale,
double* abnrm,
double* bbnrm,
double* RCONDE,
double* RCONDV,
double* WORK,
const int& lwork,
int* IWORK,
int* BWORK,
int* info)
const;
1547 void GGEVX(
const char& BALANC,
const char& JOBVL,
const char& JOBVR,
const char& SENSE,
const int& n,
double* A,
const int& lda,
double* B,
const int& ldb,
double* ALPHAR,
double* ALPHAI,
double* BETA,
double* VL,
const int& ldvl,
double* VR,
const int& ldvr,
int* ilo,
int* ihi,
double* lscale,
double* rscale,
double* abnrm,
double* bbnrm,
double* RCONDE,
double* RCONDV,
double* WORK,
const int& lwork,
double* rwork,
int* IWORK,
int* BWORK,
int* info)
const;
1548 void GGEV(
const char& JOBVL,
const char& JOBVR,
const int& n,
double* A,
const int& lda,
double* B,
const int& ldb,
double* ALPHAR,
double* ALPHAI,
double* BETA,
double* VL,
const int& ldvl,
double* VR,
const int& ldvr,
double* WORK,
const int& lwork,
int* info)
const;
1549 void TRSEN(
const char& JOB,
const char& COMPQ,
const int* SELECT,
const int& n,
double* T,
const int& ldt,
double* Q,
const int& ldq,
double* WR,
double* WI,
int* M,
double* S,
double* SEP,
double* WORK,
const int& lwork,
int* IWORK,
const int& liwork,
int* info )
const;
1550 void TGSEN(
const int& ijob,
const int& wantq,
const int& wantz,
const int* SELECT,
const int& n,
double* A,
const int& lda,
double* B,
const int& ldb,
double* ALPHAR,
double* ALPHAI,
double* BETA,
double* Q,
const int& ldq,
double* Z,
const int& ldz,
int* M,
double* PL,
double* PR,
double* DIF,
double* WORK,
const int& lwork,
int* IWORK,
const int& liwork,
int* info )
const;
1551 void GGES(
const char& JOBVL,
const char& JOBVR,
const char& SORT,
int (*ptr2func)(
double*,
double*,
double*),
const int& n,
double* A,
const int& lda,
double* B,
const int& ldb,
int* sdim,
double* ALPHAR,
double* ALPHAI,
double* BETA,
double* VL,
const int& ldvl,
double* VR,
const int& ldvr,
double* WORK,
const int& lwork,
int* bwork,
int* info )
const;
1555 void GESVD(
const char& JOBU,
const char& JOBVT,
const int& m,
const int& n,
double* A,
const int& lda,
double* S,
double* U,
const int& ldu,
double* V,
const int& ldv,
double* WORK,
const int& lwork,
double* RWORK,
int* info)
const;
1558 void ORMQR(
const char& SIDE,
const char&
TRANS,
const int& m,
const int& n,
const int& k,
const double* A,
const int& lda,
const double* TAU,
double* C,
const int& ldc,
double* WORK,
const int& lwork,
int* info)
const;
1559 void ORM2R(
const char& SIDE,
const char&
TRANS,
const int& m,
const int& n,
const int& k,
const double* A,
const int& lda,
const double* TAU,
double* C,
const int& ldc,
double* WORK,
int*
const info)
const;
1560 void UNMQR(
const char& SIDE,
const char&
TRANS,
const int& m,
const int& n,
const int& k,
const double* A,
const int& lda,
const double* TAU,
double* C,
const int& ldc,
double* WORK,
const int& lwork,
int* info)
const;
1561 void UNM2R(
const char& SIDE,
const char&
TRANS,
const int& M,
const int& N,
const int& K,
const double* A,
const int& LDA,
const double* TAU,
double* C,
const int& LDC,
double* WORK,
int*
const INFO)
const;
1562 void ORGQR(
const int& m,
const int& n,
const int& k,
double* A,
const int& lda,
const double* TAU,
double* WORK,
const int& lwork,
int* info)
const;
1563 void UNGQR(
const int& m,
const int& n,
const int& k,
double* A,
const int& lda,
const double* TAU,
double* WORK,
const int& lwork,
int* info)
const;
1564 void ORGHR(
const int& n,
const int& ilo,
const int& ihi,
double* A,
const int& lda,
const double* TAU,
double* WORK,
const int& lwork,
int* info)
const;
1565 void ORMHR(
const char& SIDE,
const char&
TRANS,
const int& m,
const int& n,
const int& ilo,
const int& ihi,
const double* A,
const int& lda,
const double* TAU,
double* C,
const int& ldc,
double* WORK,
const int& lwork,
int* info)
const;
1568 void TREVC(
const char& SIDE,
const char& HOWMNY,
int* select,
const int& n,
const double* T,
const int& ldt,
double* VL,
const int& ldvl,
double* VR,
const int& ldvr,
const int& mm,
int* m,
double* WORK,
int* info)
const;
1569 void TREVC(
const char& SIDE,
const int& n,
const double* T,
const int& ldt,
double* VL,
const int& ldvl,
double* VR,
const int& ldvr,
const int& mm,
int* m,
double* WORK,
double* RWORK,
int* info)
const;
1571 void TREXC(
const char& COMPQ,
const int& n,
double* T,
const int& ldt,
double* Q,
const int& ldq,
int* ifst,
int* ilst,
double* WORK,
int* info)
const;
1573 void TGEVC(
const char& SIDE,
const char& HOWMNY,
const int* SELECT,
const int& n,
const double* S,
const int& lds,
const double* P,
const int& ldp,
double* VL,
const int& ldvl,
double* VR,
const int& ldvr,
const int& mm,
int* M,
double* WORK,
int* info)
const;
1576 void LARTG(
const double& f,
const double& g,
double* c,
double* s,
double* r )
const;
1577 void LARFG(
const int& n,
double* alpha,
double* x,
const int& incx,
double* tau )
const;
1581 void GEBAL(
const char& JOBZ,
const int& n,
double* A,
const int& lda,
int* ilo,
int* ihi,
double* scale,
int* info)
const;
1583 void GEBAK(
const char& JOBZ,
const char& SIDE,
const int& n,
const int& ilo,
const int& ihi,
const double* scale,
const int& m,
double* V,
const int& ldv,
int* info)
const;
1586 double LARND(
const int& idist,
int* seed )
const;
1587 void LARNV(
const int& idist,
int* seed,
const int& n,
double* v )
const;
1590 double LAMCH(
const char& CMACH)
const;
1591 int ILAENV(
const int& ispec,
const std::string& NAME,
const std::string& OPTS,
const int& N1 = -1,
const int& N2 = -1,
const int& N3 = -1,
const int& N4 = -1 )
const;
1594 double LAPY2(
const double& x,
const double& y)
const;
1600 #ifdef HAVE_TEUCHOS_COMPLEX
1605 class TEUCHOSNUMERICS_LIB_DLL_EXPORT LAPACK<int, std::complex<float> >
1608 inline LAPACK(
void) {}
1609 inline LAPACK(
const LAPACK<
int, std::complex<float> >& lapack) {}
1610 inline virtual ~LAPACK(
void) {}
1613 void PTTRF(
const int& n,
float* d, std::complex<float>* e,
int* info)
const;
1614 void PTTRS(
const char& UPLO,
const int& n,
const int& nrhs,
const float* d,
const std::complex<float>* e, std::complex<float>* B,
const int& ldb,
int* info)
const;
1615 void POTRF(
const char& UPLO,
const int& n, std::complex<float>* A,
const int& lda,
int* info)
const;
1616 void POTRS(
const char& UPLO,
const int& n,
const int& nrhs,
const std::complex<float>* A,
const int& lda, std::complex<float>* B,
const int& ldb,
int* info)
const;
1617 void POTRI(
const char& UPLO,
const int& n, std::complex<float>* A,
const int& lda,
int* info)
const;
1618 void POCON(
const char& UPLO,
const int& n,
const std::complex<float>* A,
const int& lda,
const float& anorm,
float* rcond, std::complex<float>* WORK,
float* rwork,
int* info)
const;
1619 void POSV(
const char& UPLO,
const int& n,
const int& nrhs, std::complex<float>* A,
const int& lda, std::complex<float>* B,
const int& ldb,
int* info)
const;
1620 void POEQU(
const int& n,
const std::complex<float>* A,
const int& lda,
float* S,
float* scond,
float* amax,
int* info)
const;
1621 void PORFS(
const char& UPLO,
const int& n,
const int& nrhs,
const std::complex<float>* A,
const int& lda,
const std::complex<float>* AF,
const int& ldaf,
const std::complex<float>* B,
const int& ldb, std::complex<float>* X,
const int& ldx,
float* FERR,
float* BERR, std::complex<float>* WORK,
float* RWORK,
int* info)
const;
1623 void POSVX(
const char& FACT,
const char& UPLO,
const int& n,
const int& nrhs, std::complex<float>* A,
const int& lda, std::complex<float>* AF,
const int& ldaf,
char* EQUED,
float* S, std::complex<float>* B,
const int& ldb, std::complex<float>* X,
const int& ldx,
float* rcond,
float* FERR,
float* BERR, std::complex<float>* WORK,
float* RWORK,
int* info)
const;
1626 void GELS(
const char&
TRANS,
const int& m,
const int& n,
const int& nrhs, std::complex<float>* A,
const int& lda, std::complex<float>* B,
const int& ldb, std::complex<float>* WORK,
const int& lwork,
int* info)
const;
1627 void GELSS(
const int& m,
const int& n,
const int& nrhs, std::complex<float>* A,
const int& lda, std::complex<float>* B,
const int& ldb,
float* S,
const float& rcond,
int* rank, std::complex<float>* WORK,
const int& lwork,
float* RWORK,
int* info)
const;
1628 void GEQRF(
const int& m,
const int& n, std::complex<float>* A,
const int& lda, std::complex<float>* TAU, std::complex<float>* WORK,
const int& lwork,
int* info)
const;
1629 void GEQR2(
const int& m,
const int& n, std::complex<float>* A,
const int& lda, std::complex<float>* TAU, std::complex<float>* WORK,
int*
const info)
const;
1630 void UNGQR(
const int& m,
const int& n,
const int& k, std::complex<float>* A,
const int& lda,
const std::complex<float>* TAU, std::complex<float>* WORK,
const int& lwork,
int* info)
const;
1631 void UNMQR(
const char& SIDE,
const char&
TRANS,
const int& m,
const int& n,
const int& k,
const std::complex<float>* A,
const int& lda,
const std::complex<float>* TAU, std::complex<float>* C,
const int& ldc, std::complex<float>* WORK,
const int& lwork,
int* info)
const;
1632 void UNM2R(
const char& SIDE,
const char&
TRANS,
const int& M,
const int& N,
const int& K,
const std::complex<float>* A,
const int& LDA,
const std::complex<float>* TAU, std::complex<float>* C,
const int& LDC, std::complex<float>* WORK,
int*
const INFO)
const;
1633 void GETRF(
const int& m,
const int& n, std::complex<float>* A,
const int& lda,
int* IPIV,
int* info)
const;
1634 void GETRS(
const char&
TRANS,
const int& n,
const int& nrhs,
const std::complex<float>* A,
const int& lda,
const int* IPIV, std::complex<float>* B,
const int& ldb,
int* info)
const;
1635 void LASCL(
const char& TYPE,
const int& kl,
const int& ku,
const float& cfrom,
const float& cto,
const int& m,
const int& n, std::complex<float>* A,
const int& lda,
int* info)
const;
1637 void GEQP3 (
const int& m,
const int& n, std::complex<float>* A,
const int& lda,
int* jpvt, std::complex<float>* TAU, std::complex<float>* WORK,
const int& lwork,
float* RWORK,
int* info)
const;
1638 void LASWP (
const int& N, std::complex<float>* A,
const int& LDA,
const int& K1,
const int& K2,
const int* IPIV,
const int& INCX)
const;
1640 void GBTRF(
const int& m,
const int& n,
const int& kl,
const int& ku, std::complex<float>* A,
const int& lda,
int* IPIV,
int* info)
const;
1641 void GBTRS(
const char&
TRANS,
const int& n,
const int& kl,
const int& ku,
const int& nrhs,
const std::complex<float>* A,
const int& lda,
const int* IPIV, std::complex<float>* B,
const int& ldb,
int* info)
const;
1642 void GTTRF(
const int& n, std::complex<float>* dl, std::complex<float>* d, std::complex<float>* du, std::complex<float>* du2,
int* IPIV,
int* info)
const;
1643 void GTTRS(
const char&
TRANS,
const int& n,
const int& nrhs,
const std::complex<float>* dl,
const std::complex<float>* d,
const std::complex<float>* du,
const std::complex<float>* du2,
const int* IPIV, std::complex<float>* B,
const int& ldb,
int* info)
const;
1644 void GETRI(
const int& n, std::complex<float>* A,
const int& lda,
const int* IPIV, std::complex<float>* WORK,
const int& lwork,
int* info)
const;
1645 void LATRS (
const char& UPLO,
const char&
TRANS,
const char& DIAG,
const char& NORMIN,
const int& N,
const std::complex<float>* A,
const int& LDA, std::complex<float>* X,
float* SCALE,
float* CNORM,
int* INFO)
const;
1646 void GECON(
const char& NORM,
const int& n,
const std::complex<float>* A,
const int& lda,
const float& anorm,
float* rcond, std::complex<float>* WORK,
float* RWORK,
int* info)
const;
1647 void GBCON(
const char& NORM,
const int& n,
const int& kl,
const int& ku,
const std::complex<float>* A,
const int& lda,
const int* IPIV,
const float& anorm,
float* rcond, std::complex<float>* WORK,
float* RWORK,
int* info)
const;
1648 float LANGB(
const char& NORM,
const int& n,
const int& kl,
const int& ku,
const std::complex<float>* A,
const int& lda,
float* WORK)
const;
1649 void GESV(
const int& n,
const int& nrhs, std::complex<float>* A,
const int& lda,
int* IPIV, std::complex<float>* B,
const int& ldb,
int* info)
const;
1650 void GEEQU(
const int& m,
const int& n,
const std::complex<float>* A,
const int& lda,
float* R,
float* C,
float* rowcond,
float* colcond,
float* amax,
int* info)
const;
1651 void GERFS(
const char&
TRANS,
const int& n,
const int& nrhs,
const std::complex<float>* A,
const int& lda,
const std::complex<float>* AF,
const int& ldaf,
const int* IPIV,
const std::complex<float>* B,
const int& ldb, std::complex<float>* X,
const int& ldx,
float* FERR,
float* BERR, std::complex<float>* WORK,
float* RWORK,
int* info)
const;
1652 void GBEQU(
const int& m,
const int& n,
const int& kl,
const int& ku,
const std::complex<float>* A,
const int& lda,
float* R,
float* C,
float* rowcond,
float* colcond,
float* amax,
int* info)
const;
1653 void GBRFS(
const char&
TRANS,
const int& n,
const int& kl,
const int& ku,
const int& nrhs,
const std::complex<float>* A,
const int& lda,
const std::complex<float>* AF,
const int& ldaf,
const int* IPIV,
const std::complex<float>* B,
const int& ldb, std::complex<float>* X,
const int& ldx,
float* FERR,
float* BERR, std::complex<float>* WORK,
float* RWORK,
int* info)
const;
1655 void GESVX(
const char& FACT,
const char&
TRANS,
const int& n,
const int& nrhs, std::complex<float>* A,
const int& lda, std::complex<float>* AF,
const int& ldaf,
int* IPIV,
char* EQUED,
float* R,
float* C, std::complex<float>* B,
const int& ldb, std::complex<float>* X,
const int& ldx,
float* rcond,
float* FERR,
float* BERR, std::complex<float>* WORK,
float* RWORK,
int* info)
const;
1657 void GEHRD(
const int& n,
const int& ilo,
const int& ihi, std::complex<float>* A,
const int& lda, std::complex<float>* TAU, std::complex<float>* WORK,
const int& lwork,
int* info)
const;
1658 void TRTRS(
const char& UPLO,
const char&
TRANS,
const char& DIAG,
const int& n,
const int& nrhs,
const std::complex<float>* A,
const int& lda, std::complex<float>* B,
const int& ldb,
int* info)
const;
1659 void TRTRI(
const char& UPLO,
const char& DIAG,
const int& n, std::complex<float>* A,
const int& lda,
int* info)
const;
1662 void STEQR(
const char& COMPZ,
const int& n,
float* D,
float* E, std::complex<float>* Z,
const int& ldz,
float* WORK,
int* info)
const;
1663 void PTEQR(
const char& COMPZ,
const int& n,
float* D,
float* E, std::complex<float>* Z,
const int& ldz,
float* WORK,
int* info)
const;
1664 void HEEV(
const char& JOBZ,
const char& UPLO,
const int& n, std::complex<float>* A,
const int& lda,
float* W, std::complex<float>* WORK,
const int& lwork,
float* RWORK,
int* info)
const;
1665 void HEGV(
const int& itype,
const char& JOBZ,
const char& UPLO,
const int& n, std::complex<float>* A,
const int& lda, std::complex<float>* B,
const int& ldb,
float* W, std::complex<float>* WORK,
const int& lwork,
float* RWORK,
int* info)
const;
1668 void HSEQR(
const char& JOB,
const char& COMPZ,
const int& n,
const int& ilo,
const int& ihi, std::complex<float>* H,
const int& ldh, std::complex<float>* W, std::complex<float>* Z,
const int& ldz, std::complex<float>* WORK,
const int& lwork,
int* info)
const;
1669 void GEES(
const char& JOBVS,
const char& SORT,
int (*ptr2func)(std::complex<float>*),
const int& n, std::complex<float>* A,
const int& lda,
int* sdim, std::complex<float>* W, std::complex<float>* VS,
const int& ldvs, std::complex<float>* WORK,
const int& lwork,
float* RWORK,
int* BWORK,
int* info)
const;
1670 void GEES(
const char& JOBVS,
const int& n, std::complex<float>* A,
const int& lda,
int* sdim,
float* WR,
float* WI, std::complex<float>* VS,
const int& ldvs, std::complex<float>* WORK,
const int& lwork,
float* RWORK,
int* BWORK,
int* info)
const;
1672 void GEEV(
const char& JOBVL,
const char& JOBVR,
const int& n, std::complex<float>* A,
const int& lda, std::complex<float>* W, std::complex<float>* VL,
const int& ldvl, std::complex<float>* VR,
const int& ldvr, std::complex<float>* WORK,
const int& lwork,
float* RWORK,
int* info)
const;
1673 void GEEV(
const char& JOBVL,
const char& JOBVR,
const int& n, std::complex<float>* A,
const int& lda,
float* WR,
float* WI, std::complex<float>* VL,
const int& ldvl, std::complex<float>* VR,
const int& ldvr, std::complex<float>* WORK,
const int& lwork,
float* RWORK,
int* info)
const;
1675 void GEEVX(
const char& BALANC,
const char& JOBVL,
const char& JOBVR,
const char& SENSE,
const int& n, std::complex<float>* A,
const int& lda, std::complex<float>* W, std::complex<float>* VL,
const int& ldvl, std::complex<float>* VR,
const int& ldvr,
int* ilo,
int* ihi,
float* SCALE,
float* abnrm,
float* RCONDE,
float* RCONDV, std::complex<float>* WORK,
const int& lwork,
float* RWORK,
int* info)
const;
1677 void GGEVX(
const char& BALANC,
const char& JOBVL,
const char& JOBVR,
const char& SENSE,
const int& n, std::complex<float>* A,
const int& lda, std::complex<float>* B,
const int& ldb, std::complex<float>* ALPHA, std::complex<float>* BETA, std::complex<float>* VL,
const int& ldvl, std::complex<float>* VR,
const int& ldvr,
int* ilo,
int* ihi,
float* lscale,
float* rscale,
float* abnrm,
float* bbnrm,
float* RCONDE,
float* RCONDV, std::complex<float>* WORK,
const int& lwork,
float* RWORK,
int* IWORK,
int* BWORK,
int* info)
const;
1678 void GGEVX(
const char& BALANC,
const char& JOBVL,
const char& JOBVR,
const char& SENSE,
const int& n, std::complex<float>* A,
const int& lda, std::complex<float>* B,
const int& ldb,
float* ALPHAR,
float* ALPHAI, std::complex<float>* BETA, std::complex<float>* VL,
const int& ldvl, std::complex<float>* VR,
const int& ldvr,
int* ilo,
int* ihi,
float* lscale,
float* rscale,
float* abnrm,
float* bbnrm,
float* RCONDE,
float* RCONDV, std::complex<float>* WORK,
const int& lwork,
float* RWORK,
int* IWORK,
int* BWORK,
int* info)
const;
1679 void GGEV(
const char& JOBVL,
const char& JOBVR,
const int& n, std::complex<float> *A,
const int& lda, std::complex<float> *B,
const int& ldb, std::complex<float>* ALPHA, std::complex<float>* BETA, std::complex<float>* VL,
const int& ldvl, std::complex<float>* VR,
const int& ldvr, std::complex<float> *WORK,
const int& lwork,
float* RWORK,
int* info)
const;
1680 void GGES(
const char& JOBVL,
const char& JOBVR,
const char& SORT,
int (*ptr2func)(std::complex<float>*, std::complex<float>*),
const int& n, std::complex<float>* A,
const int& lda, std::complex<float>* B,
const int& ldb,
int* sdim, std::complex<float>* ALPHA, std::complex<float>* BETA, std::complex<float>* VL,
const int& ldvl, std::complex<float>* VR,
const int& ldvr, std::complex<float>* WORK,
const int& lwork,
float* rwork,
int* bwork,
int* info )
const;
1681 void TGSEN(
const int& ijob,
const int& wantq,
const int& wantz,
const int* SELECT,
const int& n, std::complex<float>* A,
const int& lda, std::complex<float>* B,
const int& ldb, std::complex<float>* ALPHA, std::complex<float>* BETA, std::complex<float>* Q,
const int& ldq, std::complex<float>* Z,
const int& ldz,
int* M,
float* PL,
float* PR,
float* DIF, std::complex<float>* WORK,
const int& lwork,
int* IWORK,
const int& liwork,
int* info )
const;
1684 void GESVD(
const char& JOBU,
const char& JOBVT,
const int& m,
const int& n, std::complex<float>* A,
const int& lda,
float* S, std::complex<float>* U,
const int& ldu, std::complex<float>* V,
const int& ldv, std::complex<float>* WORK,
const int& lwork,
float* RWORK,
int* info)
const;
1687 void TREVC(
const char& SIDE,
const char& HOWMNY,
int* select,
const int& n,
const std::complex<float>* T,
const int& ldt, std::complex<float>* VL,
const int& ldvl, std::complex<float>* VR,
const int& ldvr,
const int& mm,
int* m, std::complex<float>* WORK,
float* RWORK,
int* info)
const;
1688 void TREVC(
const char& SIDE,
const int& n,
const std::complex<float>* T,
const int& ldt, std::complex<float>* VL,
const int& ldvl, std::complex<float>* VR,
const int& ldvr,
const int& mm,
int* m, std::complex<float>* WORK,
float* RWORK,
int* info)
const;
1690 void TREXC(
const char& COMPQ,
const int& n, std::complex<float>* T,
const int& ldt, std::complex<float>* Q,
const int& ldq,
int* ifst,
int* ilst, std::complex<float>* WORK,
int* info)
const;
1693 void LARTG(
const std::complex<float> f,
const std::complex<float> g,
float* c, std::complex<float>* s, std::complex<float>* r )
const;
1694 void LARFG(
const int& n, std::complex<float>* alpha, std::complex<float>* x,
const int& incx, std::complex<float>* tau )
const;
1698 void GEBAL(
const char& JOBZ,
const int& n, std::complex<float>* A,
const int& lda,
int* ilo,
int* ihi,
float* scale,
int* info)
const;
1700 void GEBAK(
const char& JOBZ,
const char& SIDE,
const int& n,
const int& ilo,
const int& ihi,
const float* scale,
const int& m, std::complex<float>* V,
const int& ldv,
int* info)
const;
1703 std::complex<float> LARND(
const int& idist,
int* seed )
const;
1704 void LARNV(
const int& idist,
int* seed,
const int& n, std::complex<float>* v )
const;
1707 int ILAENV(
const int& ispec,
const std::string& NAME,
const std::string& OPTS,
const int& N1 = -1,
const int& N2 = -1,
const int& N3 = -1,
const int& N4 = -1 )
const;
1716 class TEUCHOSNUMERICS_LIB_DLL_EXPORT LAPACK<int, std::complex<double> >
1719 inline LAPACK(
void) {}
1720 inline LAPACK(
const LAPACK<
int, std::complex<double> >& lapack) {}
1721 inline virtual ~LAPACK(
void) {}
1724 void PTTRF(
const int& n,
double* d, std::complex<double>* e,
int* info)
const;
1725 void PTTRS(
const char& UPLO,
const int& n,
const int& nrhs,
const double* d,
const std::complex<double>* e, std::complex<double>* B,
const int& ldb,
int* info)
const;
1726 void POTRF(
const char& UPLO,
const int& n, std::complex<double>* A,
const int& lda,
int* info)
const;
1727 void POTRS(
const char& UPLO,
const int& n,
const int& nrhs,
const std::complex<double>* A,
const int& lda, std::complex<double>* B,
const int& ldb,
int* info)
const;
1728 void POTRI(
const char& UPLO,
const int& n, std::complex<double>* A,
const int& lda,
int* info)
const;
1729 void POCON(
const char& UPLO,
const int& n,
const std::complex<double>* A,
const int& lda,
const double& anorm,
double* rcond, std::complex<double>* WORK,
double* RWORK,
int* info)
const;
1730 void POSV(
const char& UPLO,
const int& n,
const int& nrhs, std::complex<double>* A,
const int& lda, std::complex<double>* B,
const int& ldb,
int* info)
const;
1731 void POEQU(
const int& n,
const std::complex<double>* A,
const int& lda,
double* S,
double* scond,
double* amax,
int* info)
const;
1732 void PORFS(
const char& UPLO,
const int& n,
const int& nrhs,
const std::complex<double>* A,
const int& lda,
const std::complex<double>* AF,
const int& ldaf,
const std::complex<double>* B,
const int& ldb, std::complex<double>* X,
const int& ldx,
double* FERR,
double* BERR, std::complex<double>* WORK,
double* RWORK,
int* info)
const;
1734 void POSVX(
const char& FACT,
const char& UPLO,
const int& n,
const int& nrhs, std::complex<double>* A,
const int& lda, std::complex<double>* AF,
const int& ldaf,
char* EQUED,
double* S, std::complex<double>* B,
const int& ldb, std::complex<double>* X,
const int& ldx,
double* rcond,
double* FERR,
double* BERR, std::complex<double>* WORK,
double* RWORK,
int* info)
const;
1737 void GELS(
const char&
TRANS,
const int& m,
const int& n,
const int& nrhs, std::complex<double>* A,
const int& lda, std::complex<double>* B,
const int& ldb, std::complex<double>* WORK,
const int& lwork,
int* info)
const;
1738 void GELSS(
const int& m,
const int& n,
const int& nrhs, std::complex<double>* A,
const int& lda, std::complex<double>* B,
const int& ldb,
double* S,
const double& rcond,
int* rank, std::complex<double>* WORK,
const int& lwork,
double* RWORK,
int* info)
const;
1739 void GEQRF(
const int& m,
const int& n, std::complex<double>* A,
const int& lda, std::complex<double>* TAU, std::complex<double>* WORK,
const int& lwork,
int* info)
const;
1740 void GEQR2(
const int& m,
const int& n, std::complex<double>* A,
const int& lda, std::complex<double>* TAU, std::complex<double>* WORK,
int*
const info)
const;
1741 void UNGQR(
const int& m,
const int& n,
const int& k, std::complex<double>* A,
const int& lda,
const std::complex<double>* TAU, std::complex<double>* WORK,
const int& lwork,
int* info)
const;
1742 void UNMQR(
const char& SIDE,
const char&
TRANS,
const int& m,
const int& n,
const int& k,
const std::complex<double>* A,
const int& lda,
const std::complex<double>* TAU, std::complex<double>* C,
const int& ldc, std::complex<double>* WORK,
const int& lwork,
int* info)
const;
1743 void UNM2R(
const char& SIDE,
const char&
TRANS,
const int& M,
const int& N,
const int& K,
const std::complex<double>* A,
const int& LDA,
const std::complex<double>* TAU, std::complex<double>* C,
const int& LDC, std::complex<double>* WORK,
int*
const INFO)
const;
1745 void GETRF(
const int& m,
const int& n, std::complex<double>* A,
const int& lda,
int* IPIV,
int* info)
const;
1746 void GETRS(
const char&
TRANS,
const int& n,
const int& nrhs,
const std::complex<double>* A,
const int& lda,
const int* IPIV, std::complex<double>* B,
const int& ldb,
int* info)
const;
1747 void LASCL(
const char& TYPE,
const int& kl,
const int& ku,
const double& cfrom,
const double& cto,
const int& m,
const int& n, std::complex<double>* A,
const int& lda,
int* info)
const;
1749 void GEQP3 (
const int& m,
const int& n, std::complex<double>* A,
const int& lda,
int* jpvt, std::complex<double>* TAU, std::complex<double>* WORK,
const int& lwork,
double* RWORK,
int* info)
const;
1750 void LASWP (
const int& N, std::complex<double>* A,
const int& LDA,
const int& K1,
const int& K2,
const int* IPIV,
const int& INCX)
const;
1752 void GBTRF(
const int& m,
const int& n,
const int& kl,
const int& ku, std::complex<double>* A,
const int& lda,
int* IPIV,
int* info)
const;
1753 void GBTRS(
const char&
TRANS,
const int& n,
const int& kl,
const int& ku,
const int& nrhs,
const std::complex<double>* A,
const int& lda,
const int* IPIV, std::complex<double>* B,
const int& ldb,
int* info)
const;
1754 void GTTRF(
const int& n, std::complex<double>* dl, std::complex<double>* d, std::complex<double>* du, std::complex<double>* du2,
int* IPIV,
int* info)
const;
1755 void GTTRS(
const char&
TRANS,
const int& n,
const int& nrhs,
const std::complex<double>* dl,
const std::complex<double>* d,
const std::complex<double>* du,
const std::complex<double>* du2,
const int* IPIV, std::complex<double>* B,
const int& ldb,
int* info)
const;
1756 void GETRI(
const int& n, std::complex<double>* A,
const int& lda,
const int* IPIV, std::complex<double>* WORK,
const int& lwork,
int* info)
const;
1757 void LATRS (
const char& UPLO,
const char&
TRANS,
const char& DIAG,
const char& NORMIN,
const int& N,
const std::complex<double>* A,
const int& LDA, std::complex<double>* X,
double* SCALE,
double* CNORM,
int* INFO)
const;
1758 void GECON(
const char& NORM,
const int& n,
const std::complex<double>* A,
const int& lda,
const double& anorm,
double* rcond, std::complex<double>* WORK,
double* RWORK,
int* info)
const;
1759 void GBCON(
const char& NORM,
const int& n,
const int& kl,
const int& ku,
const std::complex<double>* A,
const int& lda,
const int* IPIV,
const double& anorm,
double* rcond, std::complex<double>* WORK,
double* RWORK,
int* info)
const;
1760 double LANGB(
const char& NORM,
const int& n,
const int& kl,
const int& ku,
const std::complex<double>* A,
const int& lda,
double* WORK)
const;
1761 void GESV(
const int& n,
const int& nrhs, std::complex<double>* A,
const int& lda,
int* IPIV, std::complex<double>* B,
const int& ldb,
int* info)
const;
1762 void GEEQU(
const int& m,
const int& n,
const std::complex<double>* A,
const int& lda,
double* R,
double* C,
double* rowcond,
double* colcond,
double* amax,
int* info)
const;
1763 void GERFS(
const char&
TRANS,
const int& n,
const int& nrhs,
const std::complex<double>* A,
const int& lda,
const std::complex<double>* AF,
const int& ldaf,
const int* IPIV,
const std::complex<double>* B,
const int& ldb, std::complex<double>* X,
const int& ldx,
double* FERR,
double* BERR, std::complex<double>* WORK,
double* RWORK,
int* info)
const;
1764 void GBEQU(
const int& m,
const int& n,
const int& kl,
const int& ku,
const std::complex<double>* A,
const int& lda,
double* R,
double* C,
double* rowcond,
double* colcond,
double* amax,
int* info)
const;
1765 void GBRFS(
const char&
TRANS,
const int& n,
const int& kl,
const int& ku,
const int& nrhs,
const std::complex<double>* A,
const int& lda,
const std::complex<double>* AF,
const int& ldaf,
const int* IPIV,
const std::complex<double>* B,
const int& ldb, std::complex<double>* X,
const int& ldx,
double* FERR,
double* BERR, std::complex<double>* WORK,
double* RWORK,
int* info)
const;
1767 void GESVX(
const char& FACT,
const char&
TRANS,
const int& n,
const int& nrhs, std::complex<double>* A,
const int& lda, std::complex<double>* AF,
const int& ldaf,
int* IPIV,
char* EQUED,
double* R,
double* C, std::complex<double>* B,
const int& ldb, std::complex<double>* X,
const int& ldx,
double* rcond,
double* FERR,
double* BERR, std::complex<double>* WORK,
double* RWORK,
int* info)
const;
1769 void GEHRD(
const int& n,
const int& ilo,
const int& ihi, std::complex<double>* A,
const int& lda, std::complex<double>* TAU, std::complex<double>* WORK,
const int& lwork,
int* info)
const;
1770 void TRTRS(
const char& UPLO,
const char&
TRANS,
const char& DIAG,
const int& n,
const int& nrhs,
const std::complex<double>* A,
const int& lda, std::complex<double>* B,
const int& ldb,
int* info)
const;
1771 void TRTRI(
const char& UPLO,
const char& DIAG,
const int& n, std::complex<double>* A,
const int& lda,
int* info)
const;
1774 void STEQR(
const char& COMPZ,
const int& n,
double* D,
double* E, std::complex<double>* Z,
const int& ldz,
double* WORK,
int* info)
const;
1775 void PTEQR(
const char& COMPZ,
const int& n,
double* D,
double* E, std::complex<double>* Z,
const int& ldz,
double* WORK,
int* info)
const;
1776 void HEEV(
const char& JOBZ,
const char& UPLO,
const int& n, std::complex<double>* A,
const int& lda,
double* W, std::complex<double>* WORK,
const int& lwork,
double* RWORK,
int* info)
const;
1777 void HEGV(
const int& itype,
const char& JOBZ,
const char& UPLO,
const int& n, std::complex<double>* A,
const int& lda, std::complex<double>* B,
const int& ldb,
double* W, std::complex<double>* WORK,
const int& lwork,
double* RWORK,
int* info)
const;
1780 void HSEQR(
const char& JOB,
const char& COMPZ,
const int& n,
const int& ilo,
const int& ihi, std::complex<double>* H,
const int& ldh, std::complex<double>* W, std::complex<double>* Z,
const int& ldz, std::complex<double>* WORK,
const int& lwork,
int* info)
const;
1781 void GEES(
const char& JOBVS,
const char& SORT,
int (*ptr2func)(std::complex<double>*),
const int& n, std::complex<double>* A,
const int& lda,
int* sdim, std::complex<double>* W, std::complex<double>* VS,
const int& ldvs, std::complex<double>* WORK,
const int& lwork,
double* RWORK,
int* BWORK,
int* info)
const;
1782 void GEES(
const char& JOBVS,
const int& n, std::complex<double>* A,
const int& lda,
int* sdim,
double* WR,
double* WI, std::complex<double>* VS,
const int& ldvs, std::complex<double>* WORK,
const int& lwork,
double* RWORK,
int* BWORK,
int* info)
const;
1784 void GEEV(
const char& JOBVL,
const char& JOBVR,
const int& n, std::complex<double>* A,
const int& lda, std::complex<double>* W, std::complex<double>* VL,
const int& ldvl, std::complex<double>* VR,
const int& ldvr, std::complex<double>* WORK,
const int& lwork,
double* RWORK,
int* info)
const;
1785 void GEEV(
const char& JOBVL,
const char& JOBVR,
const int& n, std::complex<double>* A,
const int& lda,
double* WR,
double* WI, std::complex<double>* VL,
const int& ldvl, std::complex<double>* VR,
const int& ldvr, std::complex<double>* WORK,
const int& lwork,
double* RWORK,
int* info)
const;
1787 void GEEVX(
const char& BALANC,
const char& JOBVL,
const char& JOBVR,
const char& SENSE,
const int& n, std::complex<double>* A,
const int& lda, std::complex<double>* W, std::complex<double>* VL,
const int& ldvl, std::complex<double>* VR,
const int& ldvr,
int* ilo,
int* ihi,
double* SCALE,
double* abnrm,
double* RCONDE,
double* RCONDV, std::complex<double>* WORK,
const int& lwork,
double* RWORK,
int* info)
const;
1788 void GGEVX(
const char& BALANC,
const char& JOBVL,
const char& JOBVR,
const char& SENSE,
const int& n, std::complex<double>* A,
const int& lda, std::complex<double>* B,
const int& ldb, std::complex<double>* ALPHA, std::complex<double>* BETA, std::complex<double>* VL,
const int& ldvl, std::complex<double>* VR,
const int& ldvr,
int* ilo,
int* ihi,
double* lscale,
double* rscale,
double* abnrm,
double* bbnrm,
double* RCONDE,
double* RCONDV, std::complex<double>* work,
const int& lwork,
double* RWORK,
int* IWORK,
int* BWORK,
int* info)
const;
1789 void GGEVX(
const char& BALANC,
const char& JOBVL,
const char& JOBVR,
const char& SENSE,
const int& n, std::complex<double>* A,
const int& lda, std::complex<double>* B,
const int& ldb,
double* ALPHAR,
double* ALPHAI, std::complex<double>* BETA, std::complex<double>* VL,
const int& ldvl, std::complex<double>* VR,
const int& ldvr,
int* ilo,
int* ihi,
double* lscale,
double* rscale,
double* abnrm,
double* bbnrm,
double* RCONDE,
double* RCONDV, std::complex<double>* work,
const int& lwork,
double* RWORK,
int* IWORK,
int* BWORK,
int* info)
const;
1790 void GGEV(
const char& JOBVL,
const char& JOBVR,
const int& n, std::complex<double> *A,
const int& lda, std::complex<double> *B,
const int& ldb, std::complex<double>* ALPHA, std::complex<double>* BETA, std::complex<double>* VL,
const int& ldvl, std::complex<double>*VR,
const int& ldvr, std::complex<double> *WORK,
const int& lwork,
double* RWORK,
int* info)
const;
1791 void GGES(
const char& JOBVL,
const char& JOBVR,
const char& SORT,
int (*ptr2func)(std::complex<double>*, std::complex<double>*),
const int& n, std::complex<double>* A,
const int& lda, std::complex<double>* B,
const int& ldb,
int* sdim, std::complex<double>* ALPHA, std::complex<double>* BETA, std::complex<double>* VL,
const int& ldvl, std::complex<double>* VR,
const int& ldvr, std::complex<double>* WORK,
const int& lwork,
double* rwork,
int* bwork,
int* info )
const;
1792 void TGSEN(
const int& ijob,
const int& wantq,
const int& wantz,
const int* SELECT,
const int& n, std::complex<double>* A,
const int& lda, std::complex<double>* B,
const int& ldb, std::complex<double>* ALPHA, std::complex<double>* BETA, std::complex<double>* Q,
const int& ldq, std::complex<double>* Z,
const int& ldz,
int* M,
double* PL,
double* PR,
double* DIF, std::complex<double>* WORK,
const int& lwork,
int* IWORK,
const int& liwork,
int* info )
const;
1795 void GESVD(
const char& JOBU,
const char& JOBVT,
const int& m,
const int& n, std::complex<double>* A,
const int& lda,
double* S, std::complex<double>* U,
const int& ldu, std::complex<double>* V,
const int& ldv, std::complex<double>* WORK,
const int& lwork,
double* RWORK,
int* info)
const;
1798 void TREVC(
const char& SIDE,
const char& HOWMNY,
int* select,
const int& n,
const std::complex<double>* T,
const int& ldt, std::complex<double>* VL,
const int& ldvl, std::complex<double>* VR,
const int& ldvr,
const int& mm,
int* m, std::complex<double>* WORK,
double* RWORK,
int* info)
const;
1799 void TREVC(
const char& SIDE,
const int& n,
const std::complex<double>* T,
const int& ldt, std::complex<double>* VL,
const int& ldvl, std::complex<double>* VR,
const int& ldvr,
const int& mm,
int* m, std::complex<double>* WORK,
double* RWORK,
int* info)
const;
1801 void TREXC(
const char& COMPQ,
const int& n, std::complex<double>* T,
const int& ldt, std::complex<double>* Q,
const int& ldq,
int* ifst,
int* ilst, std::complex<double>* WORK,
int* info)
const;
1804 void LARTG(
const std::complex<double> f,
const std::complex<double> g,
double* c, std::complex<double>* s, std::complex<double>* r )
const;
1805 void LARFG(
const int& n, std::complex<double>* alpha, std::complex<double>* x,
const int& incx, std::complex<double>* tau )
const;
1809 void GEBAL(
const char& JOBZ,
const int& n, std::complex<double>* A,
const int& lda,
int* ilo,
int* ihi,
double* scale,
int* info)
const;
1811 void GEBAK(
const char& JOBZ,
const char& SIDE,
const int& n,
const int& ilo,
const int& ihi,
const double* scale,
const int& m, std::complex<double>* V,
const int& ldv,
int* info)
const;
1814 std::complex<double> LARND(
const int& idist,
int* seed )
const;
1815 void LARNV(
const int& idist,
int* seed,
const int& n, std::complex<double>* v )
const;
1818 int ILAENV(
const int& ispec,
const std::string& NAME,
const std::string& OPTS,
const int& N1 = -1,
const int& N2 = -1,
const int& N3 = -1,
const int& N4 = -1 )
const;
1824 #endif // HAVE_TEUCHOS_COMPLEX
1826 #ifdef HAVE_TEUCHOSCORE_QUADMATH
1833 class TEUCHOSNUMERICS_LIB_DLL_EXPORT LAPACK<int, __float128>
1836 inline LAPACK(
void) {}
1837 inline LAPACK(
const LAPACK<int, __float128>& lapack) {}
1838 inline virtual ~LAPACK(
void) {}
1840 void GEQRF(
const int& m,
const int& n, __float128* A,
const int& lda, __float128* TAU, __float128* WORK,
const int& lwork,
int* info)
const;
1841 void GEQR2(
const int& m,
const int& n, __float128* A,
const int& lda, __float128* TAU, __float128* WORK,
int*
const info)
const;
1842 void GETRF(
const int& m,
const int& n, __float128* A,
const int& lda,
int* IPIV,
int* info)
const;
1843 void GETRS(
const char&
TRANS,
const int& n,
const int& nrhs,
const __float128* A,
const int& lda,
const int* IPIV, __float128* B,
const int& ldb,
int* info)
const;
1844 void GETRI(
const int& n, __float128* A,
const int& lda,
const int* IPIV, __float128* WORK,
const int& lwork,
int* info)
const;
1845 void LASWP (
const int& N, __float128* A,
const int& LDA,
const int& K1,
const int& K2,
const int* IPIV,
const int& INCX)
const;
1847 void ORM2R(
const char& SIDE,
const char&
TRANS,
const int& m,
const int& n,
const int& k,
const __float128* A,
const int& lda,
const __float128* TAU, __float128* C,
const int& ldc, __float128* WORK,
int*
const info)
const;
1848 void ORGQR(
const int& m,
const int& n,
const int& k, __float128* A,
const int& lda,
const __float128* TAU, __float128* WORK,
const int& lwork,
int* info)
const;
1849 void UNGQR(
const int& m,
const int& n,
const int& k, __float128* A,
const int& lda,
const __float128* TAU, __float128* WORK,
const int& lwork,
int* info)
const;
1851 void LARFG(
const int& n, __float128* alpha, __float128* x,
const int& incx, __float128* tau )
const;
1853 __float128 LAPY2 (
const __float128 x,
const __float128 y)
const;
1854 void LASCL (
const char& TYPE,
const int& kl,
const int& ku,
const __float128 cfrom,
const __float128 cto,
const int& m,
const int& n, __float128* A,
const int& lda,
int* info)
const;
1856 void GBTRF (
const int& m,
const int& n,
const int& kl,
const int& ku, __float128* A,
const int& lda,
int* IPIV,
int* info)
const;
1857 void GBTRS (
const char&
TRANS,
const int& n,
const int& kl,
const int& ku,
const int& nrhs,
const __float128* A,
const int& lda,
const int* IPIV, __float128* B,
const int& ldb,
int* info)
const;
1862 #endif // HAVE_TEUCHOSCORE_QUADMATH
1864 #ifdef HAVE_TEUCHOS_LONG_DOUBLE
1871 class TEUCHOSNUMERICS_LIB_DLL_EXPORT LAPACK<int, long double>
1874 inline LAPACK(
void) {}
1875 inline LAPACK(
const LAPACK<int, long double>& lapack) {}
1876 inline virtual ~LAPACK(
void) {}
1878 void GESV(
const int& n,
const int& nrhs,
long double* A,
const int& lda,
int* IPIV,
long double* B,
const int& ldb,
int* info)
const;
1879 void GTTRS(
const char&
TRANS,
const int& n,
const int& nrhs,
const long double* dl,
const long double* d,
const long double* du,
const long double* du2,
const int* IPIV,
long double* B,
const int& ldb,
int* info)
const;
1880 void GTTRF(
const int& n,
long double* dl,
long double* d,
long double* du,
long double* du2,
int* IPIV,
int* info)
const;
1881 void SYGV(
const int& itype,
const char& JOBZ,
const char& UPLO,
const int& n,
long double* A,
const int& lda,
long double* B,
const int& ldb,
long double* W,
long double* WORK,
const int& lwork,
int* info)
const;
1882 void GEEV(
const char& JOBVL,
const char& JOBVR,
const int& n,
long double* A,
const int& lda,
long double* WR,
long double* WI,
long double* VL,
const int& ldvl,
long double* VR,
const int& ldvr,
long double* WORK,
const int& lwork,
int* info)
const;
1883 void GEEV(
const char& JOBVL,
const char& JOBVR,
const int& n,
long double* A,
const int& lda,
long double* WR,
long double* WI,
long double* VL,
const int& ldvl,
long double* VR,
const int& ldvr,
long double* WORK,
const int& lwork,
long double* ,
int* info)
const;
1884 void GGEVX(
const char& BALANC,
const char& JOBVL,
const char& JOBVR,
const char& SENSE,
const int& n,
long double* A,
const int& lda,
long double* B,
const int& ldb,
long double* ALPHAR,
long double* ALPHAI,
long double* BETA,
long double* VL,
const int& ldvl,
long double* VR,
const int& ldvr,
int* ilo,
int* ihi,
long double* lscale,
long double* rscale,
long double* abnrm,
long double* bbnrm,
long double* RCONDE,
long double* RCONDV,
long double* WORK,
const int& lwork,
int* IWORK,
int* BWORK,
int* info)
const;
1885 void GGEVX(
const char& BALANC,
const char& JOBVL,
const char& JOBVR,
const char& SENSE,
const int& n,
long double* A,
const int& lda,
long double* B,
const int& ldb,
long double* ALPHAR,
long double* ALPHAI,
long double* BETA,
long double* VL,
const int& ldvl,
long double* VR,
const int& ldvr,
int* ilo,
int* ihi,
long double* lscale,
long double* rscale,
long double* abnrm,
long double* bbnrm,
long double* RCONDE,
long double* RCONDV,
long double* WORK,
const int& lwork,
long double* ,
int* IWORK,
int* BWORK,
int* info)
const;
1886 void PORFS(
const char& UPLO,
const int& n,
const int& nrhs,
const long double* A,
const int& lda,
const long double* AF,
const int& ldaf,
const long double* B,
const int& ldb,
long double* X,
const int& ldx,
long double* FERR,
long double* BERR,
long double* WORK,
int* IWORK,
int* info)
const;
1887 void PTEQR(
const char& COMPZ,
const int& n,
long double* D,
long double* E,
long double* Z,
const int& ldz,
long double* WORK,
int* info)
const;
1888 void POTRF(
const char& UPLO,
const int& n,
long double* A,
const int& lda,
int* info)
const;
1889 void POTRS(
const char& UPLO,
const int& n,
const int& nrhs,
const long double* A,
const int& lda,
long double* B,
const int& ldb,
int* info)
const;
1890 void POEQU(
const int& n,
const long double* A,
const int& lda,
long double* S,
long double* scond,
long double* amax,
int* info)
const;
1891 void GEQRF(
const int& m,
const int& n,
long double* A,
const int& lda,
long double* TAU,
long double* WORK,
const int& lwork,
int* info)
const;
1892 void GEQR2(
const int& m,
const int& n,
long double* A,
const int& lda,
long double* TAU,
long double* WORK,
int*
const info)
const;
1893 void GETRF(
const int& m,
const int& n,
long double* A,
const int& lda,
int* IPIV,
int* info)
const;
1894 void GETRS(
const char&
TRANS,
const int& n,
const int& nrhs,
const long double* A,
const int& lda,
const int* IPIV,
long double* B,
const int& ldb,
int* info)
const;
1895 void GETRI(
const int& n,
long double* A,
const int& lda,
const int* IPIV,
long double* WORK,
const int& lwork,
int* info)
const;
1896 void LASWP (
const int& N,
long double* A,
const int& LDA,
const int& K1,
const int& K2,
const int* IPIV,
const int& INCX)
const;
1898 void ORM2R(
const char& SIDE,
const char&
TRANS,
const int& m,
const int& n,
const int& k,
const long double* A,
const int& lda,
const long double* TAU,
long double* C,
const int& ldc,
long double* WORK,
int*
const info)
const;
1899 void ORGQR(
const int& m,
const int& n,
const int& k,
long double* A,
const int& lda,
const long double* TAU,
long double* WORK,
const int& lwork,
int* info)
const;
1900 void UNGQR(
const int& m,
const int& n,
const int& k,
long double* A,
const int& lda,
const long double* TAU,
long double* WORK,
const int& lwork,
int* info)
const;
1902 void LARFG(
const int& n,
long double* alpha,
long double* x,
const int& incx,
long double* tau )
const;
1904 long double LAPY2 (
const long double x,
const long double y)
const;
1905 void LASCL (
const char& TYPE,
const int& kl,
const int& ku,
const long double cfrom,
const long double cto,
const int& m,
const int& n,
long double* A,
const int& lda,
int* info)
const;
1907 void GBTRF (
const int& m,
const int& n,
const int& kl,
const int& ku,
long double* A,
const int& lda,
int* IPIV,
int* info)
const;
1908 void GBTRS (
const char&
TRANS,
const int& n,
const int& kl,
const int& ku,
const int& nrhs,
const long double* A,
const int& lda,
const int* IPIV,
long double* B,
const int& ldb,
int* info)
const;
1913 #endif // HAVE_TEUCHOS_LONG_DOUBLE
1915 #endif // DOXYGEN_SHOULD_SKIP_THIS
1919 #endif // _TEUCHOS_LAPACK_HPP_
void TGEVC(const char &SIDE, const char &HOWMNY, const OrdinalType *SELECT, const OrdinalType &n, const ScalarType *S, const OrdinalType &lds, const ScalarType *P, const OrdinalType &ldp, ScalarType *VL, const OrdinalType &ldvl, ScalarType *VR, const OrdinalType &ldvr, const OrdinalType &mm, OrdinalType *M, ScalarType *WORK, OrdinalType *info) const
void ORM2R(const char &SIDE, const char &TRANS, const OrdinalType &m, const OrdinalType &n, const OrdinalType &k, const ScalarType *A, const OrdinalType &lda, const ScalarType *TAU, ScalarType *C, const OrdinalType &ldc, ScalarType *WORK, OrdinalType *const info) const
BLAS 2 version of ORMQR; known workspace size.
void ORGQR(const OrdinalType &m, const OrdinalType &n, const OrdinalType &k, ScalarType *A, const OrdinalType &lda, const ScalarType *TAU, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *info) const
Compute explicit Q factor from QR factorization (GEQRF) (real case).
void GEEQU(const OrdinalType &m, const OrdinalType &n, const ScalarType *A, const OrdinalType &lda, ScalarType *R, ScalarType *C, ScalarType *rowcond, ScalarType *colcond, ScalarType *amax, OrdinalType *info) const
Computes row and column scalings intended to equilibrate an m by n matrix A and reduce its condition ...
void LARTG(const ScalarType &f, const ScalarType &g, MagnitudeType *c, ScalarType *s, ScalarType *r) const
Gnerates a plane rotation that zeros out the second component of the input vector.
void GGEVX(const char &BALANC, const char &JOBVL, const char &JOBVR, const char &SENSE, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, ScalarType *B, const OrdinalType &ldb, MagnitudeType *ALPHAR, MagnitudeType *ALPHAI, ScalarType *BETA, ScalarType *VL, const OrdinalType &ldvl, ScalarType *VR, const OrdinalType &ldvr, OrdinalType *ilo, OrdinalType *ihi, MagnitudeType *lscale, MagnitudeType *rscale, MagnitudeType *abnrm, MagnitudeType *bbnrm, MagnitudeType *RCONDE, MagnitudeType *RCONDV, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *IWORK, OrdinalType *BWORK, OrdinalType *info) const
T magnitudeType
Mandatory typedef for result of magnitude.
void SYTRD(const char &UPLO, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, ScalarType *D, ScalarType *E, ScalarType *TAU, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *info) const
Reduces a real symmetric matrix A to tridiagonal form by orthogonal similarity transformations.
OrdinalType ILAENV(const OrdinalType &ispec, const std::string &NAME, const std::string &OPTS, const OrdinalType &N1=-1, const OrdinalType &N2=-1, const OrdinalType &N3=-1, const OrdinalType &N4=-1) const
Chooses problem-dependent parameters for the local environment.
void GBCON(const char &NORM, const OrdinalType &n, const OrdinalType &kl, const OrdinalType &ku, const ScalarType *A, const OrdinalType &lda, const OrdinalType *IPIV, const ScalarType &anorm, ScalarType *rcond, ScalarType *WORK, OrdinalType *IWORK, OrdinalType *info) const
Estimates the reciprocal of the condition number of a general banded real matrix A, in either the 1-norm or the infinity-norm, using the LU factorization computed by GETRF.
void GGLSE(const OrdinalType &m, const OrdinalType &n, const OrdinalType &p, ScalarType *A, const OrdinalType &lda, ScalarType *B, const OrdinalType &ldb, ScalarType *C, ScalarType *D, ScalarType *X, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *info) const
Solves the linear equality-constrained least squares (LSE) problem where A is an m by n matrix...
void LARFG(const OrdinalType &n, ScalarType *alpha, ScalarType *x, const OrdinalType &incx, ScalarType *tau) const
Generates an elementary reflector of order n that zeros out the last n-1 components of the input vect...
void GEBAK(const char &JOBZ, const char &SIDE, const OrdinalType &n, const OrdinalType &ilo, const OrdinalType &ihi, const MagnitudeType *scale, const OrdinalType &m, ScalarType *V, const OrdinalType &ldv, OrdinalType *info) const
Forms the left or right eigenvectors of a general matrix that has been balanced by GEBAL by backward ...
void POTRI(const char &UPLO, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, OrdinalType *info) const
Computes the inverse of a real symmetric positive definite matrix A using the Cholesky factorization ...
void GTTRS(const char &TRANS, const OrdinalType &n, const OrdinalType &nrhs, const ScalarType *dl, const ScalarType *d, const ScalarType *du, const ScalarType *du2, const OrdinalType *IPIV, ScalarType *B, const OrdinalType &ldb, OrdinalType *info) const
Solves a system of linear equations A*X=B or A'*X=B or A^H*X=B with a tridiagonal matrix A using the ...
static T squareroot(T x)
Returns a number of magnitudeType that is the square root of this scalar type x.
void SPEV(const char &JOBZ, const char &UPLO, const OrdinalType &n, ScalarType *AP, ScalarType *W, ScalarType *Z, const OrdinalType &ldz, ScalarType *WORK, OrdinalType *info) const
Computes the eigenvalues and, optionally, eigenvectors of a symmetric n by n matrix A in packed stora...
void GTTRF(const OrdinalType &n, ScalarType *dl, ScalarType *d, ScalarType *du, ScalarType *du2, OrdinalType *IPIV, OrdinalType *info) const
Computes an LU factorization of a n by n tridiagonal matrix A using partial pivoting with row interch...
void TRTRS(const char &UPLO, const char &TRANS, const char &DIAG, const OrdinalType &n, const OrdinalType &nrhs, const ScalarType *A, const OrdinalType &lda, ScalarType *B, const OrdinalType &ldb, OrdinalType *info) const
Solves a triangular linear system of the form A*X=B or A**T*X=B, where A is a triangular matrix...
void POTRS(const char &UPLO, const OrdinalType &n, const OrdinalType &nrhs, const ScalarType *A, const OrdinalType &lda, ScalarType *B, const OrdinalType &ldb, OrdinalType *info) const
Solves a system of linear equations A*X=B, where A is a symmetric positive definite matrix factored b...
void POSV(const char &UPLO, const OrdinalType &n, const OrdinalType &nrhs, ScalarType *A, const OrdinalType &lda, ScalarType *B, const OrdinalType &ldb, OrdinalType *info) const
Computes the solution to a real system of linear equations A*X=B, where A is a symmetric positive def...
Teuchos header file which uses auto-configuration information to include necessary C++ headers...
void PTTRF(const OrdinalType &n, MagnitudeType *d, ScalarType *e, OrdinalType *info) const
Computes the L*D*L' factorization of a Hermitian/symmetric positive definite tridiagonal matrix A...
void GBTRF(const OrdinalType &m, const OrdinalType &n, const OrdinalType &kl, const OrdinalType &ku, ScalarType *A, const OrdinalType &lda, OrdinalType *IPIV, OrdinalType *info) const
Computes an LU factorization of a general banded m by n matrix A using partial pivoting with row inte...
void UNM2R(const char &SIDE, const char &TRANS, const OrdinalType &M, const OrdinalType &N, const OrdinalType &K, const ScalarType *A, const OrdinalType &LDA, const ScalarType *TAU, ScalarType *C, const OrdinalType &LDC, ScalarType *WORK, OrdinalType *const INFO) const
BLAS 2 version of UNMQR; known workspace size.
LAPACK(void)
Default Constructor.
void PTEQR(const char &COMPZ, const OrdinalType &n, MagnitudeType *D, MagnitudeType *E, ScalarType *Z, const OrdinalType &ldz, MagnitudeType *WORK, OrdinalType *info) const
Computes the eigenvalues and, optionally, eigenvectors of a symmetric positive-definite tridiagonal n...
void TREVC(const char &SIDE, const char &HOWMNY, OrdinalType *select, const OrdinalType &n, const ScalarType *T, const OrdinalType &ldt, ScalarType *VL, const OrdinalType &ldvl, ScalarType *VR, const OrdinalType &ldvr, const OrdinalType &mm, OrdinalType *m, ScalarType *WORK, OrdinalType *info) const
void STEQR(const char &COMPZ, const OrdinalType &n, MagnitudeType *D, MagnitudeType *E, ScalarType *Z, const OrdinalType &ldz, MagnitudeType *WORK, OrdinalType *info) const
Computes the eigenvalues and, optionally, eigenvectors of a symmetric tridiagonal n by n matrix A usi...
static magnitudeType sfmin()
Returns safe minimum (sfmin), such that 1/sfmin does not overflow.
This structure defines some basic traits for a scalar field type.
void GETRF(const OrdinalType &m, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, OrdinalType *IPIV, OrdinalType *info) const
Computes an LU factorization of a general m by n matrix A using partial pivoting with row interchange...
ScalarTraits< ScalarType >::magnitudeType LANGB(const char &NORM, const OrdinalType &n, const OrdinalType &kl, const OrdinalType &ku, const ScalarType *A, const OrdinalType &lda, MagnitudeType *WORK) const
Returns the value of the one norm, or the Frobenius norm, or the infinity norm, or the element of lar...
void GGES(const char &JOBVL, const char &JOBVR, const char &SORT, OrdinalType &(*ptr2func)(ScalarType *, ScalarType *, ScalarType *), const OrdinalType &n, ScalarType *A, const OrdinalType &lda, ScalarType *B, const OrdinalType &ldb, OrdinalType *sdim, MagnitudeType *ALPHAR, MagnitudeType *ALPHAI, MagnitudeType *BETA, ScalarType *VL, const OrdinalType &ldvl, ScalarType *VR, const OrdinalType &ldvr, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *BWORK, OrdinalType *info) const
void ORMHR(const char &SIDE, const char &TRANS, const OrdinalType &m, const OrdinalType &n, const OrdinalType &ilo, const OrdinalType &ihi, const ScalarType *A, const OrdinalType &lda, const ScalarType *TAU, ScalarType *C, const OrdinalType &ldc, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *info) const
Overwrites the general real m by n matrix C with the product of C and Q, which is a product of ihi-il...
LAPACK(const LAPACK< OrdinalType, ScalarType > &lapack)
Copy Constructor.
void GEES(const char &JOBVS, const char &SORT, OrdinalType &(*ptr2func)(ScalarType *, ScalarType *), const OrdinalType &n, ScalarType *A, const OrdinalType &lda, OrdinalType *sdim, ScalarType *WR, ScalarType *WI, ScalarType *VS, const OrdinalType &ldvs, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *BWORK, OrdinalType *info) const
void GEBAL(const char &JOBZ, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, OrdinalType *ilo, OrdinalType *ihi, MagnitudeType *scale, OrdinalType *info) const
Balances a general matrix A, through similarity transformations to make the rows and columns as close...
void PORFS(const char &UPLO, const OrdinalType &n, const OrdinalType &nrhs, const ScalarType *A, const OrdinalType &lda, const ScalarType *AF, const OrdinalType &ldaf, const ScalarType *B, const OrdinalType &ldb, ScalarType *X, const OrdinalType &ldx, ScalarType *FERR, ScalarType *BERR, ScalarType *WORK, OrdinalType *IWORK, OrdinalType *info) const
Improves the computed solution to a system of linear equations when the coefficient matrix is symmetr...
void LATRS(const char &UPLO, const char &TRANS, const char &DIAG, const char &NORMIN, const OrdinalType &N, const ScalarType *A, const OrdinalType &LDA, ScalarType *X, MagnitudeType *SCALE, MagnitudeType *CNORM, OrdinalType *INFO) const
Robustly solve a possibly singular triangular linear system.
ScalarType LARND(const OrdinalType &idist, OrdinalType *seed) const
Returns a random number from a uniform or normal distribution.
void POTRF(const char &UPLO, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, OrdinalType *info) const
Computes Cholesky factorization of a real symmetric positive definite matrix A.
void SYEV(const char &JOBZ, const char &UPLO, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, ScalarType *W, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *info) const
Computes all the eigenvalues and, optionally, eigenvectors of a symmetric n by n matrix A...
void TREXC(const char &COMPQ, const OrdinalType &n, ScalarType *T, const OrdinalType &ldt, ScalarType *Q, const OrdinalType &ldq, OrdinalType *ifst, OrdinalType *ilst, ScalarType *WORK, OrdinalType *info) const
The Templated LAPACK Wrapper Class.
void GGEV(const char &JOBVL, const char &JOBVR, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, ScalarType *B, const OrdinalType &ldb, MagnitudeType *ALPHAR, MagnitudeType *ALPHAI, ScalarType *BETA, ScalarType *VL, const OrdinalType &ldvl, ScalarType *VR, const OrdinalType &ldvr, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *info) const
void SYGV(const OrdinalType &itype, const char &JOBZ, const char &UPLO, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, ScalarType *B, const OrdinalType &ldb, ScalarType *W, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *info) const
Computes all the eigenvalues and, optionally, eigenvectors of a symmetric n by n matrix pencil {A...
void LASCL(const char &TYPE, const OrdinalType &kl, const OrdinalType &ku, const MagnitudeType cfrom, const MagnitudeType cto, const OrdinalType &m, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, OrdinalType *info) const
Multiplies the m by n matrix A by the real scalar cto/cfrom.
void POSVX(const char &FACT, const char &UPLO, const OrdinalType &n, const OrdinalType &nrhs, ScalarType *A, const OrdinalType &lda, ScalarType *AF, const OrdinalType &ldaf, char *EQUED, ScalarType *S, ScalarType *B, const OrdinalType &ldb, ScalarType *X, const OrdinalType &ldx, ScalarType *rcond, ScalarType *FERR, ScalarType *BERR, ScalarType *WORK, OrdinalType *IWORK, OrdinalType *info) const
Uses the Cholesky factorization to compute the solution to a real system of linear equations A*X=B...
ScalarType LAPY2(const ScalarType &x, const ScalarType &y) const
Computes x^2 + y^2 safely, to avoid overflow.
void TRSEN(const char &JOB, const char &COMPQ, const OrdinalType *SELECT, const OrdinalType &n, ScalarType *T, const OrdinalType &ldt, ScalarType *Q, const OrdinalType &ldq, MagnitudeType *WR, MagnitudeType *WI, OrdinalType *M, ScalarType *S, MagnitudeType *SEP, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *IWORK, const OrdinalType &liwork, OrdinalType *info) const
void GEHRD(const OrdinalType &n, const OrdinalType &ilo, const OrdinalType &ihi, ScalarType *A, const OrdinalType &lda, ScalarType *TAU, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *info) const
Reduces a real general matrix A to upper Hessenberg form by orthogonal similarity transformations...
static magnitudeType magnitude(T a)
Returns the magnitudeType of the scalar type a.
void UNGQR(const OrdinalType &m, const OrdinalType &n, const OrdinalType &k, ScalarType *A, const OrdinalType &lda, const ScalarType *TAU, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *info) const
Compute explicit QR factor from QR factorization (GEQRF) (complex case).
void GBTRS(const char &TRANS, const OrdinalType &n, const OrdinalType &kl, const OrdinalType &ku, const OrdinalType &nrhs, const ScalarType *A, const OrdinalType &lda, const OrdinalType *IPIV, ScalarType *B, const OrdinalType &ldb, OrdinalType *info) const
Solves a system of linear equations A*X=B or A'*X=B with a general banded n by n matrix A using the L...
void GEQRF(const OrdinalType &m, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, ScalarType *TAU, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *info) const
Computes a QR factorization of a general m by n matrix A.
void GESVD(const char &JOBU, const char &JOBVT, const OrdinalType &m, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, MagnitudeType *S, ScalarType *U, const OrdinalType &ldu, ScalarType *V, const OrdinalType &ldv, ScalarType *WORK, const OrdinalType &lwork, MagnitudeType *RWORK, OrdinalType *info) const
Computes the singular values (and optionally, vectors) of a real matrix A.
void GESV(const OrdinalType &n, const OrdinalType &nrhs, ScalarType *A, const OrdinalType &lda, OrdinalType *IPIV, ScalarType *B, const OrdinalType &ldb, OrdinalType *info) const
Computes the solution to a real system of linear equations A*X=B, where A is factored through GETRF a...
virtual ~LAPACK(void)
Destructor.
void GEEV(const char &JOBVL, const char &JOBVR, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, MagnitudeType *WR, MagnitudeType *WI, ScalarType *VL, const OrdinalType &ldvl, ScalarType *VR, const OrdinalType &ldvr, ScalarType *WORK, const OrdinalType &lwork, MagnitudeType *RWORK, OrdinalType *info) const
Computes for an n by n real nonsymmetric matrix A, the eigenvalues and, optionally, the left and/or right eigenvectors.
void GETRS(const char &TRANS, const OrdinalType &n, const OrdinalType &nrhs, const ScalarType *A, const OrdinalType &lda, const OrdinalType *IPIV, ScalarType *B, const OrdinalType &ldb, OrdinalType *info) const
Solves a system of linear equations A*X=B or A'*X=B with a general n by n matrix A using the LU facto...
void GEEVX(const char &BALANC, const char &JOBVL, const char &JOBVR, const char &SENSE, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, ScalarType *WR, ScalarType *WI, ScalarType *VL, const OrdinalType &ldvl, ScalarType *VR, const OrdinalType &ldvr, OrdinalType *ilo, OrdinalType *ihi, MagnitudeType *SCALE, MagnitudeType *abnrm, MagnitudeType *RCONDE, MagnitudeType *RCONDV, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *IWORK, OrdinalType *info) const
void GESVX(const char &FACT, const char &TRANS, const OrdinalType &n, const OrdinalType &nrhs, ScalarType *A, const OrdinalType &lda, ScalarType *AF, const OrdinalType &ldaf, OrdinalType *IPIV, char *EQUED, ScalarType *R, ScalarType *C, ScalarType *B, const OrdinalType &ldb, ScalarType *X, const OrdinalType &ldx, ScalarType *rcond, ScalarType *FERR, ScalarType *BERR, ScalarType *WORK, OrdinalType *IWORK, OrdinalType *info) const
Uses the LU factorization to compute the solution to a real system of linear equations A*X=B...
Defines basic traits for the scalar field type.
static T zero()
Returns representation of zero for this scalar type.
void GETRI(const OrdinalType &n, ScalarType *A, const OrdinalType &lda, const OrdinalType *IPIV, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *info) const
Computes the inverse of a matrix A using the LU factorization computed by GETRF.
void GELS(const char &TRANS, const OrdinalType &m, const OrdinalType &n, const OrdinalType &nrhs, ScalarType *A, const OrdinalType &lda, ScalarType *B, const OrdinalType &ldb, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *info) const
Solves an over/underdetermined real m by n linear system A using QR or LQ factorization of A...
void GEQR2(const OrdinalType &m, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, ScalarType *TAU, ScalarType *WORK, OrdinalType *const info) const
BLAS 2 version of GEQRF, with known workspace size.
void ORMQR(const char &SIDE, const char &TRANS, const OrdinalType &m, const OrdinalType &n, const OrdinalType &k, const ScalarType *A, const OrdinalType &lda, const ScalarType *TAU, ScalarType *C, const OrdinalType &ldc, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *info) const
void GELSS(const OrdinalType &m, const OrdinalType &n, const OrdinalType &nrhs, ScalarType *A, const OrdinalType &lda, ScalarType *B, const OrdinalType &ldb, MagnitudeType *S, const MagnitudeType rcond, OrdinalType *rank, ScalarType *WORK, const OrdinalType &lwork, MagnitudeType *RWORK, OrdinalType *info) const
Use the SVD to solve a possibly rank-deficient linear least-squares problem.
void HEEV(const char &JOBZ, const char &UPLO, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, MagnitudeType *W, ScalarType *WORK, const OrdinalType &lwork, MagnitudeType *RWORK, OrdinalType *info) const
Computes all the eigenvalues and, optionally, eigenvectors of a Hermitian n by n matrix A...
void LARNV(const OrdinalType &idist, OrdinalType *seed, const OrdinalType &n, ScalarType *v) const
Returns a vector of random numbers from a chosen distribution.
void TGSEN(const OrdinalType &ijob, const OrdinalType &wantq, const OrdinalType &wantz, const OrdinalType *SELECT, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, ScalarType *B, const OrdinalType &ldb, MagnitudeType *ALPHAR, MagnitudeType *ALPHAI, MagnitudeType *BETA, ScalarType *Q, const OrdinalType &ldq, ScalarType *Z, const OrdinalType &ldz, OrdinalType *M, MagnitudeType *PL, MagnitudeType *PR, MagnitudeType *DIF, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *IWORK, const OrdinalType &liwork, OrdinalType *info) const
void POCON(const char &UPLO, const OrdinalType &n, const ScalarType *A, const OrdinalType &lda, const ScalarType &anorm, ScalarType *rcond, ScalarType *WORK, OrdinalType *IWORK, OrdinalType *info) const
Estimates the reciprocal of the condition number (1-norm) of a real symmetric positive definite matri...
void PTTRS(const OrdinalType &n, const OrdinalType &nrhs, const MagnitudeType *d, const ScalarType *e, ScalarType *B, const OrdinalType &ldb, OrdinalType *info) const
Solves a tridiagonal system A*X=B using the *D*L' factorization of A computed by PTTRF.
void GERFS(const char &TRANS, const OrdinalType &n, const OrdinalType &nrhs, const ScalarType *A, const OrdinalType &lda, const ScalarType *AF, const OrdinalType &ldaf, const OrdinalType *IPIV, const ScalarType *B, const OrdinalType &ldb, ScalarType *X, const OrdinalType &ldx, ScalarType *FERR, ScalarType *BERR, ScalarType *WORK, OrdinalType *IWORK, OrdinalType *info) const
Improves the computed solution to a system of linear equations and provides error bounds and backward...
void GEQP3(const OrdinalType &m, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, OrdinalType *jpvt, ScalarType *TAU, ScalarType *WORK, const OrdinalType &lwork, MagnitudeType *RWORK, OrdinalType *info) const
Computes a QR factorization with column pivoting of a matrix A: A*P = Q*R using Level 3 BLAS...
ScalarType LAMCH(const char &CMACH) const
Determines machine parameters for floating point characteristics.
void HEGV(const OrdinalType &itype, const char &JOBZ, const char &UPLO, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, ScalarType *B, const OrdinalType &ldb, MagnitudeType *W, ScalarType *WORK, const OrdinalType &lwork, MagnitudeType *RWORK, OrdinalType *info) const
Computes all the eigenvalues and, optionally, eigenvectors of a generalized Hermitian-definite n by n...
void ORGHR(const OrdinalType &n, const OrdinalType &ilo, const OrdinalType &ihi, ScalarType *A, const OrdinalType &lda, const ScalarType *TAU, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *info) const
Generates a real orthogonal matrix Q which is the product of ihi-ilo elementary reflectors of order n...
void TRTRI(const char &UPLO, const char &DIAG, const OrdinalType &n, ScalarType *A, const OrdinalType &lda, OrdinalType *info) const
Computes the inverse of an upper or lower triangular matrix A.
void GECON(const char &NORM, const OrdinalType &n, const ScalarType *A, const OrdinalType &lda, const ScalarType &anorm, ScalarType *rcond, ScalarType *WORK, OrdinalType *IWORK, OrdinalType *info) const
Estimates the reciprocal of the condition number of a general real matrix A, in either the 1-norm or ...
static T one()
Returns representation of one for this scalar type.
void HSEQR(const char &JOB, const char &COMPZ, const OrdinalType &n, const OrdinalType &ilo, const OrdinalType &ihi, ScalarType *H, const OrdinalType &ldh, ScalarType *WR, ScalarType *WI, ScalarType *Z, const OrdinalType &ldz, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *info) const
Computes the eigenvalues of a real upper Hessenberg matrix H and, optionally, the matrices T and Z fr...
void LASWP(const OrdinalType &N, ScalarType *A, const OrdinalType &LDA, const OrdinalType &K1, const OrdinalType &K2, const OrdinalType *IPIV, const OrdinalType &INCX) const
Apply a series of row interchanges to the matrix A.
void GBRFS(const char &TRANS, const OrdinalType &n, const OrdinalType &kl, const OrdinalType &ku, const OrdinalType &nrhs, const ScalarType *A, const OrdinalType &lda, const ScalarType *AF, const OrdinalType &ldaf, const OrdinalType *IPIV, const ScalarType *B, const OrdinalType &ldb, ScalarType *X, const OrdinalType &ldx, ScalarType *FERR, ScalarType *BERR, ScalarType *WORK, OrdinalType *IWORK, OrdinalType *info) const
Improves the computed solution to a banded system of linear equations and provides error bounds and b...
void POEQU(const OrdinalType &n, const ScalarType *A, const OrdinalType &lda, MagnitudeType *S, MagnitudeType *scond, MagnitudeType *amax, OrdinalType *info) const
Computes row and column scalings intended to equilibrate a symmetric positive definite matrix A and r...
void UNMQR(const char &SIDE, const char &TRANS, const OrdinalType &m, const OrdinalType &n, const OrdinalType &k, const ScalarType *A, const OrdinalType &lda, const ScalarType *TAU, ScalarType *C, const OrdinalType &ldc, ScalarType *WORK, const OrdinalType &lwork, OrdinalType *info) const
Apply Householder reflectors (complex case).
void GBEQU(const OrdinalType &m, const OrdinalType &n, const OrdinalType &kl, const OrdinalType &ku, const ScalarType *A, const OrdinalType &lda, MagnitudeType *R, MagnitudeType *C, MagnitudeType *rowcond, MagnitudeType *colcond, MagnitudeType *amax, OrdinalType *info) const
Computes row and column scalings intended to equilibrate an m by n banded matrix A and reduce its con...