Teuchos Package Browser (Single Doxygen Collection) Version of the Day
Teuchos_LAPACK.hpp
Go to the documentation of this file.
00001 // @HEADER
00002 // ***********************************************************************
00003 //
00004 //                    Teuchos: Common Tools Package
00005 //                 Copyright (2004) Sandia Corporation
00006 //
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 //
00010 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
00038 //
00039 // ***********************************************************************
00040 // @HEADER
00041 
00042 #ifndef _TEUCHOS_LAPACK_HPP_
00043 #define _TEUCHOS_LAPACK_HPP_
00044 
00052 #include "Teuchos_ConfigDefs.hpp"
00053 #include "Teuchos_ScalarTraits.hpp"
00054 
00085 namespace Teuchos
00086 {
00087 
00088   template<class T>
00089   struct UndefinedLAPACKRoutine
00090   {
00091     // This function should not compile if there is an attempt to instantiate!
00092     static inline T notDefined() { return T::LAPACK_routine_not_defined_for_this_type(); }
00093   };
00094 
00095   template<typename OrdinalType, typename ScalarType>
00096   class LAPACK
00097   {
00098   public:
00099 
00100     typedef typename Teuchos::ScalarTraits<ScalarType>::magnitudeType MagnitudeType;
00101 
00103 
00104 
00106     inline LAPACK(void) {}
00107 
00109     inline LAPACK(const LAPACK<OrdinalType, ScalarType>& lapack) {}
00110 
00112     inline virtual ~LAPACK(void) {}
00114 
00116 
00117 
00119     void PTTRF(const OrdinalType n, ScalarType* d, ScalarType* e, OrdinalType* info) const;
00120 
00122     void PTTRS(const OrdinalType n, const OrdinalType nrhs, const ScalarType* d, const ScalarType* e, ScalarType* B, const OrdinalType ldb, OrdinalType* info) const;
00123 
00125     void POTRF(const char UPLO, const OrdinalType n, ScalarType* A, const OrdinalType lda, OrdinalType* info) const;
00126 
00128     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;
00129 
00131     void POTRI(const char UPLO, const OrdinalType n, ScalarType* A, const OrdinalType lda, OrdinalType* info) const;
00132 
00134 
00135     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;
00136 
00138     void POSV(const char UPLO, const OrdinalType n, const OrdinalType nrhs, ScalarType* A, const OrdinalType lda, ScalarType* B, const OrdinalType ldb, OrdinalType* info) const;
00139 
00141     void POEQU(const OrdinalType n, const ScalarType* A, const OrdinalType lda, ScalarType* S, ScalarType* scond, ScalarType* amax, OrdinalType* info) const;
00142 
00144     void PORFS(const char UPLO, const OrdinalType n, const OrdinalType nrhs, 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;
00145 
00147     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;
00149 
00151 
00152 
00154     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;
00155 
00189     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;
00190 
00192     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;
00193 
00195     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;
00196 
00198     void GEQRF( const OrdinalType m, const OrdinalType n, ScalarType* A, const OrdinalType lda, ScalarType* TAU, ScalarType* WORK, const OrdinalType lwork, OrdinalType* info) const;
00199 
00201     void GETRF(const OrdinalType m, const OrdinalType n, ScalarType* A, const OrdinalType lda, OrdinalType* IPIV, OrdinalType* info) const;
00202 
00204     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;
00205 
00207     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;
00208 
00210     void
00211     LASWP (const OrdinalType N,
00212      ScalarType A[],
00213      const OrdinalType LDA,
00214      const OrdinalType K1,
00215      const OrdinalType K2,
00216      const OrdinalType IPIV[],
00217      const OrdinalType INCX) const;
00218 
00220     void GBTRF(const OrdinalType m, const OrdinalType n, const OrdinalType kl, const OrdinalType ku, ScalarType* A, const OrdinalType lda, OrdinalType* IPIV, OrdinalType* info) const;
00221 
00223     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;
00224 
00226     void GTTRF(const OrdinalType n, ScalarType* dl, ScalarType* d, ScalarType* du, ScalarType* du2, OrdinalType* IPIV, OrdinalType* info) const;
00227 
00229     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;
00230 
00232     void GETRI(const OrdinalType n, ScalarType* A, const OrdinalType lda, const OrdinalType* IPIV, ScalarType* WORK, const OrdinalType lwork, OrdinalType* info) const;
00233 
00238     void
00239     LATRS (const char UPLO,
00240      const char TRANS,
00241      const char DIAG,
00242      const char NORMIN,
00243      const OrdinalType N,
00244      ScalarType* A,
00245      const OrdinalType LDA,
00246      ScalarType* X,
00247      MagnitudeType* SCALE,
00248      MagnitudeType* CNORM,
00249      OrdinalType* INFO) const;
00250 
00252     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;
00253 
00255     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;
00256 
00258     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;
00259 
00261     void GESV(const OrdinalType n, const OrdinalType nrhs, ScalarType* A, const OrdinalType lda, OrdinalType* IPIV, ScalarType* B, const OrdinalType ldb, OrdinalType* info) const;
00262 
00264     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;
00265 
00267     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;
00268 
00270     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;
00271 
00273     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;
00274 
00276     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;
00277 
00281     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;
00282 
00284     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;
00285 
00287     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;
00288 
00290     void TRTRI(const char UPLO, const char DIAG, const OrdinalType n, const ScalarType* A, const OrdinalType lda, OrdinalType* info) const;
00292 
00294 
00295 
00298     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;
00299 
00303     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;
00304 
00308     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;
00309 
00313     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;
00314 
00318     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;
00319 
00321     void STEQR(const char COMPZ, const OrdinalType n, ScalarType* D, ScalarType* E, ScalarType* Z, const OrdinalType ldz, ScalarType* WORK, OrdinalType* info) const;
00323 
00325 
00326 
00327     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;
00328 
00332     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;
00333 
00337     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;
00338 
00342     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;
00343 
00345     void GEEV(const char JOBVL, const char JOBVR, const OrdinalType n, ScalarType* A, const OrdinalType lda, ScalarType* WR, ScalarType* WI, ScalarType* VL, const OrdinalType ldvl, ScalarType* VR, const OrdinalType ldvr, ScalarType* WORK, const OrdinalType lwork, OrdinalType* info) const;
00346 
00351     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;
00352 
00357     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;
00358 
00362     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;
00363 
00365 
00366 
00368 
00369 
00370     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;
00372 
00373 
00375 
00376 
00386     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;
00387 
00396     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;
00397 
00407     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;
00408 
00417     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;
00418 
00422     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;
00423 
00427     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;
00429 
00431 
00432 
00435     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;
00436 
00440     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;
00441 
00445     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;
00446 
00448 
00450 
00451 
00453     void LARTG( const ScalarType f, const ScalarType g, MagnitudeType* c, ScalarType* s, ScalarType* r ) const;
00454 
00456     void LARFG( const OrdinalType n, ScalarType* alpha, ScalarType* x, const OrdinalType incx, ScalarType* tau ) const;
00457 
00459 
00461 
00462 
00464     void GEBAL(const char JOBZ, const OrdinalType n, ScalarType* A, const OrdinalType lda, OrdinalType ilo, OrdinalType ihi, MagnitudeType* scale, OrdinalType* info) const;
00465 
00467     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;
00468 
00470 
00472 
00473 
00474     ScalarType LARND( const OrdinalType idist, OrdinalType* seed ) const;
00475 
00477     void LARNV( const OrdinalType idist, OrdinalType* seed, const OrdinalType n, ScalarType* v ) const;
00479 
00481 
00482 
00485     ScalarType LAMCH(const char CMACH) const;
00486 
00491     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;
00493 
00495 
00496 
00499     ScalarType LAPY2(const ScalarType x, const ScalarType y) const;
00501   };
00502 
00503   // END GENERAL TEMPLATE DECLARATION //
00504 
00505   // BEGIN GENERAL TEMPLATE IMPLEMENTATION //
00506 
00507 
00508   template<typename OrdinalType, typename ScalarType>
00509   void LAPACK<OrdinalType, ScalarType>::PTTRF(const OrdinalType n, ScalarType* d, ScalarType* e, OrdinalType* info) const
00510   {
00511     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00512   }
00513 
00514   template<typename OrdinalType, typename ScalarType>
00515   void LAPACK<OrdinalType, ScalarType>::PTTRS(const OrdinalType n, const OrdinalType nrhs, const ScalarType* d, const ScalarType* e, ScalarType* B, const OrdinalType ldb, OrdinalType* info) const
00516   {
00517     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00518   }
00519 
00520   template<typename OrdinalType, typename ScalarType>
00521   void LAPACK<OrdinalType, ScalarType>::POTRF(const char UPLO, const OrdinalType n, ScalarType* A, const OrdinalType lda, OrdinalType* info) const
00522   {
00523     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00524   }
00525 
00526   template<typename OrdinalType, typename ScalarType>
00527   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
00528   {
00529     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00530   }
00531 
00532   template<typename OrdinalType, typename ScalarType>
00533   void LAPACK<OrdinalType, ScalarType>::POTRI(const char UPLO, const OrdinalType n, ScalarType* A, const OrdinalType lda, OrdinalType* info) const
00534   {
00535     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00536   }
00537 
00538   template<typename OrdinalType, typename ScalarType>
00539   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
00540   {
00541     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00542   }
00543 
00544   template<typename OrdinalType, typename ScalarType>
00545   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
00546   {
00547     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00548   }
00549 
00550   template<typename OrdinalType, typename ScalarType>
00551   void LAPACK<OrdinalType, ScalarType>::POEQU(const OrdinalType n, const ScalarType* A, const OrdinalType lda, ScalarType* S, ScalarType* scond, ScalarType* amax, OrdinalType* info) const
00552   {
00553     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00554   }
00555 
00556   template<typename OrdinalType, typename ScalarType>
00557   void LAPACK<OrdinalType, ScalarType>::PORFS(const char UPLO, const OrdinalType n, const OrdinalType nrhs, 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
00558   {
00559     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00560   }
00561 
00562   template<typename OrdinalType, typename ScalarType>
00563   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
00564   {
00565     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00566   }
00567 
00568   template<typename OrdinalType, typename ScalarType>
00569   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
00570   {
00571     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00572   }
00573 
00574   template<typename OrdinalType, typename ScalarType>
00575   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
00576   {
00577     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00578   }
00579 
00580   template<typename OrdinalType, typename ScalarType>
00581   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
00582   {
00583     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00584   }
00585 
00586   template<typename OrdinalType, typename ScalarType>
00587   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
00588   {
00589     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00590   }
00591 
00592   template<typename OrdinalType, typename ScalarType>
00593   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
00594   {
00595     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00596   }
00597 
00598   template<typename OrdinalType, typename ScalarType>
00599   void LAPACK<OrdinalType,ScalarType>::GETRF(const OrdinalType m, const OrdinalType n, ScalarType* A, const OrdinalType lda, OrdinalType* IPIV, OrdinalType* info) const
00600   {
00601     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00602   }
00603 
00604   template<typename OrdinalType, typename ScalarType>
00605   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
00606   {
00607     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00608   }
00609 
00610   template<typename OrdinalType, typename ScalarType>
00611   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
00612   {
00613     MagnitudeType safeMin = ScalarTraits<ScalarType>::sfmin();
00614     ScalarType sZero = ScalarTraits<ScalarType>::zero();
00615     ScalarType sOne  = ScalarTraits<ScalarType>::one();
00616     MagnitudeType mZero = ScalarTraits<ScalarType>::magnitude(sZero);
00617     MagnitudeType mOne = ScalarTraits<ScalarType>::magnitude(sOne);
00618 
00619     MagnitudeType smlnum = ScalarTraits<ScalarType>::magnitude(safeMin);
00620     MagnitudeType bignum = ScalarTraits<ScalarType>::magnitude(sOne/smlnum);
00621 
00622     OrdinalType i, j;
00623     ScalarType* ptr;
00624     MagnitudeType mul;
00625     bool done = false;
00626 
00627     MagnitudeType cfromc = cfrom;
00628     MagnitudeType ctoc = cto;
00629     MagnitudeType cfrom1;
00630     MagnitudeType cto1;
00631 
00632     while (!done) {
00633 
00634       cfrom1 = cfromc*smlnum;
00635       if (cfrom1 == cfromc) {
00636   // cfromc is an inf. Multiply by a correctly signed zero for finite ctoc, or a NaN if ctoc is infinite.
00637   mul = ctoc / cfromc;
00638   done = true;
00639   cto1 = ctoc;
00640     } else {
00641   cto1 = ctoc / bignum;
00642   if (cto1 == ctoc) {
00643     // ctoc is either 0 or an inf. In both cases, ctoc itself serves as the correct multiplication factor.
00644     mul = ctoc;
00645     done = true;
00646     cfromc = mOne;
00647   } else if (ScalarTraits<ScalarType>::magnitude(cfrom1) > ScalarTraits<ScalarType>::magnitude(ctoc) && ctoc != mZero) {
00648     mul = smlnum;
00649     done = false;
00650     cfromc = cfrom1;
00651   } else if (ScalarTraits<ScalarType>::magnitude(cto1) > ScalarTraits<ScalarType>::magnitude(cfromc)) {
00652     mul = bignum;
00653     done = false;
00654     ctoc = cto1;
00655   } else {
00656     mul = ctoc / cfromc;
00657     done = true;
00658   }
00659       }
00660 
00661       for (j=0; j<n; j++) {
00662   ptr = A + j*lda;
00663   for (i=0; i<m; i++) { *ptr = mul * (*ptr); ptr++; }
00664       }
00665     }
00666 
00667   }
00668 
00669   template<typename OrdinalType, typename ScalarType>
00670   void
00671   LAPACK<OrdinalType, ScalarType>::
00672   LASWP (const OrdinalType N,
00673    ScalarType A[],
00674    const OrdinalType LDA,
00675    const OrdinalType K1,
00676    const OrdinalType K2,
00677    const OrdinalType IPIV[],
00678    const OrdinalType INCX) const
00679   {
00680     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00681   }
00682 
00683   template<typename OrdinalType, typename ScalarType>
00684   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
00685   {
00686     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00687   }
00688 
00689   template<typename OrdinalType, typename ScalarType>
00690   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
00691   {
00692     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00693   }
00694 
00695   template<typename OrdinalType, typename ScalarType>
00696   void LAPACK<OrdinalType,ScalarType>::GTTRF(const OrdinalType n, ScalarType* dl, ScalarType* d, ScalarType* du, ScalarType* du2, OrdinalType* IPIV, OrdinalType* info) const
00697   {
00698     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00699   }
00700 
00701   template<typename OrdinalType, typename ScalarType>
00702   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
00703   {
00704     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00705   }
00706 
00707   template<typename OrdinalType, typename ScalarType>
00708   void LAPACK<OrdinalType,ScalarType>::GETRI(const OrdinalType n, ScalarType* A, const OrdinalType lda, const OrdinalType* IPIV, ScalarType* WORK, const OrdinalType lwork, OrdinalType* info) const
00709   {
00710     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00711   }
00712 
00713   template<typename OrdinalType, typename ScalarType>
00714   void
00715   LAPACK<OrdinalType,ScalarType>::
00716   LATRS (const char UPLO,
00717    const char TRANS,
00718    const char DIAG,
00719    const char NORMIN,
00720    const OrdinalType N,
00721    ScalarType* A,
00722    const OrdinalType LDA,
00723    ScalarType* X,
00724    MagnitudeType* SCALE,
00725    MagnitudeType* CNORM,
00726    OrdinalType* INFO) const
00727   {
00728     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00729   }
00730 
00731   template<typename OrdinalType, typename ScalarType>
00732   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
00733   {
00734     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00735   }
00736 
00737   template<typename OrdinalType, typename ScalarType>
00738   void LAPACK<OrdinalType,ScalarType>::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
00739   {
00740     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00741   }
00742 
00743   template<typename OrdinalType, typename ScalarType>
00744   typename ScalarTraits<ScalarType>::magnitudeType LAPACK<OrdinalType,ScalarType>::LANGB(const char NORM, const OrdinalType n, const OrdinalType kl, const OrdinalType ku, const ScalarType* A, const OrdinalType lda, MagnitudeType* WORK) const
00745   {
00746     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00747   }
00748 
00749   template<typename OrdinalType, typename ScalarType>
00750   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
00751   {
00752     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00753   }
00754 
00755   template<typename OrdinalType, typename ScalarType>
00756   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
00757   {
00758 
00759     // Test the input parameters
00760     info = 0;
00761     if (m < 0) {
00762       info = -1;
00763     } else if (n < 0) {
00764       info = -2;
00765     } else if (lda < TEUCHOS_MAX(1, m)) {
00766       info = -4;
00767     }
00768     if (info != 0) {
00769       return;
00770     }
00771 
00772     ScalarType sZero = ScalarTraits<ScalarType>::zero();
00773     ScalarType sOne  = ScalarTraits<ScalarType>::one();
00774     MagnitudeType mZero = ScalarTraits<ScalarType>::magnitude(sZero);
00775     MagnitudeType mOne = ScalarTraits<ScalarType>::magnitude(sOne);
00776 
00777     // Quick return
00778     if (m == 0 || n == 0) {
00779       *rowcond = mOne;
00780       *colcond = mOne;
00781       *amax = mZero;
00782       return;
00783     }
00784 
00785     MagnitudeType safeMin = ScalarTraits<ScalarType>::sfmin();
00786     MagnitudeType smlnum = ScalarTraits<ScalarType>::magnitude(safeMin);
00787     MagnitudeType bignum = ScalarTraits<ScalarType>::magnitude(sOne/smlnum);
00788 
00789     // Compute the row scale factors
00790     for (OrdinalType i=0; i<m; i++) {
00791       R[i] = mZero;
00792     }
00793 
00794     // Find the maximum element in each row
00795     for (OrdinalType j=0; j<n; j++) {
00796       for (OrdinalType i=0; i<m; i++) {
00797   R[i] = TEUCHOS_MAX( R[i], ScalarTraits<ScalarType>::magnitude( A[j*lda + i] ) );
00798       }
00799     }
00800 
00801     // Find the maximum and minimum scale factors
00802     MagnitudeType rcmin = bignum;
00803     MagnitudeType rcmax = mZero;
00804     for (OrdinalType i=0; i<m; i++) {
00805       rcmax = TEUCHOS_MAX( rcmax, R[i] );
00806       rcmin = TEUCHOS_MIN( rcmin, R[i] );
00807     }
00808     *amax = rcmax;
00809 
00810     if (rcmin == mZero) {
00811       // Find the first zero scale factor and return an error code
00812       for (OrdinalType i=0; i<m; i++) {
00813   if (R[i] == mZero)
00814     *info = i;
00815       }
00816     } else {
00817       // Invert the scale factors
00818       for (OrdinalType i=0; i<m; i++) {
00819   R[i] = mOne / TEUCHOS_MIN( TEUCHOS_MAX( R[i], smlnum ), bignum );
00820       }
00821       // Compute rowcond = min(R(i)) / max(R(i))
00822       *rowcond = TEUCHOS_MAX( rcmin, smlnum ) / TEUCHOS_MIN( rcmax, bignum );
00823     }
00824 
00825     // Compute the column scale factors
00826     for (OrdinalType j=0; j<n; j++) {
00827       C[j] = mZero;
00828     }
00829 
00830     // Find the maximum element in each column, assuming the row scaling computed above
00831     for (OrdinalType j=0; j<n; j++) {
00832       for (OrdinalType i=0; i<m; i++) {
00833   C[j] = TEUCHOS_MAX( C[j], R[i]*ScalarTraits<ScalarType>::magnitude( A[j*lda + i] ) );
00834       }
00835     }
00836 
00837     // Find the maximum and minimum scale factors
00838     rcmin = bignum;
00839     rcmax = mZero;
00840     for (OrdinalType j=0; j<n; j++) {
00841       rcmax = TEUCHOS_MAX( rcmax, C[j] );
00842       rcmin = TEUCHOS_MIN( rcmin, C[j] );
00843     }
00844 
00845     if (rcmin == mZero) {
00846       // Find the first zero scale factor and return an error code
00847       for (OrdinalType j=0; j<n; j++) {
00848   if (C[j] == mZero)
00849     *info = m+j;
00850       }
00851     } else {
00852       // Invert the scale factors
00853       for (OrdinalType j=0; j<n; j++) {
00854   C[j] = mOne / TEUCHOS_MIN( TEUCHOS_MAX( C[j], smlnum ), bignum );
00855       }
00856       // Compute colcond = min(C(j)) / max(C(j))
00857       *colcond = TEUCHOS_MAX( rcmin, smlnum ) / TEUCHOS_MIN( rcmax, bignum );
00858     }
00859   }
00860 
00861   template<typename OrdinalType, typename ScalarType>
00862   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
00863   {
00864     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00865   }
00866 
00867   template<typename OrdinalType, typename ScalarType>
00868   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
00869   {
00870 
00871     // Test the input parameters
00872     info = 0;
00873     if (m < 0) {
00874       info = -1;
00875     } else if (n < 0) {
00876       info = -2;
00877     } else if (kl < 0) {
00878       info = -3;
00879     } else if (ku < 0) {
00880       info = -4;
00881     } else if (lda < kl+ku+1) {
00882       info = -6;
00883     }
00884     if (info != 0) {
00885       return;
00886     }
00887 
00888     ScalarType sZero = ScalarTraits<ScalarType>::zero();
00889     ScalarType sOne  = ScalarTraits<ScalarType>::one();
00890     MagnitudeType mZero = ScalarTraits<ScalarType>::magnitude(sZero);
00891     MagnitudeType mOne = ScalarTraits<ScalarType>::magnitude(sOne);
00892 
00893     // Quick return
00894     if (m == 0 || n == 0) {
00895       *rowcond = mOne;
00896       *colcond = mOne;
00897       *amax = mZero;
00898       return;
00899     }
00900 
00901     MagnitudeType safeMin = ScalarTraits<ScalarType>::sfmin();
00902     MagnitudeType smlnum = ScalarTraits<ScalarType>::magnitude(safeMin);
00903     MagnitudeType bignum = ScalarTraits<ScalarType>::magnitude(sOne/smlnum);
00904 
00905     // Compute the row scale factors
00906     for (OrdinalType i=0; i<m; i++) {
00907       R[i] = mZero;
00908     }
00909 
00910     // Find the maximum element in each row
00911     for (OrdinalType j=0; j<n; j++) {
00912       for (OrdinalType i=TEUCHOS_MAX(j-ku,0); i<TEUCHOS_MIN(j+kl,m-1); i++) {
00913   R[i] = TEUCHOS_MAX( R[i], ScalarTraits<ScalarType>::magnitude( A[j*lda + ku+i-j] ) );
00914       }
00915     }
00916 
00917     // Find the maximum and minimum scale factors
00918     MagnitudeType rcmin = bignum;
00919     MagnitudeType rcmax = mZero;
00920     for (OrdinalType i=0; i<m; i++) {
00921       rcmax = TEUCHOS_MAX( rcmax, R[i] );
00922       rcmin = TEUCHOS_MIN( rcmin, R[i] );
00923     }
00924     *amax = rcmax;
00925 
00926     if (rcmin == mZero) {
00927       // Find the first zero scale factor and return an error code
00928       for (OrdinalType i=0; i<m; i++) {
00929   if (R[i] == mZero)
00930     *info = i;
00931       }
00932     } else {
00933       // Invert the scale factors
00934       for (OrdinalType i=0; i<m; i++) {
00935   R[i] = mOne / TEUCHOS_MIN( TEUCHOS_MAX( R[i], smlnum ), bignum );
00936       }
00937       // Compute rowcond = min(R(i)) / max(R(i))
00938       *rowcond = TEUCHOS_MAX( rcmin, smlnum ) / TEUCHOS_MIN( rcmax, bignum );
00939     }
00940 
00941     // Compute the column scale factors
00942     for (OrdinalType j=0; j<n; j++) {
00943       C[j] = mZero;
00944     }
00945 
00946     // Find the maximum element in each column, assuming the row scaling computed above
00947     for (OrdinalType j=0; j<n; j++) {
00948       for (OrdinalType i=TEUCHOS_MAX(j-ku,0); i<TEUCHOS_MIN(j+kl,m-1); i++) {
00949   C[j] = TEUCHOS_MAX( C[j], R[i]*ScalarTraits<ScalarType>::magnitude( A[j*lda + ku+i-j] ) );
00950       }
00951     }
00952 
00953     // Find the maximum and minimum scale factors
00954     rcmin = bignum;
00955     rcmax = mZero;
00956     for (OrdinalType j=0; j<n; j++) {
00957       rcmax = TEUCHOS_MAX( rcmax, C[j] );
00958       rcmin = TEUCHOS_MIN( rcmin, C[j] );
00959     }
00960 
00961     if (rcmin == mZero) {
00962       // Find the first zero scale factor and return an error code
00963       for (OrdinalType j=0; j<n; j++) {
00964   if (C[j] == mZero)
00965     *info = m+j;
00966       }
00967     } else {
00968       // Invert the scale factors
00969       for (OrdinalType j=0; j<n; j++) {
00970   C[j] = mOne / TEUCHOS_MIN( TEUCHOS_MAX( C[j], smlnum ), bignum );
00971       }
00972       // Compute colcond = min(C(j)) / max(C(j))
00973       *colcond = TEUCHOS_MAX( rcmin, smlnum ) / TEUCHOS_MIN( rcmax, bignum );
00974     }
00975   }
00976 
00977   template<typename OrdinalType, typename ScalarType>
00978   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
00979   {
00980     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00981   }
00982 
00983   template<typename OrdinalType, typename ScalarType>
00984   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
00985   {
00986     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00987   }
00988 
00989   template<typename OrdinalType, typename ScalarType>
00990   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
00991   {
00992     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00993   }
00994 
00995   template<typename OrdinalType, typename ScalarType>
00996   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
00997   {
00998     UndefinedLAPACKRoutine<ScalarType>::notDefined();
00999   }
01000 
01001   template<typename OrdinalType, typename ScalarType>
01002   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
01003   {
01004     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01005   }
01006 
01007   template<typename OrdinalType, typename ScalarType>
01008   void LAPACK<OrdinalType,ScalarType>::TRTRI(const char UPLO, const char DIAG, const OrdinalType n, const ScalarType* A, const OrdinalType lda, OrdinalType* info) const
01009   {
01010     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01011   }
01012 
01013   template<typename OrdinalType, typename ScalarType>
01014   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
01015   {
01016     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01017   }
01018 
01019   template<typename OrdinalType, typename ScalarType>
01020   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
01021   {
01022     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01023   }
01024 
01025   template<typename OrdinalType, typename ScalarType>
01026   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
01027   {
01028     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01029   }
01030 
01031   template<typename OrdinalType, typename ScalarType>
01032   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
01033   {
01034     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01035   }
01036 
01037   template<typename OrdinalType, typename ScalarType>
01038   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
01039   {
01040     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01041   }
01042 
01043   template<typename OrdinalType, typename ScalarType>
01044   void LAPACK<OrdinalType,ScalarType>::STEQR(const char COMPZ, const OrdinalType n, ScalarType* D, ScalarType* E, ScalarType* Z, const OrdinalType ldz, ScalarType* WORK, OrdinalType* info) const
01045   {
01046     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01047   }
01048 
01049   template<typename OrdinalType, typename ScalarType>
01050   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
01051   {
01052     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01053   }
01054 
01055   template<typename OrdinalType, typename ScalarType>
01056   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
01057   {
01058     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01059   }
01060 
01061   template<typename OrdinalType, typename ScalarType>
01062   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
01063   {
01064     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01065   }
01066 
01067   template<typename OrdinalType, typename ScalarType>
01068   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
01069   {
01070     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01071   }
01072 
01073   template<typename OrdinalType, typename ScalarType>
01074   void LAPACK<OrdinalType, ScalarType>::GEEV(const char JOBVL, const char JOBVR, const OrdinalType n, ScalarType* A, const OrdinalType lda, ScalarType* WR, ScalarType* WI, ScalarType* VL, const OrdinalType ldvl, ScalarType* VR, const OrdinalType ldvr, ScalarType* WORK, const OrdinalType lwork, OrdinalType* info) const
01075   {
01076     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01077   }
01078 
01079   template<typename OrdinalType, typename ScalarType>
01080   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
01081   {
01082     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01083   }
01084 
01085   template<typename OrdinalType, typename ScalarType>
01086   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
01087   {
01088     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01089   }
01090 
01091   template<typename OrdinalType, typename ScalarType>
01092   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
01093   {
01094     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01095   }
01096 
01097   template<typename OrdinalType, typename ScalarType>
01098   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
01099   {
01100     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01101   }
01102 
01103   template<typename OrdinalType, typename ScalarType>
01104   void LAPACK<OrdinalType, ScalarType>::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
01105   {
01106     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01107   }
01108 
01109   template<typename OrdinalType, typename ScalarType>
01110   void LAPACK<OrdinalType, ScalarType>::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
01111   {
01112     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01113   }
01114 
01115   template<typename OrdinalType, typename ScalarType>
01116   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
01117   {
01118     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01119   }
01120 
01121   template<typename OrdinalType, typename ScalarType>
01122   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
01123   {
01124     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01125   }
01126 
01127   template<typename OrdinalType, typename ScalarType>
01128   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
01129   {
01130     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01131   }
01132 
01133   template<typename OrdinalType, typename ScalarType>
01134   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
01135   {
01136     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01137   }
01138 
01139   template<typename OrdinalType, typename ScalarType>
01140   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
01141   {
01142     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01143   }
01144 
01145   template<typename OrdinalType, typename ScalarType>
01146   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
01147   {
01148     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01149   }
01150 
01151   template<typename OrdinalType, typename ScalarType>
01152   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
01153   {
01154     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01155   }
01156 
01157   template<typename OrdinalType, typename ScalarType>
01158   ScalarType LAPACK<OrdinalType, ScalarType>::LAMCH(const char CMACH) const
01159   {
01160     return UndefinedLAPACKRoutine<ScalarType>::notDefined();
01161   }
01162 
01163   template<typename OrdinalType, typename ScalarType>
01164   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
01165   {
01166     return UndefinedLAPACKRoutine<OrdinalType>::notDefined();
01167   }
01168 
01169   template<typename OrdinalType, typename ScalarType>
01170   ScalarType LAPACK<OrdinalType, ScalarType>::LAPY2(const ScalarType x, const ScalarType y) const
01171   {
01172     return UndefinedLAPACKRoutine<ScalarType>::notDefined();
01173   }
01174 
01175   template<typename OrdinalType, typename ScalarType>
01176   void LAPACK<OrdinalType, ScalarType>::LARTG( const ScalarType f, const ScalarType g, MagnitudeType* c, ScalarType* s, ScalarType* r ) const
01177   {
01178     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01179   }
01180 
01181   template<typename OrdinalType, typename ScalarType>
01182   void LAPACK<OrdinalType, ScalarType>::LARFG( const OrdinalType n, ScalarType* alpha, ScalarType* x, const OrdinalType incx, ScalarType* tau ) const
01183   {
01184     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01185   }
01186 
01187   template<typename OrdinalType, typename ScalarType>
01188   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
01189   {
01190     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01191   }
01192 
01193 
01194   template<typename OrdinalType, typename ScalarType>
01195   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
01196   {
01197     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01198   }
01199 
01200   template<typename OrdinalType, typename ScalarType>
01201   ScalarType LAPACK<OrdinalType, ScalarType>::LARND( const OrdinalType idist, OrdinalType* seed ) const
01202   {
01203     return UndefinedLAPACKRoutine<ScalarType>::notDefined();
01204   }
01205 
01206   template<typename OrdinalType, typename ScalarType>
01207   void LAPACK<OrdinalType, ScalarType>::LARNV( const OrdinalType idist, OrdinalType* seed, const OrdinalType n, ScalarType* v ) const
01208   {
01209     UndefinedLAPACKRoutine<ScalarType>::notDefined();
01210   }
01211 
01212   // END GENERAL TEMPLATE IMPLEMENTATION //
01213 
01214 #ifndef DOXYGEN_SHOULD_SKIP_THIS
01215 
01216   // BEGIN INT, FLOAT SPECIALIZATION DECLARATION //
01217 
01218   template<>
01219   class TEUCHOSNUMERICS_LIB_DLL_EXPORT LAPACK<int, float>
01220   {
01221   public:
01222     inline LAPACK(void) {}
01223     inline LAPACK(const LAPACK<int, float>& lapack) {}
01224     inline virtual ~LAPACK(void) {}
01225 
01226     // Symmetric positive definite linear system routines
01227     void POTRF(const char UPLO, const int n, float* A, const int lda, int * info) const;
01228     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;
01229     void PTTRF(const int n, float* d, float* e, int* info) const;
01230     void PTTRS(const int n, const int nrhs, const float* d, const float* e, float* B, const int ldb, int* info) const;
01231     void POTRI(const char UPLO, const int n, float* A, const int lda, int* info) const;
01232     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;
01233     void POSV(const char UPLO, const int n, const int nrhs, float* A, const int lda, float* B, const int ldb, int* info) const;
01234     void POEQU(const int n, const float* A, const int lda, float* S, float* scond, float* amax, int* info) const;
01235     void PORFS(const char UPLO, const int n, const int nrhs, 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;
01236     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;
01237 
01238     // General Linear System Routines
01239     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;
01240     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;
01241     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;
01242     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;
01243     void GEQRF( const int m, const int n, float* A, const int lda, float* TAU, float* WORK, const int lwork, int* info) const;
01244     void GETRF(const int m, const int n, float* A, const int lda, int* IPIV, int* info) const;
01245     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;
01246     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;
01247 
01248 
01249     void LASWP (const int N,
01250     float A[],
01251     const int LDA,
01252     const int K1,
01253     const int K2,
01254     const int IPIV[],
01255     const int INCX) const;
01256 
01257     void GBTRF(const int m, const int n, const int kl, const int ku, float* A, const int lda, int* IPIV, int* info) const;
01258     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;
01259     void GTTRF(const int n, float* dl, float* d, float* du, float* du2, int* IPIV, int* info) const;
01260     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;
01261 
01262 
01263     void GETRI(const int n, float* A, const int lda, const int* IPIV, float* WORK, const int lwork, int* info) const;
01264     void LATRS (const char UPLO, const char TRANS, const char DIAG, const char NORMIN, const int N, float* A, const int LDA, float* X, float* SCALE, float* CNORM, int* INFO) const;
01265     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;
01266     void GBCON(const char NORM, const int n, const int kl, const int ku, const float* A, const int lda, int* IPIV, const float anorm, float* rcond, float* WORK, int* IWORK, int* info) const;
01267     float LANGB(const char NORM, const int n, const int kl, const int ku, const float* A, const int lda, float* WORK) const;
01268     void GESV(const int n, const int nrhs, float* A, const int lda, int* IPIV, float* B, const int ldb, int* info) const;
01269     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;
01270     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;
01271     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;
01272     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;
01273     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;
01274     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;
01275     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;
01276     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;
01277     void TRTRI(const char UPLO, const char DIAG, const int n, const float* A, const int lda, int* info) const;
01278 
01279     // Symmetric eigenvalue routines.
01280     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;
01281     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;
01282     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;
01283     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;
01284     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;
01285     void STEQR(const char COMPZ, const int n, float* D, float* E, float* Z, const int ldz, float* WORK, int* info) const;
01286 
01287     // Non-Hermitian eigenvalue routines.
01288     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;
01289     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;
01290     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;
01291     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;
01292     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;
01293     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;
01294     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;
01295 
01296     // SVD routine
01297     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;
01298 
01299     // Orthogonal matrix routines.
01300     void ORMQR(const char SIDE, const char TRANS, const int m, const int n, const int k, float* A, const int lda, const float* TAU, float* C, const int ldc, float* WORK, const int lwork, int* info) const;
01301     void UNMQR(const char SIDE, const char TRANS, const int m, const int n, const int k, float* A, const int lda, const float* TAU, float* C, const int ldc, float* WORK, const int lwork, int* info) const;
01302 
01303     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;
01304     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;
01305     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;
01306     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;
01307 
01308     // Triangular matrix routines.
01309     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;
01310     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;
01311     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;
01312 
01313     // Rotation/reflection generators
01314     void LARTG( const float f, const float g, float* c, float* s, float* r ) const;
01315     void LARFG( const int n, float* alpha, float* x, const int incx, float* tau ) const;
01316 
01317     // Matrix balancing routines.
01318     void GEBAL(const char JOBZ, const int n, float* A, const int lda, int ilo, int ihi, float* scale, int* info) const;
01319     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;
01320 
01321     // Random number generators
01322     float LARND( const int idist, int* seed ) const;
01323     void LARNV( const int idist, int* seed, const int n, float* v ) const;
01324 
01325     // Machine characteristics.
01326     float LAMCH(const char CMACH) const;
01327     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;
01328 
01329     // Miscellaneous routines.
01330     float LAPY2(const float x, const float y) const;
01331 
01332   };
01333 
01334   // END INT, FLOAT SPECIALIZATION DECLARATION //
01335 
01336   // BEGIN INT, DOUBLE SPECIALIZATION DECLARATION //
01337 
01338   template<>
01339   class TEUCHOSNUMERICS_LIB_DLL_EXPORT LAPACK<int, double>
01340   {
01341   public:
01342     inline LAPACK(void) {}
01343     inline LAPACK(const LAPACK<int, double>& lapack) {}
01344     inline virtual ~LAPACK(void) {}
01345 
01346     // Symmetric positive definite linear system routines
01347     void PTTRF(const int n, double* d, double* e, int* info) const;
01348     void PTTRS(const int n, const int nrhs, const double* d, const double* e, double* B, const int ldb, int* info) const;
01349     void POTRF(const char UPLO, const int n, double* A, const int lda, int* info) const;
01350     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;
01351     void POTRI(const char UPLO, const int n, double* A, const int lda, int* info) const;
01352     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;
01353     void POSV(const char UPLO, const int n, const int nrhs, double* A, const int lda, double* B, const int ldb, int* info) const;
01354     void POEQU(const int n, const double* A, const int lda, double* S, double* scond, double* amax, int* info) const;
01355     void PORFS(const char UPLO, const int n, const int nrhs, 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;
01356     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;
01357 
01358     // General linear system routines
01359     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;
01360     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;
01361     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;
01362     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;
01363     void GEQRF( const int m, const int n, double* A, const int lda, double* TAU, double* WORK, const int lwork, int* info) const;
01364     void GETRF(const int m, const int n, double* A, const int lda, int* IPIV, int* info) const;
01365     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;
01366     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;
01367 
01368     void LASWP (const int N,
01369     double A[],
01370     const int LDA,
01371     const int K1,
01372     const int K2,
01373     const int IPIV[],
01374     const int INCX) const;
01375 
01376     void GBTRF(const int m, const int n, const int kl, const int ku, double* A, const int lda, int* IPIV, int* info) const;
01377     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;
01378     void GTTRF(const int n, double* dl, double* d, double* du, double* du2, int* IPIV, int* info) const;
01379     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;
01380     void GETRI(const int n, double* A, const int lda, const int* IPIV, double* WORK, const int lwork, int* info) const;
01381     void LATRS (const char UPLO, const char TRANS, const char DIAG, const char NORMIN, const int N, double* A, const int LDA, double* X, double* SCALE, double* CNORM, int* INFO) const;
01382     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;
01383     void GBCON(const char NORM, const int n, const int kl, const int ku, const double* A, const int lda, int* IPIV, const double anorm, double* rcond, double* WORK, int* IWORK, int* info) const;
01384     double LANGB(const char NORM, const int n, const int kl, const int ku, const double* A, const int lda, double* WORK) const;
01385     void GESV(const int n, const int nrhs, double* A, const int lda, int* IPIV, double* B, const int ldb, int* info) const;
01386     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;
01387     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;
01388     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;
01389     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;
01390     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;
01391     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;
01392     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;
01393     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;
01394     void TRTRI(const char UPLO, const char DIAG, const int n, const double* A, const int lda, int* info) const;
01395 
01396     // Symmetric eigenproblem routines.
01397     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;
01398     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;
01399     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;
01400     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;
01401     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;
01402     void STEQR(const char COMPZ, const int n, double* D, double* E, double* Z, const int ldz, double* WORK, int* info) const;
01403 
01404     // Non-Hermitian eigenproblem routines.
01405     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;
01406     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;
01407     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;
01408     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;
01409     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;
01410     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;
01411     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;
01412 
01413 
01414     // SVD routine
01415     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;
01416 
01417     // Orthogonal matrix routines.
01418     void ORMQR(const char SIDE, const char TRANS, const int m, const int n, const int k, double* A, const int lda, const double* TAU, double* C, const int ldc, double* WORK, const int lwork, int* info) const;
01419     void UNMQR(const char SIDE, const char TRANS, const int m, const int n, const int k, double* A, const int lda, const double* TAU, double* C, const int ldc, double* WORK, const int lwork, int* info) const;
01420     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;
01421     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;
01422     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;
01423     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;
01424 
01425     // Triangular matrix routines.
01426     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;
01427     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;
01428     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;
01429 
01430     // Rotation/reflection generators
01431     void LARTG( const double f, const double g, double* c, double* s, double* r ) const;
01432     void LARFG( const int n, double* alpha, double* x, const int incx, double* tau ) const;
01433 
01434     // Matrix balancing routines.
01435     void GEBAL(const char JOBZ, const int n, double* A, const int lda, int ilo, int ihi, double* scale, int* info) const;
01436     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;
01437 
01438     // Random number generators
01439     double LARND( const int idist, int* seed ) const;
01440     void LARNV( const int idist, int* seed, const int n, double* v ) const;
01441 
01442     // Machine characteristic routines.
01443     double LAMCH(const char CMACH) const;
01444     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;
01445 
01446     // Miscellaneous routines.
01447     double LAPY2(const double x, const double y) const;
01448 
01449   };
01450 
01451   // END INT, DOUBLE SPECIALIZATION DECLARATION //
01452 
01453 #ifdef HAVE_TEUCHOS_COMPLEX
01454 
01455   // BEGIN INT, COMPLEX<FLOAT> SPECIALIZATION DECLARATION //
01456 
01457   template<>
01458   class TEUCHOSNUMERICS_LIB_DLL_EXPORT LAPACK<int, std::complex<float> >
01459   {
01460   public:
01461     inline LAPACK(void) {}
01462     inline LAPACK(const LAPACK<int, std::complex<float> >& lapack) {}
01463     inline virtual ~LAPACK(void) {}
01464 
01465     // Symmetric positive definite linear system routines
01466     void PTTRF(const int n, std::complex<float>* d, std::complex<float>* e, int* info) const;
01467     void PTTRS(const int n, const int nrhs, const std::complex<float>* d, const std::complex<float>* e, std::complex<float>* B, const int ldb, int* info) const;
01468     void POTRF(const char UPLO, const int n, std::complex<float>* A, const int lda, int* info) const;
01469     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;
01470     void POTRI(const char UPLO, const int n, std::complex<float>* A, const int lda, int* info) const;
01471     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;
01472     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;
01473     void POEQU(const int n, const std::complex<float>* A, const int lda, float* S, float* scond, float* amax, int* info) const;
01474     void PORFS(const char UPLO, const int n, const int nrhs, 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;
01475     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;
01476 
01477     // General Linear System Routines
01478     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;
01479     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;
01480     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;
01481     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;
01482     void UNMQR(const char SIDE, const char TRANS, const int m, const int n, const int k, 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;
01483     void GETRF(const int m, const int n, std::complex<float>* A, const int lda, int* IPIV, int* info) const;
01484     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;
01485     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;
01486 
01487     void LASWP (const int N,
01488     std::complex<float> A[],
01489     const int LDA,
01490     const int K1,
01491     const int K2,
01492     const int IPIV[],
01493     const int INCX) const;
01494 
01495     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;
01496     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;
01497     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;
01498     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;
01499     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;
01500     void LATRS (const char UPLO, const char TRANS, const char DIAG, const char NORMIN, const int N, std::complex<float>* A, const int LDA, std::complex<float>* X, float* SCALE, float* CNORM, int* INFO) const;
01501     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;
01502     void GBCON(const char NORM, const int n, const int kl, const int ku, const std::complex<float>* A, const int lda, int* IPIV, const float anorm, float* rcond, std::complex<float>* WORK, float* RWORK, int* info) const;
01503     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;
01504     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;
01505     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;
01506     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;
01507     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;
01508     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;
01509     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;
01510     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;
01511     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;
01512     void TRTRI(const char UPLO, const char DIAG, const int n, const std::complex<float>* A, const int lda, int* info) const;
01513 
01514     // Symmetric eigenvalue routines.
01515     void STEQR(const char COMPZ, const int n, float* D, float* E, std::complex<float>* Z, const int ldz, float* WORK, int* info) const;
01516     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;
01517     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;
01518 
01519     // Non-Hermitian eigenvalue routines.
01520     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;
01521     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;
01522     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;
01523     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;
01524     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;
01525     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;
01526     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;
01527 
01528     // SVD routine
01529     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;
01530 
01531     // Triangular matrix routines.
01532     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;
01533     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;
01534     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;
01535 
01536     // Rotation/reflection generators
01537     void LARTG( const std::complex<float> f, const std::complex<float> g, float* c, std::complex<float>* s, std::complex<float>* r ) const;
01538     void LARFG( const int n, std::complex<float>* alpha, std::complex<float>* x, const int incx, std::complex<float>* tau ) const;
01539 
01540     // Matrix balancing routines.
01541     void GEBAL(const char JOBZ, const int n, std::complex<float>* A, const int lda, int ilo, int ihi, float* scale, int* info) const;
01542     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;
01543 
01544     // Random number generators
01545     std::complex<float> LARND( const int idist, int* seed ) const;
01546     void LARNV( const int idist, int* seed, const int n, std::complex<float>* v ) const;
01547 
01548     // Machine characteristics
01549     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;
01550 
01551   };
01552 
01553   // END INT, COMPLEX<FLOAT> SPECIALIZATION DECLARATION //
01554 
01555   // BEGIN INT, COMPLEX<DOUBLE> SPECIALIZATION DECLARATION //
01556 
01557   template<>
01558   class TEUCHOSNUMERICS_LIB_DLL_EXPORT LAPACK<int, std::complex<double> >
01559   {
01560   public:
01561     inline LAPACK(void) {}
01562     inline LAPACK(const LAPACK<int, std::complex<double> >& lapack) {}
01563     inline virtual ~LAPACK(void) {}
01564 
01565     // Symmetric positive definite linear system routines
01566     void PTTRF(const int n, std::complex<double>* d, std::complex<double>* e, int* info) const;
01567     void PTTRS(const int n, const int nrhs, const std::complex<double>* d, const std::complex<double>* e, std::complex<double>* B, const int ldb, int* info) const;
01568     void POTRF(const char UPLO, const int n, std::complex<double>* A, const int lda, int* info) const;
01569     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;
01570     void POTRI(const char UPLO, const int n, std::complex<double>* A, const int lda, int* info) const;
01571     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;
01572     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;
01573     void POEQU(const int n, const std::complex<double>* A, const int lda, double* S, double* scond, double* amax, int* info) const;
01574     void PORFS(const char UPLO, const int n, const int nrhs, 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;
01575     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;
01576 
01577     // General Linear System Routines
01578     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;
01579     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;
01580     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;
01581     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;
01582     void UNMQR(const char SIDE, const char TRANS, const int m, const int n, const int k, 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;
01583     void GETRF(const int m, const int n, std::complex<double>* A, const int lda, int* IPIV, int* info) const;
01584     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;
01585     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;
01586 
01587     void LASWP (const int N,
01588     std::complex<double> A[],
01589     const int LDA,
01590     const int K1,
01591     const int K2,
01592     const int IPIV[],
01593     const int INCX) const;
01594 
01595     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;
01596     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;
01597     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;
01598     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;
01599     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;
01600     void LATRS (const char UPLO, const char TRANS, const char DIAG, const char NORMIN, const int N, std::complex<double>* A, const int LDA, std::complex<double>* X, double* SCALE, double* CNORM, int* INFO) const;
01601     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;
01602     void GBCON(const char NORM, const int n, const int kl, const int ku, const std::complex<double>* A, const int lda, int* IPIV, const double anorm, double* rcond, std::complex<double>* WORK, double* RWORK, int* info) const;
01603     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;
01604     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;
01605     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;
01606     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;
01607     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;
01608     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;
01609     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;
01610     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;
01611     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;
01612     void TRTRI(const char UPLO, const char DIAG, const int n, const std::complex<double>* A, const int lda, int* info) const;
01613 
01614     // Symmetric eigenvalue routines.
01615     void STEQR(const char COMPZ, const int n, double* D, double* E, std::complex<double>* Z, const int ldz, double* WORK, int* info) const;
01616     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;
01617     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;
01618 
01619     // Non-hermitian eigenvalue routines.
01620     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;
01621     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;
01622     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;
01623     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;
01624     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;
01625     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;
01626     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;
01627 
01628     // SVD routine
01629     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;
01630 
01631     // Triangular matrix routines.
01632     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;
01633     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;
01634     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;
01635 
01636     // Rotation/reflection generators
01637     void LARTG( const std::complex<double> f, const std::complex<double> g, double* c, std::complex<double>* s, std::complex<double>* r ) const;
01638     void LARFG( const int n, std::complex<double>* alpha, std::complex<double>* x, const int incx, std::complex<double>* tau ) const;
01639 
01640     // Matrix balancing routines.
01641     void GEBAL(const char JOBZ, const int n, std::complex<double>* A, const int lda, int ilo, int ihi, double* scale, int* info) const;
01642     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;
01643 
01644     // Random number generators
01645     std::complex<double> LARND( const int idist, int* seed ) const;
01646     void LARNV( const int idist, int* seed, const int n, std::complex<double>* v ) const;
01647 
01648     // Machine characteristics
01649     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;
01650 
01651   };
01652 
01653   // END INT, COMPLEX<DOUBLE> SPECIALIZATION DECLARATION //
01654 
01655 #endif // HAVE_TEUCHOS_COMPLEX
01656 
01657 #endif // DOXYGEN_SHOULD_SKIP_THIS
01658 
01659 } // namespace Teuchos
01660 
01661 #endif // _TEUCHOS_LAPACK_HPP_
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines