Epetra Package Browser (Single Doxygen Collection)  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Epetra_LAPACK.h
Go to the documentation of this file.
1 /*
2 //@HEADER
3 // ************************************************************************
4 //
5 // Epetra: Linear Algebra Services Package
6 // Copyright 2011 Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
39 //
40 // ************************************************************************
41 //@HEADER
42 */
43 
44 #ifndef EPETRA_LAPACK_H
45 #define EPETRA_LAPACK_H
46 
48 
64 #include "Epetra_ConfigDefs.h"
65 #include "Epetra_Object.h"
66 
67 class EPETRA_LIB_DLL_EXPORT Epetra_LAPACK {
68 
69  public:
71 
72 
75  Epetra_LAPACK(void);
76 
77 
79 
81  Epetra_LAPACK(const Epetra_LAPACK& LAPACK);
82 
84  virtual ~Epetra_LAPACK(void);
86 
87 
89 
90 
92  void POTRF( const char UPLO, const int N, float * A, const int LDA, int * INFO) const;
94  void POTRF( const char UPLO, const int N, double * A, const int LDA, int * INFO) const;
95 
97  void POTRS( const char UPLO, const int N, const int NRHS, const float * A, const int LDA, float * X, const int LDX, int * INFO) const;
99  void POTRS( const char UPLO, const int N, const int NRHS, const double * A, const int LDA, double * X, const int LDX, int * INFO) const;
100 
102  void POTRI( const char UPLO, const int N, float * A, const int LDA, int * INFO) const;
104  void POTRI( const char UPLO, const int N, double * A, const int LDA, int * INFO) const;
105 
107  void POCON( const char UPLO, const int N, const float * A, const int LDA, const float ANORM,
108  float * RCOND, float * WORK, int * IWORK, int * INFO) const;
110  void POCON( const char UPLO, const int N, const double * A, const int LDA, const double ANORM,
111  double * RCOND, double * WORK, int * IWORK, int * INFO) const;
112 
114  void POSV( const char UPLO, const int N, const int NRHS, float * A, const int LDA, float * X, const int LDX, int * INFO) const;
116  void POSV( const char UPLO, const int N, const int NRHS, double * A, const int LDA, double * X, const int LDX, int * INFO) const;
117 
119  void POEQU(const int N, const float * A, const int LDA, float * S, float * SCOND, float * AMAX, int * INFO) const;
121  void POEQU(const int N, const double * A, const int LDA, double * S, double * SCOND, double * AMAX, int * INFO) const;
122 
124  void PORFS(const char UPLO, const int N, const int NRHS, const float * A, const int LDA, const float * AF, const int LDAF,
125  const float * B, const int LDB, float * X, const int LDX,
126  float * FERR, float * BERR, float * WORK, int * IWORK, int * INFO) const;
128  void PORFS(const char UPLO, const int N, const int NRHS, const double * A, const int LDA, const double * AF, const int LDAF,
129  const double * B, const int LDB, double * X, const int LDX,
130  double * FERR, double * BERR, double * WORK, int * IWORK, int * INFO) const;
131 
133  void POSVX(const char FACT, const char UPLO, const int N, const int NRHS, float * A, const int LDA, float * AF, const int LDAF,
134  const char EQUED, float * S, float * B, const int LDB, float * X, const int LDX, float * RCOND,
135  float * FERR, float * BERR, float * WORK, int * IWORK, int * INFO) const;
137  void POSVX(const char FACT, const char UPLO, const int N, const int NRHS, double * A, const int LDA, double * AF, const int LDAF,
138  const char EQUED, double * S, double * B, const int LDB, double * X, const int LDX, double * RCOND,
139  double * FERR, double * BERR, double * WORK, int * IWORK, int * INFO) const;
141 
143 
144 
146  void GELS( const char TRANS, const int M, const int N, const int NRHS, double* A, const int LDA,
147  double* B, const int LDB, double* WORK, const int LWORK, int * INFO) const;
149  void GETRF( const int M, const int N, float * A, const int LDA, int * IPIV, int * INFO) const;
151  void GETRF( const int M, const int N, double * A, const int LDA, int * IPIV, int * INFO) const;
152 
154  void GEQRF( const int M, const int N, float * A, const int LDA, float * TAU, float * WORK, const int lwork, int * INFO) const;
156  void GEQRF( const int M, const int N, double * A, const int LDA, double * TAU, double * WORK, const int lwork, int * INFO) const;
157 
159  void GETRS( const char TRANS, const int N, const int NRHS, const float * A, const int LDA, const int * IPIV, float * X, const int LDX, int * INFO) const;
161  void GETRS( const char TRANS, const int N, const int NRHS, const double * A, const int LDA, const int * IPIV, double * X, const int LDX, int * INFO) const;
162 
164  void GETRI( const int N, float * A, const int LDA, int * IPIV, float * WORK, const int * LWORK, int * INFO) const;
166  void GETRI( const int N, double * A, const int LDA, int * IPIV, double * WORK, const int * LWORK, int * INFO) const;
167 
169  void GECON( const char NORM, const int N, const float * A, const int LDA, const float ANORM,
170  float * RCOND, float * WORK, int * IWORK, int * INFO) const;
172  void GECON( const char NORM, const int N, const double * A, const int LDA, const double ANORM,
173  double * RCOND, double * WORK, int * IWORK, int * INFO) const;
174 
176  void GESV( const int N, const int NRHS, float * A, const int LDA, int * IPIV, float * X, const int LDX, int * INFO) const;
178  void GESV( const int N, const int NRHS, double * A, const int LDA, int * IPIV, double * X, const int LDX, int * INFO) const;
179 
181  void GEEQU(const int M, const int N, const float * A, const int LDA, float * R, float * C, float * ROWCND, float * COLCND, float * AMAX, int * INFO) const;
183  void GEEQU(const int M, const int N, const double * A, const int LDA, double * R, double * C, double * ROWCND, double * COLCND, double * AMAX, int * INFO) const;
184 
186  void GERFS(const char TRANS, const int N, const int NRHS, const float * A, const int LDA, const float * AF, const int LDAF,
187  const int * IPIV, const float * B, const int LDB, float * X, const int LDX,
188  float * FERR, float * BERR, float * WORK, int * IWORK, int * INFO) const;
190  void GERFS(const char TRANS, const int N, const int NRHS, const double * A, const int LDA, const double * AF, const int LDAF,
191  const int * IPIV, const double * B, const int LDB, double * X, const int LDX,
192  double * FERR, double * BERR, double * WORK, int * IWORK, int * INFO) const;
193 
195  void GESVX(const char FACT, const char TRANS, const int N, const int NRHS, float * A, const int LDA, float * AF, const int LDAF, int * IPIV,
196  const char EQUED, float * R, float * C, float * B, const int LDB, float * X, const int LDX, float * RCOND,
197  float * FERR, float * BERR, float * WORK, int * IWORK, int * INFO) const;
199  void GESVX(const char FACT, const char TRANS, const int N, const int NRHS, double * A, const int LDA, double * AF, const int LDAF, int * IPIV,
200  const char EQUED, double * R, double * C, double * B, const int LDB, double * X, const int LDX, double * RCOND,
201  double * FERR, double * BERR, double * WORK, int * IWORK, int * INFO) const;
202 
203 
205  void GEHRD(const int N, const int ILO, const int IHI, float * A, const int LDA, float * TAU, float * WORK, const int LWORK, int * INFO) const;
207  void GEHRD(const int N, const int ILO, const int IHI, double * A, const int LDA, double * TAU, double * WORK, const int LWORK, int * INFO) const;
209 
211 
212  void HSEQR( const char JOB, const char COMPZ, const int N, const int ILO, const int IHI, float * H, const int LDH, float * WR, float * WI,
214  float * Z, const int LDZ, float * WORK, const int LWORK, int * INFO) const;
216  void HSEQR( const char JOB, const char COMPZ, const int N, const int ILO, const int IHI, double * H, const int LDH, double * WR, double * WI,
217  double * Z, const int LDZ, double * WORK, const int LWORK, int * INFO) const;
219 
221 
222  void ORGQR( const int M, const int N, const int K, float * A, const int LDA, float * TAU, float * WORK, const int LWORK, int * INFO) const;
225  void ORGQR( const int M, const int N, const int K, double * A, const int LDA, double * TAU, double * WORK, const int LWORK, int * INFO) const;
226 
228  void ORGHR( const int N, const int ILO, const int IHI, float * A, const int LDA, float * TAU, float * WORK, const int LWORK, int * INFO) const;
230  void ORGHR( const int N, const int ILO, const int IHI, double * A, const int LDA, double * TAU, double * WORK, const int LWORK, int * INFO) const;
231 
233  void ORMHR( const char SIDE, const char TRANS, const int M, const int N, const int ILO, const int IHI, const float * A, const int LDA,
234  const float * TAU, float * C,
235  const int LDC, float * WORK, const int LWORK, int * INFO) const;
237  void ORMHR( const char SIDE, const char TRANS, const int M, const int N, const int ILO, const int IHI, const double * A, const int LDA,
238  const double * TAU, double * C,
239  const int LDC, double * WORK, const int LWORK, int * INFO) const;
241  void LARFT( const char DIRECT, const char STOREV, const int N, const int K, double * V, const int LDV, double * TAU, double * T, const int LDT) const;
243  void LARFT( const char DIRECT, const char STOREV, const int N, const int K, float * V, const int LDV, float * TAU, float * T, const int LDT) const;
245 
247 
248 
250 
252  void TREVC( const char SIDE, const char HOWMNY, int * SELECT, const int N, const float * T, const int LDT, float *VL, const int LDVL,
253  float * VR, const int LDVR, const int MM, int * M, float * WORK, int * INFO) const;
255 
257  void TREVC( const char SIDE, const char HOWMNY, int * SELECT, const int N, const double * T, const int LDT, double *VL, const int LDVL,
258  double * VR, const int LDVR, const int MM, int *M, double * WORK, int * INFO) const;
259 
261  void TREXC( const char COMPQ, const int N, float * T, const int LDT, float * Q, const int LDQ, int IFST, int ILST,
262  float * WORK, int * INFO) const;
264  void TREXC( const char COMPQ, const int N, double * T, const int LDT, double * Q, const int LDQ, int IFST, int ILST,
265  double * WORK, int * INFO) const;
267 
269 
270 
272  void GESVD( const char JOBU, const char JOBVT, const int M, const int N, float * A, const int LDA, float * S, float * U,
273  const int LDU, float * VT, const int LDVT, float * WORK, const int * LWORK, int * INFO) const;
275  void GESVD( const char JOBU, const char JOBVT, const int M, const int N, double * A, const int LDA, double * S, double * U,
276  const int LDU, double * VT, const int LDVT, double * WORK, const int * LWORK, int * INFO) const;
277 
279  void GGSVD(const char JOBU, const char JOBV, const char JOBQ, const int M, const int N, const int P, int * K, int * L, double* A, const int LDA, double* B, const int LDB,
280  double* ALPHA, double* BETA, double* U, const int LDU, double* V, const int LDV, double* Q, const int LDQ, double* WORK,
281  #ifdef HAVE_EPETRA_LAPACK_GSSVD3
282  const int LWORK,
283  #endif
284  int* IWORK, int* INFO) const;
286  void GGSVD(const char JOBU, const char JOBV, const char JOBQ, const int M, const int N, const int P, int * K, int * L, float* A, const int LDA, float* B, const int LDB,
287  float* ALPHA, float* BETA, float* U, const int LDU, float* V, const int LDV, float* Q, const int LDQ, float* WORK,
288  #ifdef HAVE_EPETRA_LAPACK_GSSVD3
289  const int LWORK,
290  #endif
291  int* IWORK, int* INFO) const;
293 
295 
296  void GEEV(const char JOBVL, const char JOBVR, const int N, double* A, const int LDA, double* WR, double* WI,
298  double* VL, const int LDVL, double* VR, const int LDVR, double* WORK, const int LWORK, int* INFO) const;
300  void GEEV(const char JOBVL, const char JOBVR, const int N, float* A, const int LDA, float* WR, float* WI,
301  float* VL, const int LDVL, float* VR, const int LDVR, float* WORK, const int LWORK, int* INFO) const;
302 
304  void SPEV(const char JOBZ, const char UPLO, const int N, double* AP, double* W, double* Z, int LDZ, double* WORK, int* INFO) const;
306  void SPEV(const char JOBZ, const char UPLO, const int N, float* AP, float* W, float* Z, int LDZ, float* WORK, int* INFO) const;
307 
309  void SPGV(const int ITYPE, const char JOBZ, const char UPLO, const int N, double* AP, double* BP, double* W, double* Z, const int LDZ, double* WORK, int* INFO) const;
311  void SPGV(const int ITYPE, const char JOBZ, const char UPLO, const int N, float* AP, float* BP, float* W, float* Z, const int LDZ, float* WORK, int* INFO) const;
312 
314  void SYEV(const char JOBZ, const char UPLO, const int N, double* A, const int LDA, double* W, double* WORK, const int LWORK, int* INFO) const;
316  void SYEV(const char JOBZ, const char UPLO, const int N, float* A, const int LDA, float* W, float* WORK, const int LWORK, int* INFO) const;
317 
319  void SYEVD(const char JOBZ, const char UPLO, const int N, double* A, const int LDA, double* W,
320  double* WORK, const int LWORK, int* IWORK, const int LIWORK, int* INFO) const;
322  void SYEVD(const char JOBZ, const char UPLO, const int N, float* A, const int LDA, float* W,
323  float* WORK, const int LWORK, int* IWORK, const int LIWORK, int* INFO) const;
324 
326  void SYEVX(const char JOBZ, const char RANGE, const char UPLO, const int N, double* A, const int LDA,
327  const double* VL, const double* VU, const int* IL, const int* IU,
328  const double ABSTOL, int * M, double* W, double* Z, const int LDZ, double* WORK,
329  const int LWORK, int* IWORK, int* IFAIL,
330  int* INFO) const;
332  void SYEVX(const char JOBZ, const char RANGE, const char UPLO, const int N, float* A, const int LDA,
333  const float* VL, const float* VU, const int* IL, const int* IU,
334  const float ABSTOL, int * M, float* W, float* Z, const int LDZ, float* WORK,
335  const int LWORK, int* IWORK, int* IFAIL,
336  int* INFO) const;
337 
339  void SYGV(const int ITYPE, const char JOBZ, const char UPLO, const int N, double* A, const int LDA, double* B,
340  const int LDB, double* W, double* WORK, const int LWORK, int* INFO) const;
342  void SYGV(const int ITYPE, const char JOBZ, const char UPLO, const int N, float* A, const int LDA, float* B,
343  const int LDB, float* W, float* WORK, const int LWORK, int* INFO) const;
344 
346  void SYGVX(const int ITYPE, const char JOBZ, const char RANGE, const char UPLO, const int N,
347  double* A, const int LDA, double* B, const int LDB, const double* VL, const double* VU,
348  const int* IL, const int* IU, const double ABSTOL, int* M, double* W, double* Z,
349  const int LDZ, double* WORK, const int LWORK, int* IWORK,
350  int* IFAIL, int* INFO) const;
352  void SYGVX(const int ITYPE, const char JOBZ, const char RANGE, const char UPLO, const int N,
353  float* A, const int LDA, float* B, const int LDB, const float* VL, const float* VU,
354  const int* IL, const int* IU, const float ABSTOL, int* M, float* W, float* Z,
355  const int LDZ, float* WORK, const int LWORK, int* IWORK,
356  int* IFAIL, int* INFO) const;
357 
359  void SYEVR(const char JOBZ, const char RANGE, const char UPLO, const int N, double* A, const int LDA, const double* VL, const double* VU, const int *IL, const int *IU,
360  const double ABSTOL, int* M, double* W, double* Z, const int LDZ, int* ISUPPZ, double* WORK, const int LWORK, int* IWORK,
361  const int LIWORK, int* INFO) const;
363  void SYEVR(const char JOBZ, const char RANGE, const char UPLO, const int N, float* A, const int LDA,
364  const float* VL, const float* VU, const int *IL, const int *IU,
365  const float ABSTOL, int* M, float* W, float* Z, const int LDZ, int* ISUPPZ,
366  float* WORK, const int LWORK, int* IWORK,
367  const int LIWORK, int* INFO) const;
368 
370  void GEEVX(const char BALANC, const char JOBVL, const char JOBVR, const char SENSE, const int N, double* A, const int LDA, double* WR, double* WI, double* VL,
371  const int LDVL, double* VR, const int LDVR, int* ILO, int* IHI, double* SCALE, double* ABNRM, double* RCONDE,
372  double* RCONDV, double* WORK, const int LWORK, int* IWORK, int* INFO) const;
374  void GEEVX(const char BALANC, const char JOBVL, const char JOBVR, const char SENSE, const int N, float* A, const int LDA, float* WR, float* WI, float* VL,
375  const int LDVL, float* VR, const int LDVR, int* ILO, int* IHI, float* SCALE, float* ABNRM, float* RCONDE,
376  float* RCONDV, float* WORK, const int LWORK, int* IWORK, int* INFO) const;
377 
379  void GESDD(const char JOBZ, const int M, const int N, double* A, const int LDA, double* S, double* U, const int LDU, double* VT, const int LDVT, double* WORK,
380  const int LWORK, int* IWORK, int* INFO) const;
382  void GESDD(const char JOBZ, const int M, const int N, float* A, const int LDA, float* S, float* U, const int LDU, float* VT, const int LDVT, float* WORK,
383  const int LWORK, int* IWORK, int* INFO) const;
385 
386  void GGEV(const char JOBVL, const char JOBVR, const int N, double* A, const int LDA, double* B, const int LDB, double* ALPHAR, double* ALPHAI,
387  double* BETA, double* VL, const int LDVL, double* VR, const int LDVR, double* WORK, const int LWORK, int* INFO) const;
389  void GGEV(const char JOBVL, const char JOBVR, const int N, float* A, const int LDA, float* B, const int LDB, float* ALPHAR, float* ALPHAI,
390  float* BETA, float* VL, const int LDVL, float* VR, const int LDVR, float* WORK, const int LWORK, int* INFO) const;
391 
393 
395 
396  void GGLSE(const int M, const int N, const int P, double* A, const int LDA, double* B, const int LDB,
398  double* C, double* D, double* X, double* WORK, const int LWORK, int* INFO) const;
400  void GGLSE(const int M, const int N, const int P, float* A, const int LDA, float* B, const int LDB,
401  float* C, float* D, float* X, float* WORK, const int LWORK, int* INFO) const;
403 
405 
406  void LAMCH ( const char CMACH, float & T) const;
409  void LAMCH ( const char CMACH, double & T) const;
411 
413 
415 
416  void TRTRS(const char UPLO, const char TRANS, const char DIAG, const int N, const int NRHS, const float *A,
418  const int LDA, float *B, const int LDB, int *INFO) const;
420  void TRTRS(const char UPLO, const char TRANS, const char DIAG, const int N, const int NRHS, const double *A,
421  const int LDA, double *B, const int LDB, int *INFO) const;
422 };
423 
424 // Epetra_LAPACK constructor
426 // Epetra_LAPACK constructor
427 inline Epetra_LAPACK::Epetra_LAPACK(const Epetra_LAPACK& LAPACK){(void)LAPACK;}
428 // Epetra_LAPACK destructor
430 
431 #endif /* EPETRA_LAPACK_H */
Epetra_LAPACK: The Epetra LAPACK Wrapper Class.
Definition: Epetra_LAPACK.h:67
virtual ~Epetra_LAPACK(void)
Epetra_LAPACK Destructor.
Epetra_LAPACK(void)
Epetra_LAPACK Constructor.