Teuchos - Trilinos Tools Package Version of the Day
Teuchos_LAPACK.cpp
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 #include "Teuchos_ConfigDefs.hpp"
00043 #include "Teuchos_ScalarTraits.hpp"
00044 #include "Teuchos_LAPACK.hpp"
00045 #include "Teuchos_LAPACK_wrappers.hpp"
00046 
00047 /* for INTEL_CXML, the second arg may need to be changed to 'one'.  If so
00048 the appropriate declaration of one will need to be added back into
00049 functions that include the macro:
00050 */
00051 
00052 #ifdef CHAR_MACRO
00053 #undef CHAR_MACRO
00054 #endif
00055 #if defined (INTEL_CXML)
00056 #define CHAR_MACRO(char_var) &char_var, one
00057 #else
00058 #define CHAR_MACRO(char_var) &char_var
00059 #endif
00060 
00061 namespace {
00062 
00063 #if defined (INTEL_CXML)
00064         unsigned int one=1;
00065 #endif
00066 
00067 // Use a warpper function to handle claling ILAENV().  This removes
00068 // duplicaiton and avoid name lookup problems with member functions called
00069 // ILAENV() trying to call nonmember functions called ILAENV() (which does not
00070 // work on Intel compiler on Windows, see Trilinos bug 5762).
00071 inline
00072 int ilaenv_wrapper(
00073   const int* ispec, const char* name, unsigned int name_length,
00074   const char* opts, unsigned int opts_length,
00075   const int* N1, const int* N2, const int* N3, const int* N4 )
00076 {
00077 #if defined (INTEL_CXML)
00078     return ILAENV_F77(ispec, name, name_length, opts, opts_length, N1, N2, N3, N4 );
00079 #else
00080     return ILAENV_F77(ispec, name, opts, N1, N2, N3, N4, name_length, opts_length );
00081 #endif
00082 }
00083 
00084 } // namespace
00085 
00086 
00087 
00088 extern "C" {
00089 
00090 
00091 typedef int (*gees_nullfptr_t)(double*,double*);
00092 
00093 
00094 } // extern "C"
00095 
00096 
00097 namespace Teuchos
00098 {
00099   // BEGIN INT, FLOAT SPECIALIZATION IMPLEMENTATION //
00100 
00101   void LAPACK<int, float>::PTTRF(const int n, float* d, float* e, int* info) const
00102   { SPTTRF_F77(&n,d,e,info); }
00103 
00104 
00105   void LAPACK<int, float>::PTTRS(const int n, const int nrhs, const float* d, const float* e, float* B, const int ldb, int* info) const
00106   { SPTTRS_F77(&n,&nrhs,d,e,B,&ldb,info); }
00107 
00108 
00109   void LAPACK<int, float>::POTRF(const char UPLO, const int n, float* A, const int lda, int* info) const
00110   { SPOTRF_F77(CHAR_MACRO(UPLO), &n, A, &lda, info); }
00111 
00112 
00113   void LAPACK<int, float>::POTRS(const char UPLO, const int n, const int nrhs, const float* A, const int lda, float* B, const int ldb, int* info) const
00114   { SPOTRS_F77(CHAR_MACRO(UPLO), &n, &nrhs, A, &lda, B, &ldb, info); }
00115 
00116 
00117   void LAPACK<int, float>::POTRI(const char UPLO, const int n, float* A, const int lda, int* info) const
00118   { SPOTRI_F77(CHAR_MACRO(UPLO), &n, A, &lda, info); }
00119 
00120 
00121   void LAPACK<int, float>::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
00122   { SPOCON_F77(CHAR_MACRO(UPLO), &n, A, &lda, &anorm, rcond, WORK, IWORK, info); }
00123 
00124 
00125   void LAPACK<int, float>::POSV(const char UPLO, const int n, const int nrhs, float* A, const int lda, float* B, const int ldb, int* info) const
00126   { SPOSV_F77(CHAR_MACRO(UPLO), &n, &nrhs, A, &lda, B, &ldb, info); }
00127 
00128 
00129   void LAPACK<int, float>::POEQU(const int n, const float* A, const int lda, float* S, float* scond, float* amax, int* info) const
00130   { SPOEQU_F77(&n, A, &lda, S, scond, amax, info); }
00131 
00132 
00133   void LAPACK<int, float>::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
00134   { SPORFS_F77(CHAR_MACRO(UPLO), &n, &nrhs, A, &lda, AF, &ldaf, B, &ldb, X, &ldx, FERR, BERR, WORK, IWORK, info); }
00135 
00136 
00137   void LAPACK<int, float>::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
00138   { SPOSVX_F77(CHAR_MACRO(FACT), CHAR_MACRO(UPLO), &n, &nrhs, A, &lda, AF, &ldaf, CHAR_MACRO(EQUED), S, B, &ldb, X, &ldx, rcond, FERR, BERR, WORK, IWORK, info); }
00139 
00140 
00141   void LAPACK<int,float>::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
00142   { SGELS_F77(CHAR_MACRO(TRANS), &m, &n, &nrhs, A, &lda, B, &ldb, WORK, &lwork, info); }
00143 
00144   void LAPACK<int,float>::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
00145   {
00146     (void) rwork;
00147     SGELSS_F77(&m, &n, &nrhs, A, &lda, B, &ldb, S, &rcond, rank, WORK, &lwork, info);
00148   }
00149 
00150   void LAPACK<int,float>::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
00151   { SGELSS_F77(&m, &n, &nrhs, A, &lda, B, &ldb, S, &rcond, rank, WORK, &lwork, info); }
00152 
00153 
00154   void LAPACK<int,float>::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
00155   { SGGLSE_F77(&m, &n, &p, A, &lda, B, &ldb, C, D, X, WORK, &lwork, info); }
00156 
00157 
00158   void LAPACK<int,float>::GEQRF( const int m, const int n, float* A, const int lda, float* TAU, float* WORK, const int lwork, int* info) const
00159   { SGEQRF_F77(&m, &n, A, &lda, TAU, WORK, &lwork, info); }
00160 
00161 
00162   void LAPACK<int,float>::GETRF(const int m, const int n, float* A, const int lda, int* IPIV, int* info) const
00163   { SGETRF_F77(&m, &n, A, &lda, IPIV, info); }
00164 
00165 
00166   void LAPACK<int,float>::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
00167   { SGETRS_F77(CHAR_MACRO(TRANS), &n, &nrhs, A, &lda, IPIV, B, &ldb, info); }
00168 
00169 
00170   void LAPACK<int,float>::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
00171   { SLASCL_F77(CHAR_MACRO(TYPE), &kl, &ku, &cfrom, &cto, &m, &n, A, &lda, info); }
00172 
00173   void
00174   LAPACK<int,float>::
00175   GEQP3 (const int m,
00176          const int n,
00177          float* A,
00178          const int lda,
00179          int *jpvt,
00180          float* TAU,
00181          float* WORK,
00182          const int lwork,
00183          float* RWORK,
00184          int* info) const
00185   {
00186     (void) RWORK;
00187     SGEQP3_F77(&m, &n, A, &lda, jpvt, TAU, WORK, &lwork, info);
00188   }
00189 
00190   void LAPACK<int, float>::
00191   LASWP (const int N,
00192          float A[],
00193          const int LDA,
00194          const int K1,
00195          const int K2,
00196          const int IPIV[],
00197          const int INCX) const
00198   {
00199     SLASWP_F77(&N, A, &LDA, &K1, &K2, IPIV, &INCX);
00200   }
00201 
00202   void LAPACK<int,float>::GBTRF(const int m, const int n, const int kl, const int ku, float* A, const int lda, int* IPIV, int* info) const
00203   { SGBTRF_F77(&m, &n, &kl, &ku, A, &lda, IPIV, info); }
00204 
00205 
00206   void LAPACK<int,float>::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
00207   { SGBTRS_F77(CHAR_MACRO(TRANS), &n, &kl, &ku, &nrhs, A, &lda, IPIV, B, &ldb, info); }
00208 
00209 
00210   void LAPACK<int,float>::GTTRF(const int n, float* dl, float* d, float* du, float* du2, int* IPIV, int* info) const
00211   { SGTTRF_F77(&n, dl, d, du, du2, IPIV, info); }
00212 
00213 
00214   void LAPACK<int,float>::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
00215   { SGTTRS_F77(CHAR_MACRO(TRANS), &n, &nrhs, dl, d, du, du2, IPIV, B, &ldb, info); }
00216 
00217 
00218   void LAPACK<int,float>::GETRI(const int n, float* A, const int lda, const int* IPIV, float* WORK, const int lwork, int* info) const
00219   { SGETRI_F77(&n, A, &lda, IPIV, WORK, &lwork, info); }
00220 
00221   void
00222   LAPACK<int, float>::LATRS (const char UPLO,
00223                              const char TRANS,
00224                              const char DIAG,
00225                              const char NORMIN,
00226                              const int N,
00227                              float* A,
00228                              const int LDA,
00229                              float* X,
00230                              float* SCALE,
00231                              float* CNORM,
00232                              int* INFO) const
00233   {
00234     SLATRS_F77(CHAR_MACRO(UPLO),
00235                CHAR_MACRO(TRANS),
00236                CHAR_MACRO(DIAG),
00237                CHAR_MACRO(NORMIN),
00238                &N,
00239                A,
00240                &LDA,
00241                X,
00242                SCALE,
00243                CNORM,
00244                INFO);
00245   }
00246 
00247 
00248   void LAPACK<int,float>::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
00249   { SGECON_F77(CHAR_MACRO(NORM), &n, A, &lda, &anorm, rcond, WORK, IWORK, info); }
00250 
00251 
00252   void LAPACK<int,float>::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
00253   { SGBCON_F77(CHAR_MACRO(NORM), &n, &kl, &ku, A, &lda, IPIV, &anorm, rcond, WORK, IWORK, info); }
00254 
00255 
00256   float LAPACK<int,float>::LANGB(const char NORM, const int n, const int kl, const int ku, const float* A, const int lda, float* WORK) const
00257   { return( SLANGB_F77(CHAR_MACRO(NORM), &n, &kl, &ku, A, &lda, WORK) ); }
00258 
00259 
00260   void LAPACK<int,float>::GESV(const int n, const int nrhs, float* A, const int lda, int* IPIV, float* B, const int ldb, int* info) const
00261   { SGESV_F77(&n, &nrhs, A, &lda, IPIV, B, &ldb, info); }
00262 
00263 
00264   void LAPACK<int,float>::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
00265   { SGEEQU_F77(&m, &n, A, &lda, R, C, rowcond, colcond, amax, info); }
00266 
00267 
00268   void LAPACK<int,float>::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
00269   { SGERFS_F77(CHAR_MACRO(TRANS), &n, &nrhs, A, &lda, AF, &ldaf, IPIV, B, &ldb, X, &ldx, FERR, BERR, WORK, IWORK, info); }
00270 
00271 
00272   void LAPACK<int,float>::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
00273   { SGBEQU_F77(&m, &n, &kl, &ku, A, &lda, R, C, rowcond, colcond, amax, info); }
00274 
00275 
00276   void LAPACK<int,float>::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
00277   { SGBRFS_F77(CHAR_MACRO(TRANS), &n, &kl, &ku, &nrhs, A, &lda, AF, &ldaf, IPIV, B, &ldb, X, &ldx, FERR, BERR, WORK, IWORK, info); }
00278 
00279 
00280   void LAPACK<int,float>::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
00281   { SGESVX_F77(CHAR_MACRO(FACT), CHAR_MACRO(TRANS), &n, &nrhs, A, &lda, AF, &ldaf, IPIV, CHAR_MACRO(EQUED), R, C, B, &ldb, X, &ldx, rcond, FERR, BERR, WORK, IWORK, info); }
00282 
00283 
00284   void LAPACK<int,float>::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
00285   { SSYTRD_F77(CHAR_MACRO(UPLO), &n, A, &lda, D, E, TAU, WORK, &lwork, info); }
00286 
00287 
00288   void LAPACK<int,float>::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
00289   { SGEHRD_F77(&n, &ilo, &ihi, A, &lda, TAU, WORK, &lwork, info); }
00290 
00291 
00292   void LAPACK<int,float>::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
00293   { STRTRS_F77(CHAR_MACRO(UPLO), CHAR_MACRO(TRANS), CHAR_MACRO(DIAG), &n, &nrhs, A, &lda, B, &ldb, info); }
00294 
00295 
00296   void LAPACK<int,float>::TRTRI(const char UPLO, const char DIAG, const int n, const float* A, const int lda, int* info) const
00297   { STRTRI_F77(CHAR_MACRO(UPLO), CHAR_MACRO(DIAG), &n, A, &lda, info); }
00298 
00299 
00300   void LAPACK<int,float>::SPEV(const char JOBZ, const char UPLO, const int n, float* AP, float* W, float* Z, const int ldz, float* WORK, int* info) const
00301   { SSPEV_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &n, AP, W, Z, &ldz, WORK, info); }
00302 
00303 
00304   void LAPACK<int,float>::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
00305   { SSYEV_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &n, A, &lda, W, WORK, &lwork, info); }
00306 
00307 
00308   void LAPACK<int,float>::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
00309   { SSYGV_F77(&itype, CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &n, A, &lda, B, &ldb, W, WORK, &lwork, info); }
00310 
00311 
00312   void LAPACK<int,float>::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
00313   { SSYEV_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &n, A, &lda, W, WORK, &lwork, info); }
00314 
00315 
00316   void LAPACK<int,float>::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
00317   { SSYGV_F77(&itype, CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &n, A, &lda, B, &ldb, W, WORK, &lwork, info); }
00318 
00319 
00320   void LAPACK<int,float>::STEQR(const char COMPZ, const int n, float* D, float* E, float* Z, const int ldz, float* WORK, int* info) const
00321   { SSTEQR_F77(CHAR_MACRO(COMPZ), &n, D, E, Z, &ldz, WORK, info); }
00322 
00323 
00324   void LAPACK<int, float>::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
00325   { SHSEQR_F77(CHAR_MACRO(JOB), CHAR_MACRO(COMPZ), &n, &ilo, &ihi, H, &ldh, WR, WI, Z, &ldz, WORK, &lwork, info); }
00326 
00327 
00328   void LAPACK<int, float>::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
00329   { SGEES_F77(CHAR_MACRO(JOBVS), CHAR_MACRO(SORT), ptr2func, &n, A, &lda, sdim, WR, WI, VS, &ldvs, WORK, &lwork, BWORK, info); }
00330 
00331 
00332   void LAPACK<int, float>::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
00333   {
00334     int (*nullfptr)(float*,float*) = NULL;
00335     const char sort = 'N';
00336     SGEES_F77(CHAR_MACRO(JOBVS), CHAR_MACRO(sort), nullfptr, &n, A, &lda, sdim, WR, WI, VS, &ldvs, WORK, &lwork, BWORK, info);
00337   }
00338 
00339 
00340   void LAPACK<int, float>::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
00341   { SGEEV_F77(CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), &n, A, &lda, WR, WI, VL, &ldvl, VR, &ldvr, WORK, &lwork, info); }
00342 
00343   void LAPACK<int, float>::GEEV(const char JOBVL, const char JOBVR, const int n, float* A, const int lda, float* WR, float* WI, float* VL, const int ldvl, float* VR, const int ldvr, float* WORK, const int lwork, float* /* RWORK */, int* info) const
00344   {
00345     GEEV (JOBVL, JOBVR, n, A, lda, WR, WI, VL, ldvl, VR, ldvr, WORK, lwork, info);
00346   }
00347 
00348 
00349   void LAPACK<int, float>::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
00350   { SGESVD_F77(CHAR_MACRO(JOBU), CHAR_MACRO(JOBVT), &m, &n, A, &lda, S, U, &ldu, V, &ldv, WORK, &lwork, info); }
00351 
00352 
00353   void LAPACK<int,float>::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
00354   { SGEEVX_F77(CHAR_MACRO(BALANC), CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), CHAR_MACRO(SENSE), &n, A, &lda, WR, WI, VL, &ldvl, VR, &ldvr, ilo, ihi, SCALE, abnrm, RCONDE, RCONDV, WORK, &lwork, IWORK, info); }
00355 
00356 
00357   void LAPACK<int,float>::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
00358   { SGGEVX_F77(CHAR_MACRO(BALANC), CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), CHAR_MACRO(SENSE), &n, A, &lda, B, &ldb, ALPHAR, ALPHAI, BETA, VL, &ldvl, VR, &ldvr, ilo, ihi, lscale, rscale, abnrm, bbnrm, RCONDE, RCONDV, WORK, &lwork, IWORK, BWORK, info); }
00359 
00360   void LAPACK<int,float>::GGEVX(const char BALANC, const char JOBVL, const char JOBVR, const char SENSE, const int n, float* A, const int lda, float* B, const int ldb, float* ALPHAR, float* ALPHAI, float* BETA, float* VL, const int ldvl, float* VR, const int ldvr, int* ilo, int* ihi, float* lscale, float* rscale, float* abnrm, float* bbnrm, float* RCONDE, float* RCONDV, float* WORK, const int lwork, float* /* RWORK */, int* IWORK, int* BWORK, int* info) const
00361   {
00362     GGEVX(BALANC, JOBVL, JOBVR, SENSE, n, A, lda, B, ldb, ALPHAR, ALPHAI, BETA, VL, ldvl, VR, ldvr, ilo, ihi, lscale, rscale, abnrm, bbnrm, RCONDE, RCONDV, WORK, lwork, IWORK, BWORK, info);
00363   }
00364 
00365   void LAPACK<int, float>::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
00366   { SGGEV_F77(CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), &n, A, &lda, B, &ldb, ALPHAR, ALPHAI, BETA, VL, &ldvl, VR, &ldvr, WORK, &lwork, info); }
00367 
00368 
00369   void LAPACK<int, float>::TRSEN(const char JOB, const char COMPQ, const int *SELECT, const int n, float *T, const int ldt, float *Q, const int ldq, float *WR, float *WI, int *M, float *S, float *SEP, float *WORK, const int lwork, int *IWORK, const int liwork, int *info ) const
00370   { STRSEN_F77(CHAR_MACRO(JOB), CHAR_MACRO(COMPQ), SELECT, &n, T, &ldt, Q, &ldq, WR, WI, M, S, SEP, WORK, &lwork, IWORK, &liwork, info); }
00371 
00372 
00373   void LAPACK<int, float>::TGSEN(const int ijob, const int wantq, const int wantz, const int *SELECT, const int n, float *A, const int lda, float *B, const int ldb, float *ALPHAR, float *ALPHAI, float *BETA, float *Q, const int ldq, float *Z, const int ldz, int *M, float *PL, float *PR, float *DIF, float *WORK, const int lwork, int *IWORK, const int liwork, int *info ) const
00374   { STGSEN_F77(&ijob, &wantq, &wantz, SELECT, &n, A, &lda, B, &ldb, ALPHAR, ALPHAI, BETA, Q, &ldq, Z, &ldz, M, PL, PR, DIF, WORK, &lwork, IWORK, &liwork, info); }
00375 
00376 
00377   void LAPACK<int, float>::GGES(const char JOBVL, const char JOBVR, const char SORT, int (*ptr2func)(float *, float *, float *), const int n, float *A, const int lda, float *B, const int ldb, int *sdim, float *ALPHAR, float *ALPHAI, float *BETA, float *VL, const int ldvl, float *VR, const int ldvr, float *WORK, const int lwork, int *BWORK, int *info ) const
00378   { SGGES_F77(CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), CHAR_MACRO(SORT), ptr2func, &n, A, &lda, B, &ldb, sdim, ALPHAR, ALPHAI, BETA, VL, &ldvl, VR, &ldvr, WORK, &lwork, BWORK, info); }
00379 
00380 
00381   void LAPACK<int, float>::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
00382   { SORMQR_F77(CHAR_MACRO(SIDE), CHAR_MACRO(TRANS), &m, &n, &k, A, &lda, TAU, C, &ldc, WORK, &lwork, info); }
00383 
00384   void LAPACK<int, float>::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
00385   {
00386     // LAPACK only defines UNMQR for Z and C (complex*8
00387     // resp. complex*16), but logically, UNMQR means the same thing as
00388     // ORMQR for real arithmetic.
00389     ORMQR (SIDE, TRANS, m, n, k, A, lda, TAU, C, ldc, WORK, lwork, info);
00390   }
00391 
00392   void LAPACK<int, float>::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
00393   { SORGQR_F77( &m, &n, &k, A, &lda, TAU, WORK, &lwork, info); }
00394 
00395 
00396   void LAPACK<int, float>::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
00397   { SORGQR_F77( &m, &n, &k, A, &lda, TAU, WORK, &lwork, info); }
00398 
00399 
00400   void LAPACK<int, float>::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
00401   { SORGHR_F77(&n, &ilo, &ihi, A, &lda, TAU, WORK, &lwork, info); }
00402 
00403 
00404   void LAPACK<int, float>::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
00405   { SORMHR_F77(CHAR_MACRO(SIDE), CHAR_MACRO(TRANS), &m, &n, &ilo, &ihi, A, &lda, TAU, C, &ldc, WORK, &lwork, info); }
00406 
00407 
00408   void LAPACK<int, float>::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
00409   { STREVC_F77(CHAR_MACRO(SIDE), CHAR_MACRO(HOWMNY), select, &n, T, &ldt, VL, &ldvl, VR, &ldvr, &mm, m, WORK, info); }
00410 
00411 
00412   void LAPACK<int, float>::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
00413   {
00414     std::vector<int> select(1);
00415     const char whch = 'A';
00416     STREVC_F77(CHAR_MACRO(SIDE), CHAR_MACRO(whch), &select[0], &n, T, &ldt, VL, &ldvl, VR, &ldvr, &mm, m, WORK, info);
00417   }
00418 
00419 
00420   void LAPACK<int, float>::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
00421   { STREXC_F77(CHAR_MACRO(COMPQ), &n, T, &ldt, Q, &ldq, &ifst, &ilst, WORK, info); }
00422 
00423 
00424   void LAPACK<int, float>::TGEVC(const char SIDE, const char HOWMNY, const int *SELECT, const int n, float *S, const int lds, float *P, const int ldp, float *VL, const int ldvl, float *VR, const int ldvr, const int mm, int *M, float *WORK, int *info) const
00425   { STGEVC_F77(CHAR_MACRO(SIDE), CHAR_MACRO(HOWMNY), SELECT, &n, S, &lds, P, &ldp, VL, &ldvl, VR, &ldvr, &mm, M, WORK, info); }
00426 
00427 
00428   void LAPACK<int, float>::LARTG( const float f, const float g, float* c, float* s, float* r ) const
00429   { SLARTG_F77(&f, &g, c, s, r); }
00430 
00431 
00432   void LAPACK<int, float>::LARFG( const int n, float* alpha, float* x, const int incx, float* tau ) const
00433   { SLARFG_F77(&n, alpha, x, &incx, tau); }
00434 
00435 
00436   void LAPACK<int, float>::GEBAL(const char JOBZ, const int n, float* A, const int lda, int ilo, int ihi, float* scale, int* info) const
00437   { SGEBAL_F77(CHAR_MACRO(JOBZ),&n, A, &lda, &ilo, &ihi, scale, info); }
00438 
00439 
00440   void LAPACK<int, float>::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
00441   { SGEBAK_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(SIDE), &n, &ilo, &ihi, scale, &m, V, &ldv, info); }
00442 
00443 #ifdef HAVE_TEUCHOS_LAPACKLARND
00444   float LAPACK<int, float>::LARND( const int idist, int* seed ) const
00445   { return(SLARND_F77(&idist, seed)); }
00446 #endif
00447 
00448   void LAPACK<int, float>::LARNV( const int idist, int* seed, const int n, float* v ) const
00449   { SLARNV_F77(&idist, seed, &n, v); }
00450 
00451 
00452   float LAPACK<int, float>::LAMCH(const char CMACH) const
00453   { return(SLAMCH_F77(CHAR_MACRO(CMACH))); }
00454 
00455 
00456   int LAPACK<int, float>::ILAENV( const int ispec, const std::string& NAME, const std::string& OPTS, const int N1, const int N2, const int N3, const int N4 ) const
00457   {
00458     unsigned int opts_length = OPTS.length();
00459     // if user queries a Hermitian routine, change it to a symmetric routine
00460     std::string temp_NAME = "s" + NAME;
00461     if (temp_NAME.substr(1,2) == "he") {
00462       temp_NAME.replace(1,2,"sy");
00463     }
00464     unsigned int name_length = temp_NAME.length();
00465     return ilaenv_wrapper(&ispec, &temp_NAME[0], name_length, &OPTS[0], opts_length, &N1, &N2, &N3, &N4);
00466   }
00467 
00468 
00469   float LAPACK<int, float>::LAPY2(const float x, const float y) const
00470   {
00471 #if defined(HAVE_TEUCHOS_BLASFLOAT)
00472     return SLAPY2_F77(&x, &y);
00473 #else
00474     typedef ScalarTraits<float> ST;
00475     const float xabs = ST::magnitude(x);
00476     const float yabs = ST::magnitude(y);
00477     const float w = TEUCHOS_MAX(xabs, yabs);
00478     const float z = TEUCHOS_MIN(xabs, yabs);
00479     if ( z == 0.0 ) {
00480       return w;
00481     }
00482     const float z_over_w = z/w;
00483     return w*ST::squareroot( 1.0+(z_over_w*z_over_w));
00484 #endif
00485   }
00486 
00487   // END INT, FLOAT SPECIALIZATION IMPLEMENTATION //
00488 
00489   // BEGIN INT, DOUBLE SPECIALIZATION IMPLEMENTATION //
00490 
00491 
00492 
00493   void LAPACK<int, double>::PTTRF(const int n, double* d, double* e, int* info) const
00494   { DPTTRF_F77(&n,d,e,info); }
00495 
00496 
00497   void LAPACK<int, double>::PTTRS(const int n, const int nrhs, const double* d, const double* e, double* B, const int ldb, int* info) const
00498   { DPTTRS_F77(&n,&nrhs,d,e,B,&ldb,info); }
00499 
00500 
00501   void LAPACK<int, double>::POTRF(const char UPLO, const int n, double* A, const int lda, int* info) const
00502   { DPOTRF_F77(CHAR_MACRO(UPLO), &n, A, &lda, info); }
00503 
00504 
00505   void LAPACK<int, double>::POTRS(const char UPLO, const int n, const int nrhs, const double* A, const int lda, double* B, const int ldb, int* info) const
00506   { DPOTRS_F77(CHAR_MACRO(UPLO), &n, &nrhs, A, &lda, B, &ldb, info); }
00507 
00508 
00509   void LAPACK<int, double>::POTRI(const char UPLO, const int n, double* A, const int lda, int* info) const
00510   { DPOTRI_F77(CHAR_MACRO(UPLO), &n, A, &lda, info); }
00511 
00512 
00513   void LAPACK<int, double>::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
00514   { DPOCON_F77(CHAR_MACRO(UPLO), &n, A, &lda, &anorm, rcond, WORK, IWORK, info); }
00515 
00516 
00517   void LAPACK<int, double>::POSV(const char UPLO, const int n, const int nrhs, double* A, const int lda, double* B, const int ldb, int* info) const
00518   { DPOSV_F77(CHAR_MACRO(UPLO), &n, &nrhs, A, &lda, B, &ldb, info); }
00519 
00520 
00521   void LAPACK<int, double>::POEQU(const int n, const double* A, const int lda, double* S, double* scond, double* amax, int* info) const
00522   { DPOEQU_F77(&n, A, &lda, S, scond, amax, info); }
00523 
00524 
00525   void LAPACK<int, double>::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
00526   { DPORFS_F77(CHAR_MACRO(UPLO), &n, &nrhs, A, &lda, AF, &ldaf, B, &ldb, X, &ldx, FERR, BERR, WORK, IWORK, info); }
00527 
00528 
00529   void LAPACK<int, double>::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
00530   { DPOSVX_F77(CHAR_MACRO(FACT), CHAR_MACRO(UPLO), &n, &nrhs, A, &lda, AF, &ldaf, CHAR_MACRO(EQUED), S, B, &ldb, X, &ldx, rcond, FERR, BERR, WORK, IWORK, info); }
00531 
00532 
00533   void LAPACK<int,double>::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
00534   { DGELS_F77(CHAR_MACRO(TRANS), &m, &n, &nrhs, A, &lda, B, &ldb, WORK, &lwork, info); }
00535 
00536 
00537   void LAPACK<int,double>::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
00538   {
00539     (void) rwork;
00540     DGELSS_F77(&m, &n, &nrhs, A, &lda, B, &ldb, S, &rcond, rank, WORK, &lwork, info);
00541   }
00542 
00543 
00544   void LAPACK<int,double>::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
00545   { DGELSS_F77(&m, &n, &nrhs, A, &lda, B, &ldb, S, &rcond, rank, WORK, &lwork, info); }
00546 
00547 
00548   void LAPACK<int,double>::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
00549   { DGGLSE_F77(&m, &n, &p, A, &lda, B, &ldb, C, D, X, WORK, &lwork, info); }
00550 
00551 
00552   void LAPACK<int,double>::GEQRF( const int m, const int n, double* A, const int lda, double* TAU, double* WORK, const int lwork, int* info) const
00553   { DGEQRF_F77(&m, &n, A, &lda, TAU, WORK, &lwork, info); }
00554 
00555 
00556   void LAPACK<int,double>::GETRF(const int m, const int n, double* A, const int lda, int* IPIV, int* info) const
00557   { DGETRF_F77(&m, &n, A, &lda, IPIV, info); }
00558 
00559 
00560   void LAPACK<int,double>::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
00561   { DGETRS_F77(CHAR_MACRO(TRANS), &n, &nrhs, A, &lda, IPIV, B, &ldb, info); }
00562 
00563 
00564   void LAPACK<int,double>::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
00565   { DLASCL_F77(CHAR_MACRO(TYPE), &kl, &ku, &cfrom, &cto, &m, &n, A, &lda, info); }
00566 
00567   void LAPACK<int,double>::GEQP3(const int m, const int n, double* A, const int lda, int *jpvt, double* TAU, double* WORK, const int lwork, double* RWORK, int* info ) const
00568   {
00569     (void) RWORK;
00570     DGEQP3_F77(&m, &n, A, &lda, jpvt, TAU, WORK, &lwork, info);
00571   }
00572 
00573   void LAPACK<int, double>::
00574   LASWP (const int N,
00575          double A[],
00576          const int LDA,
00577          const int K1,
00578          const int K2,
00579          const int IPIV[],
00580          const int INCX) const
00581   {
00582     DLASWP_F77(&N, A, &LDA, &K1, &K2, IPIV, &INCX);
00583   }
00584 
00585   void LAPACK<int,double>::GBTRF(const int m, const int n, const int kl, const int ku, double* A, const int lda, int* IPIV, int* info) const
00586   { DGBTRF_F77(&m, &n, &kl, &ku, A, &lda, IPIV, info); }
00587 
00588 
00589   void LAPACK<int,double>::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
00590   { DGBTRS_F77(CHAR_MACRO(TRANS), &n, &kl, &ku, &nrhs, A, &lda, IPIV, B, &ldb, info); }
00591 
00592 
00593   void LAPACK<int,double>::GTTRF(const int n, double* dl, double* d, double* du, double* du2, int* IPIV, int* info) const
00594   { DGTTRF_F77(&n, dl, d, du, du2, IPIV, info); }
00595 
00596 
00597   void LAPACK<int,double>::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
00598   { DGTTRS_F77(CHAR_MACRO(TRANS), &n, &nrhs, dl, d, du, du2, IPIV, B, &ldb, info); }
00599 
00600 
00601   void LAPACK<int,double>::GETRI(const int n, double* A, const int lda, const int* IPIV, double* WORK, const int lwork, int* info) const
00602   { DGETRI_F77(&n, A, &lda, IPIV, WORK, &lwork, info); }
00603 
00604   void
00605   LAPACK<int, double>::LATRS (const char UPLO,
00606                               const char TRANS,
00607                               const char DIAG,
00608                               const char NORMIN,
00609                               const int N,
00610                               double* A,
00611                               const int LDA,
00612                               double* X,
00613                               double* SCALE,
00614                               double* CNORM,
00615                               int* INFO) const
00616   {
00617     DLATRS_F77(CHAR_MACRO(UPLO),
00618                CHAR_MACRO(TRANS),
00619                CHAR_MACRO(DIAG),
00620                CHAR_MACRO(NORMIN),
00621                &N,
00622                A,
00623                &LDA,
00624                X,
00625                SCALE,
00626                CNORM,
00627                INFO);
00628   }
00629 
00630   void LAPACK<int,double>::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
00631   { DGECON_F77(CHAR_MACRO(NORM), &n, A, &lda, &anorm, rcond, WORK, IWORK, info); }
00632 
00633 
00634   void LAPACK<int,double>::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
00635   { DGBCON_F77(CHAR_MACRO(NORM), &n, &kl, &ku, A, &lda, IPIV, &anorm, rcond, WORK, IWORK, info); }
00636 
00637 
00638   double LAPACK<int,double>::LANGB(const char NORM, const int n, const int kl, const int ku, const double* A, const int lda, double* WORK) const
00639   { return( DLANGB_F77(CHAR_MACRO(NORM), &n, &kl, &ku, A, &lda, WORK) ); }
00640 
00641 
00642   void LAPACK<int,double>::GESV(const int n, const int nrhs, double* A, const int lda, int* IPIV, double* B, const int ldb, int* info) const
00643   { DGESV_F77(&n, &nrhs, A, &lda, IPIV, B, &ldb, info); }
00644 
00645 
00646   void LAPACK<int,double>::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
00647   { DGEEQU_F77(&m, &n, A, &lda, R, C, rowcond, colcond, amax, info); }
00648 
00649 
00650   void LAPACK<int,double>::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
00651   { DGERFS_F77(CHAR_MACRO(TRANS), &n, &nrhs, A, &lda, AF, &ldaf, IPIV, B, &ldb, X, &ldx, FERR, BERR, WORK, IWORK, info); }
00652 
00653 
00654   void LAPACK<int,double>::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
00655   { DGBEQU_F77(&m, &n, &kl, &ku, A, &lda, R, C, rowcond, colcond, amax, info); }
00656 
00657 
00658   void LAPACK<int,double>::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
00659   { DGBRFS_F77(CHAR_MACRO(TRANS), &n, &kl, &ku, &nrhs, A, &lda, AF, &ldaf, IPIV, B, &ldb, X, &ldx, FERR, BERR, WORK, IWORK, info); }
00660 
00661 
00662   void LAPACK<int,double>::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
00663   { DGESVX_F77(CHAR_MACRO(FACT), CHAR_MACRO(TRANS), &n, &nrhs, A, &lda, AF, &ldaf, IPIV, CHAR_MACRO(EQUED), R, C, B, &ldb, X, &ldx, rcond, FERR, BERR, WORK, IWORK, info); }
00664 
00665 
00666   void LAPACK<int,double>::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
00667   { DSYTRD_F77(CHAR_MACRO(UPLO), &n, A, &lda, D, E, TAU, WORK, &lwork, info); }
00668 
00669 
00670   void LAPACK<int, double>::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
00671   { DGEHRD_F77(&n, &ilo, &ihi, A, &lda, TAU, WORK, &lwork, info); }
00672 
00673 
00674   void LAPACK<int,double>::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
00675   { DTRTRS_F77(CHAR_MACRO(UPLO), CHAR_MACRO(TRANS), CHAR_MACRO(DIAG), &n, &nrhs, A, &lda, B, &ldb, info); }
00676 
00677 
00678   void LAPACK<int,double>::TRTRI(const char UPLO, const char DIAG, const int n, const double* A, const int lda, int* info) const
00679   { DTRTRI_F77(CHAR_MACRO(UPLO), CHAR_MACRO(DIAG), &n, A, &lda, info); }
00680 
00681 
00682   void LAPACK<int,double>::SPEV(const char JOBZ, const char UPLO, const int n, double* AP, double* W, double* Z, const int ldz, double* WORK, int* info) const
00683   { DSPEV_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &n, AP, W, Z, &ldz, WORK, info); }
00684 
00685 
00686   void LAPACK<int,double>::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
00687   {
00688     DSYEV_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &n, A, &lda, W, WORK, &lwork, info);
00689   }
00690 
00691 
00692   void LAPACK<int,double>::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
00693   {
00694     DSYGV_F77(&itype, CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &n, A, &lda, B, &ldb, W, WORK, &lwork, info);
00695   }
00696 
00697 
00698   void LAPACK<int,double>::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
00699   {
00700     DSYEV_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &n, A, &lda, W, WORK, &lwork, info);
00701   }
00702 
00703 
00704   void LAPACK<int,double>::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
00705   {
00706     DSYGV_F77(&itype, CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &n, A, &lda, B, &ldb, W, WORK, &lwork, info);
00707   }
00708 
00709 
00710   void LAPACK<int,double>::STEQR(const char COMPZ, const int n, double* D, double* E, double* Z, const int ldz, double* WORK, int* info) const
00711   {
00712     DSTEQR_F77(CHAR_MACRO(COMPZ), &n, D, E, Z, &ldz, WORK, info);
00713   }
00714 
00715 
00716   void LAPACK<int, double>::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
00717   {
00718     DHSEQR_F77(CHAR_MACRO(JOB), CHAR_MACRO(COMPZ), &n, &ilo, &ihi, H, &ldh, WR, WI, Z, &ldz, WORK, &lwork, info);
00719   }
00720 
00721 
00722   void LAPACK<int, double>::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
00723   {
00724     DGEES_F77(CHAR_MACRO(JOBVS), CHAR_MACRO(SORT), ptr2func, &n, A, &lda, sdim, WR, WI, VS, &ldvs, WORK, &lwork, BWORK, info);
00725   }
00726 
00727 
00728   void LAPACK<int, double>::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
00729   {
00730     //int (*nullfptr)(double*,double*) = NULL;
00731     gees_nullfptr_t nullfptr = 0;
00732     const char sort = 'N';
00733     DGEES_F77(CHAR_MACRO(JOBVS), CHAR_MACRO(sort), nullfptr, &n, A, &lda, sdim, WR, WI, VS, &ldvs, WORK, &lwork, BWORK, info);
00734   }
00735 
00736 
00737   void LAPACK<int, double>::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
00738   {
00739     DGEEV_F77(CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), &n, A, &lda, WR, WI, VL, &ldvl, VR, &ldvr, WORK, &lwork, info);
00740   }
00741 
00742   void LAPACK<int, double>::GEEV(const char JOBVL, const char JOBVR, const int n, double* A, const int lda, double* WR, double* WI, double* VL, const int ldvl, double* VR, const int ldvr, double* WORK, const int lwork, double* /* RWORK */, int* info) const
00743   {
00744     GEEV (JOBVL, JOBVR, n, A, lda, WR, WI, VL, ldvl, VR, ldvr, WORK, lwork, info);
00745   }
00746 
00747 
00748   void LAPACK<int, double>::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 {
00749     DGESVD_F77(CHAR_MACRO(JOBU), CHAR_MACRO(JOBVT), &m, &n, A, &lda, S, U, &ldu, V, &ldv, WORK, &lwork, info);
00750   }
00751 
00752 
00753   void LAPACK<int,double>::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
00754   {
00755     DGEEVX_F77(CHAR_MACRO(BALANC), CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), CHAR_MACRO(SENSE), &n, A, &lda, WR, WI, VL, &ldvl, VR, &ldvr, ilo, ihi, SCALE, abnrm, RCONDE, RCONDV, WORK, &lwork, IWORK, info);
00756   }
00757 
00758 
00759   void LAPACK<int, double>::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
00760   {
00761     DGGEVX_F77(CHAR_MACRO(BALANC), CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), CHAR_MACRO(SENSE), &n, A, &lda, B, &ldb, ALPHAR, ALPHAI, BETA, VL, &ldvl, VR, &ldvr, ilo, ihi, lscale, rscale, abnrm, bbnrm, RCONDE, RCONDV, WORK, &lwork, IWORK, BWORK, info);
00762   }
00763 
00764   void LAPACK<int, double>::GGEVX(const char BALANC, const char JOBVL, const char JOBVR, const char SENSE, const int n, double* A, const int lda, double* B, const int ldb, double* ALPHAR, double* ALPHAI, double* BETA, double* VL, const int ldvl, double* VR, const int ldvr, int* ilo, int* ihi, double* lscale, double* rscale, double* abnrm, double* bbnrm, double* RCONDE, double* RCONDV, double* WORK, const int lwork, double* /* RWORK */, int* IWORK, int* BWORK, int* info) const
00765   {
00766     GGEVX(BALANC, JOBVL, JOBVR, SENSE, n, A, lda, B, ldb, ALPHAR, ALPHAI, BETA, VL, ldvl, VR, ldvr, ilo, ihi, lscale, rscale, abnrm, bbnrm, RCONDE, RCONDV, WORK, lwork, IWORK, BWORK, info);
00767   }
00768 
00769   void LAPACK<int, double>::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
00770   {
00771     DGGEV_F77(CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), &n, A, &lda, B, &ldb, ALPHAR, ALPHAI, BETA, VL, &ldvl, VR, &ldvr, WORK, &lwork, info);
00772   }
00773 
00774   void LAPACK<int, double>::TRSEN(const char JOB, const char COMPQ, const int *SELECT, const int n, double *T, const int ldt, double *Q, const int ldq, double *WR, double *WI, int *M, double *S, double *SEP, double *WORK, const int lwork, int *IWORK, const int liwork, int *info ) const
00775   { DTRSEN_F77(CHAR_MACRO(JOB), CHAR_MACRO(COMPQ), SELECT, &n, T, &ldt, Q, &ldq, WR, WI, M, S, SEP, WORK, &lwork, IWORK, &liwork, info); }
00776 
00777 
00778   void LAPACK<int, double>::TGSEN(const int ijob, const int wantq, const int wantz, const int *SELECT, const int n, double *A, const int lda, double *B, const int ldb, double *ALPHAR, double *ALPHAI, double *BETA, double *Q, const int ldq, double *Z, const int ldz, int *M, double *PL, double *PR, double *DIF, double *WORK, const int lwork, int *IWORK, const int liwork, int *info ) const
00779   { DTGSEN_F77(&ijob, &wantq, &wantz, SELECT, &n, A, &lda, B, &ldb, ALPHAR, ALPHAI, BETA, Q, &ldq, Z, &ldz, M, PL, PR, DIF, WORK, &lwork, IWORK, &liwork, info); }
00780 
00781 
00782   void LAPACK<int, double>::GGES(const char JOBVL, const char JOBVR, const char SORT, int (*ptr2func)(double *, double *, double *), const int n, double *A, const int lda, double *B, const int ldb, int *sdim, double *ALPHAR, double *ALPHAI, double *BETA, double *VL, const int ldvl, double *VR, const int ldvr, double *WORK, const int lwork, int *BWORK, int *info ) const
00783   { DGGES_F77(CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), CHAR_MACRO(SORT), ptr2func, &n, A, &lda, B, &ldb, sdim, ALPHAR, ALPHAI, BETA, VL, &ldvl, VR, &ldvr, WORK, &lwork, BWORK, info); }
00784 
00785 
00786   void LAPACK<int, double>::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
00787   {
00788     DORMQR_F77(CHAR_MACRO(SIDE), CHAR_MACRO(TRANS), &m, &n, &k, A, &lda, TAU, C, &ldc, WORK, &lwork, info);
00789   }
00790 
00791   void LAPACK<int, double>::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
00792   {
00793     // LAPACK only defines UNMQR for Z and C (complex*8
00794     // resp. complex*16), but logically, UNMQR means the same thing as
00795     // ORMQR for real arithmetic.
00796     ORMQR (SIDE, TRANS, m, n, k, A, lda, TAU, C, ldc, WORK, lwork, info);
00797   }
00798 
00799   void LAPACK<int, double>::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
00800   {
00801     DORGQR_F77( &m, &n, &k, A, &lda, TAU, WORK, &lwork, info);
00802   }
00803 
00804 
00805   void LAPACK<int, double>::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
00806   {
00807     DORGQR_F77( &m, &n, &k, A, &lda, TAU, WORK, &lwork, info);
00808   }
00809 
00810 
00811   void LAPACK<int, double>::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
00812   {
00813     DORGHR_F77(&n, &ilo, &ihi, A, &lda, TAU, WORK, &lwork, info);
00814   }
00815 
00816 
00817   void LAPACK<int, double>::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
00818   {
00819     DORMHR_F77(CHAR_MACRO(SIDE), CHAR_MACRO(TRANS), &m, &n, &ilo, &ihi, A, &lda, TAU, C, &ldc, WORK, &lwork, info);
00820   }
00821 
00822 
00823   void LAPACK<int, double>::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
00824   {
00825     DTREVC_F77(CHAR_MACRO(SIDE), CHAR_MACRO(HOWMNY), select, &n, T, &ldt, VL, &ldvl, VR, &ldvr, &mm, m, WORK, info);
00826   }
00827 
00828 
00829   void LAPACK<int, double>::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
00830   {
00831     std::vector<int> select(1);
00832     const char whch = 'A';
00833     DTREVC_F77(CHAR_MACRO(SIDE), CHAR_MACRO(whch), &select[0], &n, T, &ldt, VL, &ldvl, VR, &ldvr, &mm, m, WORK, info);
00834   }
00835 
00836 
00837   void LAPACK<int, double>::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
00838   {
00839     DTREXC_F77(CHAR_MACRO(COMPQ), &n, T, &ldt, Q, &ldq, &ifst, &ilst, WORK, info);
00840   }
00841 
00842 
00843   void LAPACK<int, double>::TGEVC(const char SIDE, const char HOWMNY, const int *SELECT, const int n, double *S, const int lds, double *P, const int ldp, double *VL, const int ldvl, double *VR, const int ldvr, const int mm, int *M, double *WORK, int *info) const
00844   { DTGEVC_F77(CHAR_MACRO(SIDE), CHAR_MACRO(HOWMNY), SELECT, &n, S, &lds, P, &ldp, VL, &ldvl, VR, &ldvr, &mm, M, WORK, info); }
00845 
00846 
00847   void LAPACK<int, double>::LARTG( const double f, const double g, double* c, double* s, double* r ) const
00848   {
00849     DLARTG_F77(&f, &g, c, s, r);
00850   }
00851 
00852 
00853   void LAPACK<int, double>::LARFG( const int n, double* alpha, double* x, const int incx, double* tau ) const
00854   {
00855     DLARFG_F77(&n, alpha, x, &incx, tau);
00856   }
00857 
00858 
00859   void LAPACK<int, double>::GEBAL(const char JOBZ, const int n, double* A, const int lda, int ilo, int ihi, double* scale, int* info) const
00860   {
00861     DGEBAL_F77(CHAR_MACRO(JOBZ),&n, A, &lda, &ilo, &ihi, scale, info);
00862   }
00863 
00864 
00865   void LAPACK<int, double>::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
00866   {
00867     DGEBAK_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(SIDE), &n, &ilo, &ihi, scale, &m, V, &ldv, info);
00868   }
00869 
00870 
00871 #ifdef HAVE_TEUCHOS_LAPACKLARND
00872   double LAPACK<int, double>::LARND( const int idist, int* seed ) const
00873   {
00874     return(DLARND_F77(&idist, seed));
00875   }
00876 #endif
00877 
00878   void LAPACK<int, double>::LARNV( const int idist, int* seed, const int n, double* v ) const
00879   {
00880     DLARNV_F77(&idist, seed, &n, v);
00881   }
00882 
00883 
00884   double LAPACK<int, double>::LAMCH(const char CMACH) const
00885   {
00886     return(DLAMCH_F77(CHAR_MACRO(CMACH)));
00887   }
00888 
00889 
00890   int LAPACK<int, double>::ILAENV( const int ispec, const std::string& NAME, const std::string& OPTS, const int N1, const int N2, const int N3, const int N4 ) const
00891   {
00892     unsigned int opts_length = OPTS.length();
00893     // if user queries a Hermitian routine, change it to a symmetric routine
00894     std::string temp_NAME = "d" + NAME;
00895     if (temp_NAME.substr(1,2) == "he") {
00896       temp_NAME.replace(1,2,"sy");
00897     }
00898     unsigned int name_length = temp_NAME.length();
00899     return ilaenv_wrapper(&ispec, &temp_NAME[0], name_length, &OPTS[0], opts_length, &N1, &N2, &N3, &N4);
00900   }
00901 
00902 
00903   double LAPACK<int, double>::LAPY2(const double x, const double y) const
00904   {
00905     return DLAPY2_F77(&x, &y);
00906   }
00907 
00908   // END INT, DOUBLE SPECIALIZATION IMPLEMENTATION //
00909 
00910 #ifdef HAVE_TEUCHOS_COMPLEX
00911 
00912   // BEGIN INT, COMPLEX<FLOAT> SPECIALIZATION IMPLEMENTATION //
00913 
00914 
00915   void LAPACK<int, std::complex<float> >::PTTRF(const int n, std::complex<float>* d, std::complex<float>* e, int* info) const
00916   {
00917     CPTTRF_F77(&n,d,e,info);
00918   }
00919 
00920 
00921   void LAPACK<int, std::complex<float> >::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
00922   {
00923     CPTTRS_F77(&n,&nrhs,d,e,B,&ldb,info);
00924   }
00925 
00926 
00927   void LAPACK<int, std::complex<float> >::POTRF(const char UPLO, const int n, std::complex<float>* A, const int lda, int* info) const
00928   {
00929     CPOTRF_F77(CHAR_MACRO(UPLO), &n, A, &lda, info);
00930   }
00931 
00932 
00933   void LAPACK<int, std::complex<float> >::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
00934   {
00935     CPOTRS_F77(CHAR_MACRO(UPLO), &n, &nrhs, A, &lda, B, &ldb, info);
00936   }
00937 
00938 
00939   void LAPACK<int, std::complex<float> >::POTRI(const char UPLO, const int n, std::complex<float>* A, const int lda, int* info) const
00940   {
00941     CPOTRI_F77(CHAR_MACRO(UPLO), &n, A, &lda, info);
00942   }
00943 
00944 
00945   void LAPACK<int, std::complex<float> >::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
00946   {
00947     CPOCON_F77(CHAR_MACRO(UPLO), &n, A, &lda, &anorm, rcond, WORK, RWORK, info);
00948   }
00949 
00950 
00951   void LAPACK<int, std::complex<float> >::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
00952   {
00953     CPOSV_F77(CHAR_MACRO(UPLO), &n, &nrhs, A, &lda, B, &ldb, info);
00954   }
00955 
00956 
00957   void LAPACK<int, std::complex<float> >::POEQU(const int n, const std::complex<float>* A, const int lda, float* S, float* scond, float* amax, int* info) const
00958   {
00959     CPOEQU_F77(&n, A, &lda, S, scond, amax, info);
00960   }
00961 
00962 
00963   void LAPACK<int, std::complex<float> >::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
00964   {
00965     CPORFS_F77(CHAR_MACRO(UPLO), &n, &nrhs, A, &lda, AF, &ldaf, B, &ldb, X, &ldx, FERR, BERR, WORK, RWORK, info);
00966   }
00967 
00968 
00969   void LAPACK<int, std::complex<float> >::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
00970   {
00971     CPOSVX_F77(CHAR_MACRO(FACT), CHAR_MACRO(UPLO), &n, &nrhs, A, &lda, AF, &ldaf, CHAR_MACRO(EQUED), S, B, &ldb, X, &ldx, rcond, FERR, BERR, WORK, RWORK, info);
00972   }
00973 
00974 
00975   void LAPACK<int,std::complex<float> >::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
00976   {
00977     CGELS_F77(CHAR_MACRO(TRANS), &m, &n, &nrhs, A, &lda, B, &ldb, WORK, &lwork, info);
00978   }
00979 
00980 
00981   void LAPACK<int, std::complex<float> >::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
00982   {
00983     CGELSS_F77(&m, &n, &nrhs, A, &lda, B, &ldb, S, &rcond, rank, WORK, &lwork, rwork, info);
00984   }
00985 
00986 
00987   void LAPACK<int,std::complex<float> >::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
00988   {
00989     CGEQRF_F77(&m, &n, A, &lda, TAU, WORK, &lwork, info);
00990   }
00991 
00992 
00993   void LAPACK<int,std::complex<float> >::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
00994   {
00995     CUNGQR_F77( &m, &n, &k, A, &lda, TAU, WORK, &lwork, info);
00996   }
00997 
00998 
00999   void LAPACK<int,std::complex<float> >::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
01000   {
01001     CUNMQR_F77(CHAR_MACRO(SIDE), CHAR_MACRO(TRANS), &m, &n, &k, A, &lda, TAU, C, &ldc, WORK, &lwork, info);
01002   }
01003 
01004 
01005   void LAPACK<int,std::complex<float> >::GETRF(const int m, const int n, std::complex<float>* A, const int lda, int* IPIV, int* info) const
01006   {
01007     CGETRF_F77(&m, &n, A, &lda, IPIV, info);
01008   }
01009 
01010 
01011   void LAPACK<int,std::complex<float> >::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
01012   {
01013     CGETRS_F77(CHAR_MACRO(TRANS), &n, &nrhs, A, &lda, IPIV, B, &ldb, info);
01014   }
01015 
01016 
01017   void LAPACK<int,std::complex<float> >::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
01018   { CLASCL_F77(CHAR_MACRO(TYPE), &kl, &ku, &cfrom, &cto, &m, &n, A, &lda, info); }
01019 
01020   void LAPACK<int,std::complex<float> >::GEQP3(const int m, const int n, std::complex<float>* A, const int lda, int *jpvt, std::complex<float>* TAU, std::complex<float>* WORK, const int lwork, float* RWORK, int* info ) const
01021   {
01022     CGEQP3_F77(&m, &n, A, &lda, jpvt, TAU, WORK, &lwork, RWORK, info);
01023   }
01024 
01025   void LAPACK<int, std::complex<float> >::
01026   LASWP (const int N,
01027          std::complex<float> A[],
01028          const int LDA,
01029          const int K1,
01030          const int K2,
01031          const int IPIV[],
01032          const int INCX) const
01033   {
01034     CLASWP_F77(&N, A, &LDA, &K1, &K2, IPIV, &INCX);
01035   }
01036 
01037   void LAPACK<int,std::complex<float> >::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
01038   {
01039     CGBTRF_F77(&m, &kl, &ku, &n, A, &lda, IPIV, info);
01040   }
01041 
01042 
01043   void LAPACK<int,std::complex<float> >::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
01044   {
01045     CGBTRS_F77(CHAR_MACRO(TRANS), &n, &kl, &ku, &nrhs, A, &lda, IPIV, B, &ldb, info);
01046   }
01047 
01048 
01049   void LAPACK<int,std::complex<float> >::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
01050   {
01051     CGTTRF_F77(&n, dl, d, du, du2, IPIV, info);
01052   }
01053 
01054 
01055   void LAPACK<int,std::complex<float> >::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
01056   {
01057     CGTTRS_F77(CHAR_MACRO(TRANS), &n, &nrhs, dl, d, du, du2, IPIV, B, &ldb, info);
01058   }
01059 
01060 
01061   void LAPACK<int,std::complex<float> >::GETRI(const int n, std::complex<float>* A, const int lda, const int* IPIV, std::complex<float>* WORK, const int lwork, int* info) const
01062   {
01063     CGETRI_F77(&n, A, &lda, IPIV, WORK, &lwork, info);
01064   }
01065 
01066 
01067   void
01068   LAPACK<int, std::complex<float> >::LATRS (const char UPLO,
01069                                             const char TRANS,
01070                                             const char DIAG,
01071                                             const char NORMIN,
01072                                             const int N,
01073                                             std::complex<float>* A,
01074                                             const int LDA,
01075                                             std::complex<float>* X,
01076                                             float* SCALE,
01077                                             float* CNORM,
01078                                             int* INFO) const
01079   {
01080     CLATRS_F77(CHAR_MACRO(UPLO),
01081                CHAR_MACRO(TRANS),
01082                CHAR_MACRO(DIAG),
01083                CHAR_MACRO(NORMIN),
01084                &N,
01085                A,
01086                &LDA,
01087                X,
01088                SCALE,
01089                CNORM,
01090                INFO);
01091   }
01092 
01093 
01094   void LAPACK<int,std::complex<float> >::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
01095   {
01096     CGECON_F77(CHAR_MACRO(NORM), &n, A, &lda, &anorm, rcond, WORK, RWORK, info);
01097   }
01098 
01099 
01100   void LAPACK<int,std::complex<float> >::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
01101   {
01102     CGBCON_F77(CHAR_MACRO(NORM), &n, &kl, &ku, A, &lda, IPIV, &anorm, rcond, WORK, RWORK, info);
01103   }
01104 
01105 
01106   float LAPACK<int,std::complex<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
01107   {
01108     return( CLANGB_F77(CHAR_MACRO(NORM), &n, &kl, &ku, A, &lda, WORK) );
01109   }
01110 
01111 
01112   void LAPACK<int,std::complex<float> >::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
01113   {
01114     CGESV_F77(&n, &nrhs, A, &lda, IPIV, B, &ldb, info);
01115   }
01116 
01117 
01118   void LAPACK<int,std::complex<float> >::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
01119   {
01120     CGEEQU_F77(&m, &n, A, &lda, R, C, rowcond, colcond, amax, info);
01121   }
01122 
01123 
01124   void LAPACK<int,std::complex<float> >::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
01125   {
01126     CGERFS_F77(CHAR_MACRO(TRANS), &n, &nrhs, A, &lda, AF, &ldaf, IPIV, B, &ldb, X, &ldx, FERR, BERR, WORK, RWORK, info);
01127   }
01128 
01129 
01130   void LAPACK<int,std::complex<float> >::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
01131   {
01132     CGBEQU_F77(&m, &n, &kl, &ku, A, &lda, R, C, rowcond, colcond, amax, info);
01133   }
01134 
01135 
01136   void LAPACK<int,std::complex<float> >::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
01137   {
01138     CGBRFS_F77(CHAR_MACRO(TRANS), &n, &kl, &ku, &nrhs, A, &lda, AF, &ldaf, IPIV, B, &ldb, X, &ldx, FERR, BERR, WORK, RWORK, info);
01139   }
01140 
01141 
01142   void LAPACK<int,std::complex<float> >::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
01143   {
01144     CGESVX_F77(CHAR_MACRO(FACT), CHAR_MACRO(TRANS), &n, &nrhs, A, &lda, AF, &ldaf, IPIV, CHAR_MACRO(EQUED), R, C, B, &ldb, X, &ldx, rcond, FERR, BERR, WORK, RWORK, info);
01145   }
01146 
01147 
01148   void LAPACK<int,std::complex<float> >::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
01149   {
01150     CGEHRD_F77(&n, &ilo, &ihi, A, &lda, TAU, WORK, &lwork, info);
01151   }
01152 
01153 
01154   void LAPACK<int,std::complex<float> >::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
01155   {
01156     CTRTRS_F77(CHAR_MACRO(UPLO), CHAR_MACRO(TRANS), CHAR_MACRO(DIAG), &n, &nrhs, A, &lda, B, &ldb, info);
01157   }
01158 
01159 
01160   void LAPACK<int,std::complex<float> >::TRTRI(const char UPLO, const char DIAG, const int n, const std::complex<float>* A, const int lda, int* info) const
01161   {
01162     CTRTRI_F77(CHAR_MACRO(UPLO), CHAR_MACRO(DIAG), &n, A, &lda, info);
01163   }
01164 
01165 
01166   void LAPACK<int,std::complex<float> >::STEQR(const char COMPZ, const int n, float* D, float* E, std::complex<float>* Z, const int ldz, float* WORK, int* info) const
01167   {
01168     CSTEQR_F77(CHAR_MACRO(COMPZ), &n, D, E, Z, &ldz, WORK, info);
01169   }
01170 
01171 
01172   void LAPACK<int,std::complex<float> >::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
01173   {
01174     CHEEV_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &n, A, &lda, W, WORK, &lwork, RWORK, info);
01175   }
01176 
01177 
01178   void LAPACK<int,std::complex<float> >::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
01179   {
01180     CHEGV_F77(&itype, CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &n, A, &lda, B, &ldb, W, WORK, &lwork, RWORK, info);
01181   }
01182 
01183 
01184   void LAPACK<int, std::complex<float> >::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
01185   {
01186     CHSEQR_F77(CHAR_MACRO(JOB), CHAR_MACRO(COMPZ), &n, &ilo, &ihi, H, &ldh, W, Z, &ldz, WORK, &lwork, info);
01187   }
01188 
01189 
01190   void LAPACK<int, std::complex<float> >::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
01191   {
01192     CGEES_F77(CHAR_MACRO(JOBVS), CHAR_MACRO(SORT), ptr2func, &n, A, &lda, sdim, W, VS, &ldvs, WORK, &lwork, RWORK, BWORK, info);
01193   }
01194 
01195 
01196   void LAPACK<int, std::complex<float> >::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
01197   {
01198     int (*nullfptr)(std::complex<float>*) = NULL;
01199     std::vector< std::complex<float> > W(n);
01200     const char sort = 'N';
01201     CGEES_F77(CHAR_MACRO(JOBVS), CHAR_MACRO(sort), nullfptr, &n, A, &lda, sdim, &W[0], VS, &ldvs, WORK, &lwork, RWORK, BWORK, info);
01202     for (int i=0; i<n; i++) {
01203       WR[i] = W[i].real();
01204       WI[i] = W[i].imag();
01205     }
01206   }
01207 
01208 
01209   void LAPACK<int, std::complex<float> >::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
01210   {
01211     CGEEV_F77(CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), &n, A, &lda, W, VL, &ldvl, VR, &ldvr, WORK, &lwork, RWORK, info);
01212   }
01213 
01214   void LAPACK<int, std::complex<float> >::GEEV(const char JOBVL, const char JOBVR, const int n, std::complex<float>* A, const int lda, float* WR, float* WI, std::complex<float>* VL, const int ldvl, std::complex<float>* VR, const int ldvr, std::complex<float>* WORK, const int lwork, float* RWORK, int* info) const
01215   {
01216     std::vector<std::complex<float> > w (n);
01217     std::complex<float>* w_rawPtr = (n == 0) ? NULL : &w[0];
01218     GEEV (JOBVL, JOBVR, n, A, lda, w_rawPtr, VL, ldvl, VR, ldvr, WORK, lwork, RWORK, info);
01219     if (info == 0) {
01220       // The eigenvalues are only valid on output if INFO is zero.
01221       // Otherwise, we shouldn't even write to WR or WI.
01222       for (int k = 0; k < n; ++k) {
01223   WR[k] = w[k].real ();
01224   WI[k] = w[k].imag ();
01225       }
01226     }
01227   }
01228 
01229   void LAPACK<int, std::complex<float> >::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 {
01230     CGESVD_F77(CHAR_MACRO(JOBU), CHAR_MACRO(JOBVT), &m, &n, A, &lda, S, U, &ldu, V, &ldv, WORK, &lwork, RWORK, info);
01231   }
01232 
01233 
01234   void LAPACK<int, std::complex<float> >::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
01235   {
01236     CGEEVX_F77(CHAR_MACRO(BALANC), CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), CHAR_MACRO(SENSE), &n, A, &lda, W, VL, &ldvl, VR, &ldvr, ilo, ihi, SCALE, abnrm, RCONDE, RCONDV, WORK, &lwork, RWORK, info);
01237   }
01238 
01239 
01240   void LAPACK<int, std::complex<float> >::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
01241   {
01242     CGGEVX_F77(CHAR_MACRO(BALANC), CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), CHAR_MACRO(SENSE), &n, A, &lda, B, &ldb, ALPHA, BETA, VL, &ldvl, VR, &ldvr, ilo, ihi, lscale, rscale, abnrm, bbnrm, RCONDE, RCONDV, WORK, &lwork, RWORK, IWORK, BWORK, info);
01243   }
01244 
01245   void LAPACK<int, std::complex<float> >::GGEVX(const char BALANC, const char JOBVL, const char JOBVR, const char SENSE, const int n, std::complex<float>* A, const int lda, std::complex<float>* B, const int ldb, float* ALPHAR, float* ALPHAI, std::complex<float>* BETA, std::complex<float>* VL, const int ldvl, std::complex<float>* VR, const int ldvr, int* ilo, int* ihi, float* lscale, float* rscale, float* abnrm, float* bbnrm, float* RCONDE, float* RCONDV, std::complex<float>* WORK, const int lwork, float* RWORK, int* IWORK, int* BWORK, int* info) const
01246   {
01247     std::vector<std::complex<float> > w (n);
01248     std::complex<float>* w_rawPtr = (n == 0) ? NULL : &w[0];
01249     GGEVX(BALANC, JOBVL, JOBVR, SENSE, n, A, lda, B, ldb, w_rawPtr, BETA, VL, ldvl, VR, ldvr, ilo, ihi, lscale, rscale, abnrm, bbnrm, RCONDE, RCONDV, WORK, lwork, RWORK, IWORK, BWORK, info);
01250     if (info == 0) {
01251       // The eigenvalues are only valid on output if INFO is zero.
01252       // Otherwise, we shouldn't even write to WR or WI.
01253       for (int k = 0; k < n; ++k) {
01254   ALPHAR[k] = w[k].real ();
01255   ALPHAI[k] = w[k].imag ();
01256       }
01257     }
01258   }
01259 
01260 
01261   void LAPACK<int, std::complex<float> >::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
01262   {
01263     CTREVC_F77(CHAR_MACRO(SIDE), CHAR_MACRO(HOWMNY), select, &n, T, &ldt, VL, &ldvl, VR, &ldvr, &mm, m, WORK, RWORK, info);
01264   }
01265 
01266 
01267   void LAPACK<int, std::complex<float> >::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
01268   {
01269     std::vector<int> select(1);
01270     const char whch = 'A';
01271     CTREVC_F77(CHAR_MACRO(SIDE), CHAR_MACRO(whch), &select[0], &n, T, &ldt, VL, &ldvl, VR, &ldvr, &mm, m, WORK, RWORK, info);
01272   }
01273 
01274 
01275   void LAPACK<int, std::complex<float> >::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
01276   {
01277     CTREXC_F77(CHAR_MACRO(COMPQ), &n, T, &ldt, Q, &ldq, &ifst, &ilst, info);
01278   }
01279 
01280 
01281   void LAPACK<int, std::complex<float> >::LARTG( const std::complex<float> f, const std::complex<float> g, float* c, std::complex<float>* s, std::complex<float>* r ) const
01282   {
01283     CLARTG_F77(&f, &g, c, s, r);
01284   }
01285 
01286 
01287   void LAPACK<int, std::complex<float> >::LARFG( const int n, std::complex<float>* alpha, std::complex<float>* x, const int incx, std::complex<float>* tau ) const
01288   {
01289     CLARFG_F77(&n, alpha, x, &incx, tau);
01290   }
01291 
01292 
01293   void LAPACK<int, std::complex<float> >::GEBAL(const char JOBZ, const int n, std::complex<float>* A, const int lda, int ilo, int ihi, float* scale, int* info) const
01294   {
01295     CGEBAL_F77(CHAR_MACRO(JOBZ),&n, A, &lda, &ilo, &ihi, scale, info);
01296   }
01297 
01298 
01299   void LAPACK<int, std::complex<float> >::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
01300   {
01301     CGEBAK_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(SIDE), &n, &ilo, &ihi, scale, &m, V, &ldv, info);
01302   }
01303 
01304 
01305 #ifdef HAVE_TEUCHOS_LAPACKLARND
01306   std::complex<float> LAPACK<int, std::complex<float> >::LARND( const int idist, int* seed ) const
01307   {
01308     return(CLARND_F77(&idist, seed));
01309   }
01310 #endif
01311 
01312   void LAPACK<int, std::complex<float> >::LARNV( const int idist, int* seed, const int n, std::complex<float>* v ) const
01313   {
01314     CLARNV_F77(&idist, seed, &n, v);
01315   }
01316 
01317 
01318   int LAPACK<int, std::complex<float> >::ILAENV( const int ispec, const std::string& NAME, const std::string& OPTS, const int N1, const int N2, const int N3, const int N4 ) const
01319   {
01320     unsigned int opts_length = OPTS.length();
01321     std::string temp_NAME = "c" + NAME;
01322     unsigned int name_length = temp_NAME.length();
01323     return ilaenv_wrapper(&ispec, &temp_NAME[0], name_length, &OPTS[0], opts_length, &N1, &N2, &N3, &N4);
01324   }
01325 
01326   // END INT, COMPLEX<FLOAT> SPECIALIZATION IMPLEMENTATION //
01327 
01328   // BEGIN INT, COMPLEX<DOUBLE> SPECIALIZATION IMPLEMENTATION //
01329 
01330 
01331   void LAPACK<int, std::complex<double> >::PTTRF(const int n, std::complex<double>* d, std::complex<double>* e, int* info) const
01332   {
01333     ZPTTRF_F77(&n,d,e,info);
01334   }
01335 
01336 
01337   void LAPACK<int, std::complex<double> >::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
01338   {
01339     ZPTTRS_F77(&n,&nrhs,d,e,B,&ldb,info);
01340   }
01341 
01342 
01343   void LAPACK<int, std::complex<double> >::POTRF(const char UPLO, const int n, std::complex<double>* A, const int lda, int* info) const
01344   {
01345     ZPOTRF_F77(CHAR_MACRO(UPLO), &n, A, &lda, info);
01346   }
01347 
01348 
01349   void LAPACK<int, std::complex<double> >::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
01350   {
01351     ZPOTRS_F77(CHAR_MACRO(UPLO), &n, &nrhs, A, &lda, B, &ldb, info);
01352   }
01353 
01354 
01355   void LAPACK<int, std::complex<double> >::POTRI(const char UPLO, const int n, std::complex<double>* A, const int lda, int* info) const
01356   {
01357     ZPOTRI_F77(CHAR_MACRO(UPLO), &n, A, &lda, info);
01358   }
01359 
01360 
01361   void LAPACK<int, std::complex<double> >::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
01362   {
01363     ZPOCON_F77(CHAR_MACRO(UPLO), &n, A, &lda, &anorm, rcond, WORK, RWORK, info);
01364   }
01365 
01366 
01367   void LAPACK<int, std::complex<double> >::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
01368   {
01369     ZPOSV_F77(CHAR_MACRO(UPLO), &n, &nrhs, A, &lda, B, &ldb, info);
01370   }
01371 
01372 
01373   void LAPACK<int, std::complex<double> >::POEQU(const int n, const std::complex<double>* A, const int lda, double* S, double* scond, double* amax, int* info) const
01374   {
01375     ZPOEQU_F77(&n, A, &lda, S, scond, amax, info);
01376   }
01377 
01378 
01379   void LAPACK<int, std::complex<double> >::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
01380   {
01381     ZPORFS_F77(CHAR_MACRO(UPLO), &n, &nrhs, A, &lda, AF, &ldaf, B, &ldb, X, &ldx, FERR, BERR, WORK, RWORK, info);
01382   }
01383 
01384 
01385   void LAPACK<int, std::complex<double> >::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
01386   {
01387     ZPOSVX_F77(CHAR_MACRO(FACT), CHAR_MACRO(UPLO), &n, &nrhs, A, &lda, AF, &ldaf, CHAR_MACRO(EQUED), S, B, &ldb, X, &ldx, rcond, FERR, BERR, WORK, RWORK, info);
01388   }
01389 
01390 
01391   void LAPACK<int,std::complex<double> >::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
01392   {
01393     ZGELS_F77(CHAR_MACRO(TRANS), &m, &n, &nrhs, A, &lda, B, &ldb, WORK, &lwork, info);
01394   }
01395 
01396 
01397   void LAPACK<int, std::complex<double> >::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
01398   {
01399     ZGELSS_F77(&m, &n, &nrhs, A, &lda, B, &ldb, S, &rcond, rank, WORK, &lwork, rwork, info);
01400   }
01401 
01402 
01403   void LAPACK<int,std::complex<double> >::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
01404   {
01405     ZGEQRF_F77(&m, &n, A, &lda, TAU, WORK, &lwork, info);
01406   }
01407 
01408 
01409   void LAPACK<int,std::complex<double> >::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
01410   {
01411     ZUNGQR_F77( &m, &n, &k, A, &lda, TAU, WORK, &lwork, info);
01412   }
01413 
01414 
01415   void LAPACK<int,std::complex<double> >::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
01416   {
01417     ZUNMQR_F77(CHAR_MACRO(SIDE), CHAR_MACRO(TRANS), &m, &n, &k, A, &lda, TAU, C, &ldc, WORK, &lwork, info);
01418   }
01419 
01420 
01421   void LAPACK<int,std::complex<double> >::GETRF(const int m, const int n, std::complex<double>* A, const int lda, int* IPIV, int* info) const
01422   {
01423     ZGETRF_F77(&m, &n, A, &lda, IPIV, info);
01424   }
01425 
01426 
01427   void LAPACK<int,std::complex<double> >::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
01428   {
01429     ZGETRS_F77(CHAR_MACRO(TRANS), &n, &nrhs, A, &lda, IPIV, B, &ldb, info);
01430   }
01431 
01432 
01433   void LAPACK<int,std::complex<double> >::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
01434   { ZLASCL_F77(CHAR_MACRO(TYPE), &kl, &ku, &cfrom, &cto, &m, &n, A, &lda, info); }
01435 
01436   void LAPACK<int,std::complex<double> >::GEQP3(const int m, const int n, std::complex<double>* A, const int lda, int *jpvt, std::complex<double>* TAU, std::complex<double>* WORK, const int lwork, double* RWORK, int* info ) const
01437   {
01438     ZGEQP3_F77(&m, &n, A, &lda, jpvt, TAU, WORK, &lwork, RWORK, info);
01439   }
01440 
01441   void LAPACK<int, std::complex<double> >::
01442   LASWP (const int N,
01443          std::complex<double> A[],
01444          const int LDA,
01445          const int K1,
01446          const int K2,
01447          const int IPIV[],
01448          const int INCX) const
01449   {
01450     ZLASWP_F77(&N, A, &LDA, &K1, &K2, IPIV, &INCX);
01451   }
01452 
01453   void LAPACK<int,std::complex<double> >::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
01454   {
01455     ZGBTRF_F77(&m, &n, &kl, &ku, A, &lda, IPIV, info);
01456   }
01457 
01458 
01459   void LAPACK<int,std::complex<double> >::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
01460   {
01461     ZGBTRS_F77(CHAR_MACRO(TRANS), &n, &kl, &ku, &nrhs, A, &lda, IPIV, B, &ldb, info);
01462   }
01463 
01464 
01465   void LAPACK<int,std::complex<double> >::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
01466   {
01467     ZGTTRF_F77(&n, dl, d, du, du2, IPIV, info);
01468   }
01469 
01470 
01471   void LAPACK<int,std::complex<double> >::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
01472   {
01473     ZGTTRS_F77(CHAR_MACRO(TRANS), &n, &nrhs, dl, d, du, du2, IPIV, B, &ldb, info);
01474   }
01475 
01476 
01477   void LAPACK<int,std::complex<double> >::GETRI(const int n, std::complex<double>* A, const int lda, const int* IPIV, std::complex<double>* WORK, const int lwork, int* info) const
01478   {
01479     ZGETRI_F77(&n, A, &lda, IPIV, WORK, &lwork, info);
01480   }
01481 
01482   void
01483   LAPACK<int, std::complex<double> >::LATRS (const char UPLO,
01484                                              const char TRANS,
01485                                              const char DIAG,
01486                                              const char NORMIN,
01487                                              const int N,
01488                                              std::complex<double>* A,
01489                                              const int LDA,
01490                                              std::complex<double>* X,
01491                                              double* SCALE,
01492                                              double* CNORM,
01493                                              int* INFO) const
01494   {
01495     ZLATRS_F77(CHAR_MACRO(UPLO),
01496                CHAR_MACRO(TRANS),
01497                CHAR_MACRO(DIAG),
01498                CHAR_MACRO(NORMIN),
01499                &N,
01500                A,
01501                &LDA,
01502                X,
01503                SCALE,
01504                CNORM,
01505                INFO);
01506   }
01507 
01508   void LAPACK<int,std::complex<double> >::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
01509   {
01510     ZGECON_F77(CHAR_MACRO(NORM), &n, A, &lda, &anorm, rcond, WORK, RWORK, info);
01511   }
01512 
01513 
01514   void LAPACK<int,std::complex<double> >::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
01515   {
01516     ZGBCON_F77(CHAR_MACRO(NORM), &n, &kl, &ku, A, &lda, IPIV, &anorm, rcond, WORK, RWORK, info);
01517   }
01518 
01519 
01520   double LAPACK<int,std::complex<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
01521   {
01522     return( ZLANGB_F77(CHAR_MACRO(NORM), &n, &kl, &ku, A, &lda, WORK) );
01523   }
01524 
01525 
01526   void LAPACK<int,std::complex<double> >::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
01527   {
01528     ZGESV_F77(&n, &nrhs, A, &lda, IPIV, B, &ldb, info);
01529   }
01530 
01531 
01532   void LAPACK<int,std::complex<double> >::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
01533   {
01534     ZGEEQU_F77(&m, &n, A, &lda, R, C, rowcond, colcond, amax, info);
01535   }
01536 
01537 
01538   void LAPACK<int,std::complex<double> >::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
01539   {
01540     ZGERFS_F77(CHAR_MACRO(TRANS), &n, &nrhs, A, &lda, AF, &ldaf, IPIV, B, &ldb, X, &ldx, FERR, BERR, WORK, RWORK, info);
01541   }
01542 
01543 
01544   void LAPACK<int,std::complex<double> >::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
01545   {
01546     ZGBEQU_F77(&m, &n, &kl, &ku, A, &lda, R, C, rowcond, colcond, amax, info);
01547   }
01548 
01549 
01550   void LAPACK<int,std::complex<double> >::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
01551   {
01552     ZGBRFS_F77(CHAR_MACRO(TRANS), &n, &kl, &ku, &nrhs, A, &lda, AF, &ldaf, IPIV, B, &ldb, X, &ldx, FERR, BERR, WORK, RWORK, info);
01553   }
01554 
01555 
01556   void LAPACK<int,std::complex<double> >::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
01557   {
01558     ZGESVX_F77(CHAR_MACRO(FACT), CHAR_MACRO(TRANS), &n, &nrhs, A, &lda, AF, &ldaf, IPIV, CHAR_MACRO(EQUED), R, C, B, &ldb, X, &ldx, rcond, FERR, BERR, WORK, RWORK, info);
01559   }
01560 
01561 
01562   void LAPACK<int,std::complex<double> >::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
01563   {
01564     ZGEHRD_F77(&n, &ilo, &ihi, A, &lda, TAU, WORK, &lwork, info);
01565   }
01566 
01567 
01568   void LAPACK<int,std::complex<double> >::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
01569   {
01570     ZTRTRS_F77(CHAR_MACRO(UPLO), CHAR_MACRO(TRANS), CHAR_MACRO(DIAG), &n, &nrhs, A, &lda, B, &ldb, info);
01571   }
01572 
01573 
01574   void LAPACK<int,std::complex<double> >::TRTRI(const char UPLO, const char DIAG, const int n, const std::complex<double>* A, const int lda, int* info) const
01575   {
01576     ZTRTRI_F77(CHAR_MACRO(UPLO), CHAR_MACRO(DIAG), &n, A, &lda, info);
01577   }
01578 
01579 
01580   void LAPACK<int,std::complex<double> >::STEQR(const char COMPZ, const int n, double* D, double* E, std::complex<double>* Z, const int ldz, double* WORK, int* info) const
01581   {
01582     ZSTEQR_F77(CHAR_MACRO(COMPZ), &n, D, E, Z, &ldz, WORK, info);
01583   }
01584 
01585 
01586   void LAPACK<int,std::complex<double> >::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
01587   {
01588     ZHEEV_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &n, A, &lda, W, WORK, &lwork, RWORK, info);
01589   }
01590 
01591 
01592   void LAPACK<int,std::complex<double> >::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
01593   {
01594     ZHEGV_F77(&itype, CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &n, A, &lda, B, &ldb, W, WORK, &lwork, RWORK, info);
01595   }
01596 
01597 
01598   void LAPACK<int, std::complex<double> >::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
01599   {
01600     ZHSEQR_F77(CHAR_MACRO(JOB), CHAR_MACRO(COMPZ), &n, &ilo, &ihi, H, &ldh, W, Z, &ldz, WORK, &lwork, info);
01601   }
01602 
01603 
01604   void LAPACK<int, std::complex<double> >::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
01605   {
01606     ZGEES_F77(CHAR_MACRO(JOBVS), CHAR_MACRO(SORT), ptr2func, &n, A, &lda, sdim, W, VS, &ldvs, WORK, &lwork, RWORK, BWORK, info);
01607   }
01608 
01609 
01610   void LAPACK<int, std::complex<double> >::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
01611   {
01612     int (*nullfptr)(std::complex<double>*) = NULL;
01613     std::vector< std::complex<double> > W(n);
01614     const char sort = 'N';
01615     ZGEES_F77(CHAR_MACRO(JOBVS), CHAR_MACRO(sort), nullfptr, &n, A, &lda, sdim, &W[0], VS, &ldvs, WORK, &lwork, RWORK, BWORK, info);
01616     for (int i=0; i<n; i++) {
01617       WR[i] = W[i].real();
01618       WI[i] = W[i].imag();
01619     }
01620   }
01621 
01622 
01623   void LAPACK<int, std::complex<double> >::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   {
01625     ZGEEV_F77(CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), &n, A, &lda, W, VL, &ldvl, VR, &ldvr, WORK, &lwork, RWORK, info);
01626   }
01627 
01628 
01629   void LAPACK<int, std::complex<double> >::GEEV(const char JOBVL, const char JOBVR, const int n, std::complex<double>* A, const int lda, double* WR, double* WI, std::complex<double>* VL, const int ldvl, std::complex<double>* VR, const int ldvr, std::complex<double>* WORK, const int lwork, double* RWORK, int* info) const
01630   {
01631     std::vector<std::complex<double> > w (n);
01632     std::complex<double>* w_rawPtr = (n == 0) ? NULL : &w[0];
01633     GEEV (JOBVL, JOBVR, n, A, lda, w_rawPtr, VL, ldvl, VR, ldvr, WORK, lwork, RWORK, info);
01634     if (info == 0) {
01635       // The eigenvalues are only valid on output if INFO is zero.
01636       // Otherwise, we shouldn't even write to WR or WI.
01637       for (int k = 0; k < n; ++k) {
01638   WR[k] = w[k].real ();
01639   WI[k] = w[k].imag ();
01640       }
01641     }
01642   }
01643 
01644 
01645   void LAPACK<int, std::complex<double> >::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 {
01646     ZGESVD_F77(CHAR_MACRO(JOBU), CHAR_MACRO(JOBVT), &m, &n, A, &lda, S, U, &ldu, V, &ldv, WORK, &lwork, RWORK, info);
01647   }
01648 
01649   void LAPACK<int, std::complex<double> >::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
01650   {
01651     ZGEEVX_F77(CHAR_MACRO(BALANC), CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), CHAR_MACRO(SENSE), &n, A, &lda, W, VL, &ldvl, VR, &ldvr, ilo, ihi, SCALE, abnrm, RCONDE, RCONDV, WORK, &lwork, RWORK, info);
01652   }
01653 
01654   void LAPACK<int, std::complex<double> >::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
01655   {
01656     ZGGEVX_F77(CHAR_MACRO(BALANC), CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), CHAR_MACRO(SENSE), &n, A, &lda, B, &ldb, ALPHA, BETA, VL, &ldvl, VR, &ldvr, ilo, ihi, lscale, rscale, abnrm, bbnrm, RCONDE, RCONDV, WORK, &lwork, RWORK, IWORK, BWORK, info);
01657   }
01658 
01659   void LAPACK<int, std::complex<double> >::GGEVX(const char BALANC, const char JOBVL, const char JOBVR, const char SENSE, const int n, std::complex<double>* A, const int lda, std::complex<double>* B, const int ldb, double* ALPHAR, double* ALPHAI, std::complex<double>* BETA, std::complex<double>* VL, const int ldvl, std::complex<double>* VR, const int ldvr, int* ilo, int* ihi, double* lscale, double* rscale, double* abnrm, double* bbnrm, double* RCONDE, double* RCONDV, std::complex<double>* WORK, const int lwork, double* RWORK, int* IWORK, int* BWORK, int* info) const
01660   {
01661     std::vector<std::complex<double> > w (n);
01662     std::complex<double>* w_rawPtr = (n == 0) ? NULL : &w[0];
01663     GGEVX(BALANC, JOBVL, JOBVR, SENSE, n, A, lda, B, ldb, w_rawPtr, BETA, VL, ldvl, VR, ldvr, ilo, ihi, lscale, rscale, abnrm, bbnrm, RCONDE, RCONDV, WORK, lwork, RWORK, IWORK, BWORK, info);
01664     if (info == 0) {
01665       // The eigenvalues are only valid on output if INFO is zero.
01666       // Otherwise, we shouldn't even write to WR or WI.
01667       for (int k = 0; k < n; ++k) {
01668   ALPHAR[k] = w[k].real ();
01669   ALPHAI[k] = w[k].imag ();
01670       }
01671     }
01672   }
01673 
01674   void LAPACK<int, std::complex<double> >::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
01675   {
01676     ZTREVC_F77(CHAR_MACRO(SIDE), CHAR_MACRO(HOWMNY), select, &n, T, &ldt, VL, &ldvl, VR, &ldvr, &mm, m, WORK, RWORK, info);
01677   }
01678 
01679 
01680   void LAPACK<int, std::complex<double> >::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
01681   {
01682     std::vector<int> select(1);
01683     const char whch = 'A';
01684     ZTREVC_F77(CHAR_MACRO(SIDE), CHAR_MACRO(whch), &select[0], &n, T, &ldt, VL, &ldvl, VR, &ldvr, &mm, m, WORK, RWORK, info);
01685   }
01686 
01687 
01688   void LAPACK<int, std::complex<double> >::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
01689   {
01690     ZTREXC_F77(CHAR_MACRO(COMPQ), &n, T, &ldt, Q, &ldq, &ifst, &ilst, info);
01691   }
01692 
01693 
01694   void LAPACK<int, std::complex<double> >::LARTG( const std::complex<double> f, const std::complex<double> g, double* c, std::complex<double>* s, std::complex<double>* r ) const
01695   {
01696     ZLARTG_F77(&f, &g, c, s, r);
01697   }
01698 
01699 
01700   void LAPACK<int, std::complex<double> >::LARFG( const int n, std::complex<double>* alpha, std::complex<double>* x, const int incx, std::complex<double>* tau ) const
01701   {
01702     ZLARFG_F77(&n, alpha, x, &incx, tau);
01703   }
01704 
01705 
01706   void LAPACK<int, std::complex<double> >::GEBAL(const char JOBZ, const int n, std::complex<double>* A, const int lda, int ilo, int ihi, double* scale, int* info) const
01707   {
01708     ZGEBAL_F77(CHAR_MACRO(JOBZ),&n, A, &lda, &ilo, &ihi, scale, info);
01709   }
01710 
01711 
01712   void LAPACK<int, std::complex<double> >::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
01713   {
01714     ZGEBAK_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(SIDE), &n, &ilo, &ihi, scale, &m, V, &ldv, info);
01715   }
01716 
01717 
01718 #ifdef HAVE_TEUCHOS_LAPACKLARND
01719   std::complex<double> LAPACK<int, std::complex<double> >::LARND( const int idist, int* seed ) const
01720   {
01721     return(ZLARND_F77(&idist, seed));
01722   }
01723 #endif
01724 
01725   void LAPACK<int, std::complex<double> >::LARNV( const int idist, int* seed, const int n, std::complex<double>* v ) const
01726   {
01727     ZLARNV_F77(&idist, seed, &n, v);
01728   }
01729 
01730 
01731   int LAPACK<int, std::complex<double> >::ILAENV( const int ispec, const std::string& NAME, const std::string& OPTS, const int N1, const int N2, const int N3, const int N4 ) const
01732   {
01733     unsigned int opts_length = OPTS.length();
01734     std::string temp_NAME = "z" + NAME;
01735     unsigned int name_length = temp_NAME.length();
01736     return ilaenv_wrapper(&ispec, &temp_NAME[0], name_length, &OPTS[0], opts_length, &N1, &N2, &N3, &N4);
01737   }
01738 
01739   // END INT, COMPLEX<DOUBLE> SPECIALIZATION IMPLEMENTATION //
01740 
01741 #endif // HAVE_TEUCHOS_COMPLEX
01742 
01743 } // namespace Teuchos
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines