Stokhos Package Browser (Single Doxygen Collection)  Version of the Day
Teuchos_LAPACK_UQ_PCE.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Teuchos: Common Tools Package
5 // Copyright (2004) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ***********************************************************************
40 // @HEADER
41 
42 #ifndef _TEUCHOS_LAPACK_UQ_PCE_HPP_
43 #define _TEUCHOS_LAPACK_UQ_PCE_HPP_
44 
45 #include "Teuchos_LAPACK.hpp"
46 #include "Teuchos_TestForException.hpp"
47 #include "Sacado_UQ_PCE.hpp"
48 
56 namespace Teuchos {
57 
58  template<typename OrdinalType, typename Storage>
59  class LAPACK<OrdinalType, Sacado::UQ::PCE<Storage> >
60  {
61  public:
62 
64  typedef typename Teuchos::ScalarTraits<ScalarType>::magnitudeType MagnitudeType;
65 
67 
68 
70  inline LAPACK(void) {}
71 
73  inline LAPACK(const LAPACK<OrdinalType, ScalarType>& lapack) {}
74 
76  inline virtual ~LAPACK(void) {}
78 
80 
81 
83  void PTTRF(const OrdinalType n, ScalarType* d, ScalarType* e, OrdinalType* info) const
84  { throw_error("PTTRF"); }
85 
87  void PTTRS(const OrdinalType n, const OrdinalType nrhs, const ScalarType* d, const ScalarType* e, ScalarType* B, const OrdinalType ldb, OrdinalType* info) const
88  { throw_error("PTTRS"); }
89 
91  void POTRF(const char UPLO, const OrdinalType n, ScalarType* A, const OrdinalType lda, OrdinalType* info) const
92  { throw_error("POTRF"); }
93 
95  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
96  { throw_error("POTRS"); }
97 
99  void POTRI(const char UPLO, const OrdinalType n, ScalarType* A, const OrdinalType lda, OrdinalType* info) const
100  { throw_error("POTRI"); }
101 
103 
104  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
105  { throw_error("POCON"); }
106 
108  void POSV(const char UPLO, const OrdinalType n, const OrdinalType nrhs, ScalarType* A, const OrdinalType lda, ScalarType* B, const OrdinalType ldb, OrdinalType* info) const
109  { throw_error("POSV"); }
110 
112  void POEQU(const OrdinalType n, const ScalarType* A, const OrdinalType lda, MagnitudeType* S, MagnitudeType* scond, MagnitudeType* amax, OrdinalType* info) const
113  { throw_error("POEQU"); }
114 
116  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
117  { throw_error("PORFS"); }
118 
120  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
121  { throw_error("POSVX"); }
123 
125 
126 
128  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
129  { throw_error("GELS"); }
130 
164  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
165  { throw_error("GELSS"); }
166 
168  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
169  { throw_error("GELSS"); }
170 
172  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
173  { throw_error("GGLSE"); }
174 
176  void GEQRF( const OrdinalType m, const OrdinalType n, ScalarType* A, const OrdinalType lda, ScalarType* TAU, ScalarType* WORK, const OrdinalType lwork, OrdinalType* info) const
177  { throw_error("GEQRF"); }
178 
180  void GETRF(const OrdinalType m, const OrdinalType n, ScalarType* A, const OrdinalType lda, OrdinalType* IPIV, OrdinalType* info) const
181  { throw_error("GETRF"); }
182 
184  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
185  { throw_error("GETRS"); }
186 
188  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
189  { throw_error("LASCL"); }
190 
192  void
193  GEQP3(const OrdinalType m,
194  const OrdinalType n, ScalarType* A,
195  const OrdinalType lda,
196  OrdinalType *jpvt,
197  ScalarType* TAU,
198  ScalarType* WORK,
199  const OrdinalType lwork,
200  MagnitudeType* RWORK,
201  OrdinalType* info ) const
202  { throw_error("GEQP3"); }
203 
205  void
206  LASWP (const OrdinalType N,
207  ScalarType A[],
208  const OrdinalType LDA,
209  const OrdinalType K1,
210  const OrdinalType K2,
211  const OrdinalType IPIV[],
212  const OrdinalType INCX) const
213  { throw_error("LASWP"); }
214 
216  void GBTRF(const OrdinalType m, const OrdinalType n, const OrdinalType kl, const OrdinalType ku, ScalarType* A, const OrdinalType lda, OrdinalType* IPIV, OrdinalType* info) const
217  { throw_error("GBTRF"); }
218 
220  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
221  { throw_error("GBTRS"); }
222 
224  void GTTRF(const OrdinalType n, ScalarType* dl, ScalarType* d, ScalarType* du, ScalarType* du2, OrdinalType* IPIV, OrdinalType* info) const
225  { throw_error("GTTRF"); }
226 
228  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
229  { throw_error("GTTRS"); }
230 
232  void GETRI(const OrdinalType n, ScalarType* A, const OrdinalType lda, const OrdinalType* IPIV, ScalarType* WORK, const OrdinalType lwork, OrdinalType* info) const
233  { throw_error("GETRI"); }
234 
239  void
240  LATRS (const char UPLO,
241  const char TRANS,
242  const char DIAG,
243  const char NORMIN,
244  const OrdinalType N,
245  ScalarType* A,
246  const OrdinalType LDA,
247  ScalarType* X,
248  MagnitudeType* SCALE,
249  MagnitudeType* CNORM,
250  OrdinalType* INFO) const
251  { throw_error("LATRS"); }
252 
254  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
255  { throw_error("GECON"); }
256 
258  void GBCON(const char NORM, const OrdinalType n, const OrdinalType kl, const OrdinalType ku, const ScalarType* A, const OrdinalType lda, OrdinalType* IPIV, const ScalarType anorm, ScalarType* rcond, ScalarType* WORK, OrdinalType* IWORK, OrdinalType* info) const
259  { throw_error("GBCON"); }
260 
262  typename 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
263  { throw_error("LANGB"); return MagnitudeType(0); }
264 
266  void GESV(const OrdinalType n, const OrdinalType nrhs, ScalarType* A, const OrdinalType lda, OrdinalType* IPIV, ScalarType* B, const OrdinalType ldb, OrdinalType* info) const
267  { throw_error("GESV"); }
268 
270  void GEEQU(const OrdinalType m, const OrdinalType n, const ScalarType* A, const OrdinalType lda, MagnitudeType* R, MagnitudeType* C, MagnitudeType* rowcond, MagnitudeType* colcond, MagnitudeType* amax, OrdinalType* info) const
271  { throw_error("GEEQU"); }
272 
274  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, MagnitudeType* FERR, MagnitudeType* BERR, ScalarType* WORK, OrdinalType* IWORK, OrdinalType* info) const
275  { throw_error("GERFS"); }
276 
278  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
279  { throw_error("GBEQU"); }
280 
282  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
283  { throw_error("GBRFS"); }
284 
286  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
287  { throw_error("GESVX"); }
288 
292  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
293  { throw_error("SYTRD"); }
294 
296  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
297  { throw_error("GEHRD"); }
298 
300  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
301  { throw_error("TRTRS"); }
302 
304  void TRTRI(const char UPLO, const char DIAG, const OrdinalType n, const ScalarType* A, const OrdinalType lda, OrdinalType* info) const
305  { throw_error("TRTRI"); }
307 
309 
310 
313  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
314  { throw_error("SPEV"); }
315 
319  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
320  { throw_error("SYEV"); }
321 
325  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
326  { throw_error("SYGV"); }
327 
331  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
332  { throw_error("HEEV"); }
333 
337  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
338  { throw_error("HEGV"); }
339 
341  void STEQR(const char COMPZ, const OrdinalType n, ScalarType* D, ScalarType* E, ScalarType* Z, const OrdinalType ldz, ScalarType* WORK, OrdinalType* info) const
342  { throw_error("STEQR"); }
344 
346 
347  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
349  { throw_error("HSEQR"); }
350 
354  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
355  { throw_error("GEES"); }
356 
360  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
361  { throw_error("GEES"); }
362 
366  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
367  { throw_error("GEES"); }
368 
374  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
375  { throw_error("GEEV"); }
376 
381  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
382  { throw_error("GEEVX"); }
383 
388  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
389  { throw_error("GGEVX"); }
390 
394  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
395  { throw_error("GGEV"); }
396 
397 
401  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
402  { throw_error("TRSEN"); }
403 
404 
408  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
409  { throw_error("TGSEN"); }
410 
411 
415  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
416  { throw_error("GGES"); }
417 
419 
420 
422 
423  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
425  { throw_error("GESVD"); }
427 
428 
430 
431 
441  void ORMQR(const char SIDE, const char TRANS, const OrdinalType m, const OrdinalType n, const OrdinalType k, ScalarType* A, const OrdinalType lda, const ScalarType* TAU, ScalarType* C, const OrdinalType ldc, ScalarType* WORK, const OrdinalType lwork, OrdinalType* info) const
442  { throw_error("ORMQR"); }
443 
452  void UNMQR(const char SIDE, const char TRANS, const OrdinalType m, const OrdinalType n, const OrdinalType k, ScalarType* A, const OrdinalType lda, const ScalarType* TAU, ScalarType* C, const OrdinalType ldc, ScalarType* WORK, const OrdinalType lwork, OrdinalType* info) const
453  { throw_error("UNMQR"); }
454 
464  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
465  { throw_error("ORGQR"); }
466 
475  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
476  { throw_error("UNGQR"); }
477 
481  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
482  { throw_error("ORGHR"); }
483 
487  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
488  { throw_error("ORMHR"); }
490 
492 
493 
496  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
497  { throw_error("TREVC"); }
498 
502  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
503  { throw_error("TREVC"); }
504 
508  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
509  { throw_error("TREXC"); }
510 
514  void TGEVC(const char SIDE, const char HOWMNY, const OrdinalType *SELECT, const OrdinalType n, ScalarType *S, const OrdinalType lds, 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
515  { throw_error("TGEVC"); }
516 
517 
519 
521 
522 
524  void LARTG( const ScalarType f, const ScalarType g, MagnitudeType* c, ScalarType* s, ScalarType* r ) const
525  { throw_error("LARTG"); }
526 
528  void LARFG( const OrdinalType n, ScalarType* alpha, ScalarType* x, const OrdinalType incx, ScalarType* tau ) const
529  { throw_error("LARFG"); }
530 
532 
534 
535 
537  void GEBAL(const char JOBZ, const OrdinalType n, ScalarType* A, const OrdinalType lda, OrdinalType ilo, OrdinalType ihi, MagnitudeType* scale, OrdinalType* info) const
538  { throw_error("GEBAL"); }
539 
541  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
542  { throw_error("GEBAK"); }
543 
545 
547 
548  ScalarType LARND( const OrdinalType idist, OrdinalType* seed ) const
550  { throw_error("LARND"); return ScalarType(0); }
551 
553  void LARNV( const OrdinalType idist, OrdinalType* seed, const OrdinalType n, ScalarType* v ) const
554  { throw_error("LARNV"); }
556 
558 
559 
562  ScalarType LAMCH(const char CMACH) const
563  { throw_error("LAMCH"); return ScalarType(0); }
564 
569  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
570  { throw_error("ILAENV"); return OrdinalType(0); }
572 
574 
575 
578  ScalarType LAPY2(const ScalarType x, const ScalarType y) const
579  { throw_error("LAPY2"); return ScalarType(0); }
581 
582  private:
583 
584  void throw_error(const char *func) const {
585  TEUCHOS_TEST_FOR_EXCEPTION(
586  true, std::logic_error,
587  func << ": Not implemented for Sacado::UQ::PCE scalar type!");
588  }
589 
590  };
591 
592 } // namespace Teuchos
593 
594 #endif // _TEUCHOS_LAPACK__UQ_PCE_HPP_
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 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 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
Legacy GELSS interface for real-valued ScalarType.
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...
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, MagnitudeType *FERR, MagnitudeType *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 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 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
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 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 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 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...
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...
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 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 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&#39;*X=B or A^H*X=B with a tridiagonal matrix A using the ...
void STEQR(const char COMPZ, const OrdinalType n, ScalarType *D, ScalarType *E, ScalarType *Z, const OrdinalType ldz, ScalarType *WORK, OrdinalType *info) const
Computes the eigenvalues and, optionally, eigenvectors of a symmetric tridiagonal n by n matrix A usi...
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 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...
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 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...
Teuchos::ScalarTraits< ScalarType >::magnitudeType MagnitudeType
void LARNV(const OrdinalType idist, OrdinalType *seed, const OrdinalType n, ScalarType *v) const
Returns a vector of random numbers from a chosen distribution.
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...
ScalarType LAMCH(const char CMACH) const
Determines machine parameters for floating point characteristics.
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 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 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 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
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 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 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 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
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 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 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...
const IndexType const IndexType const IndexType const IndexType const ValueType const ValueType * x
Definition: csr_vector.h:260
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 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 TGEVC(const char SIDE, const char HOWMNY, const OrdinalType *SELECT, const OrdinalType n, ScalarType *S, const OrdinalType lds, 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 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
void PTTRS(const OrdinalType n, const OrdinalType nrhs, const ScalarType *d, const ScalarType *e, ScalarType *B, const OrdinalType ldb, OrdinalType *info) const
Solves a tridiagonal system A*X=B using the *D*L&#39; factorization of A computed by PTTRF.
void UNMQR(const char SIDE, const char TRANS, const OrdinalType m, const OrdinalType n, const OrdinalType k, 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 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.
Specialization for Sacado::UQ::PCE< Storage<...> >
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 LATRS(const char UPLO, const char TRANS, const char DIAG, const char NORMIN, const OrdinalType N, ScalarType *A, const OrdinalType LDA, ScalarType *X, MagnitudeType *SCALE, MagnitudeType *CNORM, OrdinalType *INFO) const
Robustly solve a possibly singular triangular linear system.
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 TRTRI(const char UPLO, const char DIAG, const OrdinalType n, const ScalarType *A, const OrdinalType lda, OrdinalType *info) const
Computes the inverse of an upper or lower triangular matrix A.
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...
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 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 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 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 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 ...
void GEEQU(const OrdinalType m, const OrdinalType n, 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 matrix A and reduce its condition ...
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
void ORMQR(const char SIDE, const char TRANS, const OrdinalType m, const OrdinalType n, const OrdinalType k, ScalarType *A, const OrdinalType lda, const ScalarType *TAU, ScalarType *C, const OrdinalType ldc, ScalarType *WORK, const OrdinalType lwork, OrdinalType *info) const
ScalarType f(const Teuchos::Array< ScalarType > &x, double a, double b)
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 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&#39;*X=B with a general banded n by n matrix A using the L...
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 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 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...
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 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 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 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 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.
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 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 GBCON(const char NORM, const OrdinalType n, const OrdinalType kl, const OrdinalType ku, const ScalarType *A, const OrdinalType lda, 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 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&#39;*X=B with a general n by n matrix A using the LU facto...
const IndexType const IndexType const IndexType const IndexType const ValueType const ValueType ValueType * y
Definition: csr_vector.h:267
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...
ScalarType g(const Teuchos::Array< ScalarType > &x, const ScalarType &y)
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...
void PTTRF(const OrdinalType n, ScalarType *d, ScalarType *e, OrdinalType *info) const
Computes the L*D*L&#39; factorization of a Hermitian/symmetric positive definite tridiagonal matrix 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).