Epetra Package Browser (Single Doxygen Collection)  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Epetra_LAPACK.cpp
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 #ifdef CHAR_MACRO
44 #undef CHAR_MACRO
45 #endif
46 #if defined (INTEL_CXML)
47 #define CHAR_MACRO(char_var) &char_var, 1
48 #else
49 #define CHAR_MACRO(char_var) &char_var
50 #endif
51 
52 #include "Epetra_LAPACK.h"
53 #include "Epetra_LAPACK_wrappers.h"
54 
55 
56 // Symmetric positive definite linear systems
57 
58 //=============================================================================
59 void Epetra_LAPACK::POTRF( const char UPLO, const int N, float * A, const int LDA, int * INFO) const {
60  SPOTRF_F77(CHAR_MACRO(UPLO), &N, A, &LDA, INFO);
61 }
62 //=============================================================================
63 void Epetra_LAPACK::POTRF( const char UPLO, const int N, double * A, const int LDA, int * INFO) const {
64  DPOTRF_F77(CHAR_MACRO(UPLO), &N, A, &LDA, INFO);
65 }
66 //=============================================================================
67 void Epetra_LAPACK::POTRS( const char UPLO, const int N, const int NRHS, const float * A, const int LDA,
68  float * X, const int LDX, int * INFO) const {
69  SPOTRS_F77(CHAR_MACRO(UPLO), &N, &NRHS, A, &LDA, X, &LDX, INFO);
70 }
71 //=============================================================================
72 void Epetra_LAPACK::POTRS( const char UPLO, const int N, const int NRHS, const double * A, const int LDA,
73  double * X, const int LDX, int * INFO) const {
74  DPOTRS_F77(CHAR_MACRO(UPLO), &N, &NRHS, A, &LDA, X, &LDX, INFO);
75 }
76 //=============================================================================
77 void Epetra_LAPACK::POTRI( const char UPLO, const int N, float * A, const int LDA, int * INFO) const {
78  SPOTRI_F77(CHAR_MACRO(UPLO), &N, A, &LDA, INFO);
79 }
80 //=============================================================================
81 void Epetra_LAPACK::POTRI( const char UPLO, const int N, double * A, const int LDA, int * INFO) const {
82  DPOTRI_F77(CHAR_MACRO(UPLO), &N, A, &LDA, INFO);
83 }
84 //=============================================================================
85 void Epetra_LAPACK::POCON( const char UPLO, const int N, const float * A, const int LDA, const float ANORM,
86  float * RCOND, float * WORK, int * IWORK,
87  int * INFO) const {
88  SPOCON_F77(CHAR_MACRO(UPLO), &N, A, &LDA, &ANORM, RCOND, WORK, IWORK, INFO);
89 }
90 //=============================================================================
91 void Epetra_LAPACK::POCON( const char UPLO, const int N, const double * A, const int LDA, const double ANORM,
92  double * RCOND, double * WORK, int * IWORK,
93  int * INFO) const {
94  DPOCON_F77(CHAR_MACRO(UPLO), &N, A, &LDA, &ANORM, RCOND, WORK, IWORK, INFO);
95 }
96 //=============================================================================
97 void Epetra_LAPACK::POSV( const char UPLO, const int N, const int NRHS, float * A, const int LDA,
98  float * X, const int LDX, int * INFO) const {
99  SPOSV_F77(CHAR_MACRO(UPLO), &N, &NRHS, A, &LDA, X, &LDX, INFO);
100 }
101 //=============================================================================
102 void Epetra_LAPACK::POSV( const char UPLO, const int N, const int NRHS, double * A, const int LDA,
103  double * X, const int LDX, int * INFO) const {
104  DPOSV_F77(CHAR_MACRO(UPLO), &N, &NRHS, A, &LDA, X, &LDX, INFO);
105 }
106 //=============================================================================
107 void Epetra_LAPACK::POEQU(const int N, const float * A, const int LDA, float * S, float * SCOND,
108  float * AMAX, int * INFO) const {
109  SPOEQU_F77(&N, A, &LDA, S, SCOND, AMAX, INFO);
110 }
111 //=============================================================================
112 void Epetra_LAPACK::POEQU(const int N, const double * A, const int LDA, double * S, double * SCOND,
113  double * AMAX, int * INFO) const {
114  DPOEQU_F77(&N, A, &LDA, S, SCOND, AMAX, INFO);
115 }
116 //=============================================================================
117 void Epetra_LAPACK::PORFS(const char UPLO, const int N, const int NRHS, const float * A, const int LDA, const float * AF, const int LDAF,
118  const float * B, const int LDB, float * X, const int LDX,
119  float * FERR, float * BERR, float * WORK, int * IWORK, int * INFO) const {
120  SPORFS_F77(CHAR_MACRO(UPLO), &N, &NRHS, A, &LDA, AF, &LDAF, B, &LDB, X, &LDX,
121  FERR, BERR, WORK, IWORK, INFO);
122 }
123 //=============================================================================
124 void Epetra_LAPACK::PORFS(const char UPLO, const int N, const int NRHS, const double * A, const int LDA, const double * AF, const int LDAF,
125  const double * B, const int LDB, double * X, const int LDX,
126  double * FERR, double * BERR, double * WORK, int * IWORK, int * INFO) const {
127  DPORFS_F77( CHAR_MACRO(UPLO), &N, &NRHS, A, &LDA, AF, &LDAF,B, &LDB, X, &LDX,
128  FERR, BERR, WORK, IWORK, INFO);
129 }
130 //=============================================================================
131 void Epetra_LAPACK::POSVX(const char FACT, const char UPLO, const int N, const int NRHS, float * A, const int LDA, float * AF, const int LDAF,
132  const char EQUED, float * S, float * B, const int LDB, float * X, const int LDX, float * RCOND,
133  float * FERR, float * BERR, float * WORK, int * IWORK, int * INFO) const {
134  SPOSVX_F77(CHAR_MACRO(FACT), CHAR_MACRO(UPLO), &N, &NRHS, A, &LDA, AF, &LDAF, CHAR_MACRO(EQUED), S, B, &LDB, X, &LDX,
135  RCOND, FERR, BERR, WORK, IWORK, INFO);
136 }
137 //=============================================================================
138 void Epetra_LAPACK::POSVX(const char FACT, const char UPLO, const int N, const int NRHS, double * A, const int LDA, double * AF, const int LDAF,
139  const char EQUED, double * S, double * B, const int LDB, double * X, const int LDX, double * RCOND,
140  double * FERR, double * BERR, double * WORK, int * IWORK, int * INFO) const {
141  DPOSVX_F77(CHAR_MACRO(FACT), CHAR_MACRO(UPLO), &N, &NRHS, A, &LDA, AF, &LDAF, CHAR_MACRO(EQUED), S, B, &LDB, X, &LDX, RCOND,
142  FERR, BERR, WORK, IWORK, INFO);
143 }
144 
145 // General linear systems
146 //=============================================================================
147 void Epetra_LAPACK::GELS( const char TRANS, const int M, const int N, const int NRHS,
148  double* A, const int LDA, double* B, const int LDB, double* WORK, const int LWORK,
149  int * INFO) const {
150  DGELS_F77 (CHAR_MACRO(TRANS), &M, &N, &NRHS, A, &LDA, B, &LDB, WORK, &LWORK, INFO );
151 
152 }
153 //=============================================================================
154 void Epetra_LAPACK::GETRF( const int M, const int N, float * A, const int LDA, int * IPIV, int * INFO) const {
155  SGETRF_F77(&M, &N, A, &LDA, IPIV, INFO);
156 }
157 //=============================================================================
158 void Epetra_LAPACK::GETRF( const int M, const int N, double * A, const int LDA, int * IPIV, int * INFO) const {
159  DGETRF_F77(&M, &N, A, &LDA, IPIV, INFO);
160 }
161 //=============================================================================
162 void Epetra_LAPACK::GEQRF( const int M, const int N, float * A, const int LDA, float * TAU, float * WORK, const int LWORK, int * INFO) const {
163  SGEQRF_F77(&M, &N, A, &LDA, TAU, WORK, &LWORK, INFO);
164 }
165 //=============================================================================
166 void Epetra_LAPACK::GEQRF( const int M, const int N, double * A, const int LDA, double * TAU, double * WORK, const int LWORK, int * INFO) const {
167  DGEQRF_F77(&M, &N, A, &LDA, TAU, WORK, &LWORK, INFO);
168 }
169 //=============================================================================
170 void Epetra_LAPACK::GESVD(const char JOBU, const char JOBVT, const int M, const int N, float * A,
171  const int LDA, float * S, float * U,
172  const int LDU, float * VT, const int LDVT, float * WORK,
173  const int * LWORK, int * INFO) const {
174  SGESVD_F77(CHAR_MACRO(JOBU), CHAR_MACRO(JOBVT), &M, &N, A, &LDA, S, U, &LDU,
175  VT, &LDVT, WORK, LWORK, INFO);
176 }
177 //=============================================================================
178 void Epetra_LAPACK::GESVD(const char JOBU, const char JOBVT, const int M, const int N, double * A,
179  const int LDA, double * S, double * U,
180  const int LDU, double * VT, const int LDVT, double * WORK,
181  const int * LWORK, int * INFO) const {
182  DGESVD_F77(CHAR_MACRO(JOBU), CHAR_MACRO(JOBVT), &M, &N, A, &LDA, S, U, &LDU,
183  VT, &LDVT, WORK, LWORK, INFO);
184 }
185 //=============================================================================
186 void Epetra_LAPACK::GETRS( const char TRANS, const int N, const int NRHS, const float * A, const int LDA,
187  const int * IPIV, float * X, const int LDX, int * INFO) const {
188  SGETRS_F77(CHAR_MACRO(TRANS), &N, &NRHS, A, &LDA, IPIV, X, &LDX, INFO);
189 }
190 //=============================================================================
191 void Epetra_LAPACK::GETRS( const char TRANS, const int N, const int NRHS, const double * A, const int LDA,
192  const int * IPIV, double * X, const int LDX, int * INFO) const {
193  DGETRS_F77(CHAR_MACRO(TRANS), &N, &NRHS, A, &LDA, IPIV, X, &LDX, INFO);
194 }
195 //=============================================================================
196 void Epetra_LAPACK::GETRI( const int N, float * A, const int LDA, int * IPIV,
197  float * WORK, const int * LWORK, int * INFO) const {
198  SGETRI_F77(&N, A, &LDA, IPIV, WORK, LWORK, INFO);
199 }
200 //=============================================================================
201 void Epetra_LAPACK::GETRI( const int N, double * A, const int LDA, int * IPIV,
202  double * WORK, const int * LWORK, int * INFO) const {
203  DGETRI_F77(&N, A, &LDA, IPIV, WORK, LWORK, INFO);
204 }
205 //=============================================================================
206 void Epetra_LAPACK::GECON( const char NORM, const int N, const float * A, const int LDA, const float ANORM,
207  float * RCOND, float * WORK, int * IWORK,
208  int * INFO) const {
209  SGECON_F77(CHAR_MACRO(NORM), &N, A, &LDA, &ANORM, RCOND, WORK, IWORK, INFO);
210 }
211 //=============================================================================
212 void Epetra_LAPACK::GECON( const char NORM, const int N, const double * A, const int LDA, const double ANORM,
213  double * RCOND, double * WORK, int * IWORK,
214  int * INFO) const {
215  DGECON_F77(CHAR_MACRO(NORM), &N, A, &LDA, &ANORM, RCOND, WORK, IWORK, INFO);
216 }
217 //=============================================================================
218 void Epetra_LAPACK::GESV( const int N, const int NRHS, float * A, const int LDA, int * IPIV,
219  float * X, const int LDX, int * INFO) const {
220  SGESV_F77(&N, &NRHS, A, &LDA, IPIV, X, &LDX, INFO);
221 }
222 //=============================================================================
223 void Epetra_LAPACK::GESV( const int N, const int NRHS, double * A, const int LDA, int * IPIV,
224  double * X, const int LDX, int * INFO) const {
225  DGESV_F77(&N, &NRHS, A, &LDA, IPIV, X, &LDX, INFO);
226 }
227 //=============================================================================
228 void Epetra_LAPACK::GEEQU(const int M, const int N, const float * A, const int LDA, float * R, float * C,
229  float * ROWCND, float * COLCND, float * AMAX, int * INFO) const {
230  SGEEQU_F77(&M, &N, A, &LDA, R, C, ROWCND, COLCND, AMAX, INFO);
231 }
232 //=============================================================================
233 void Epetra_LAPACK::GEEQU(const int M, const int N, const double * A, const int LDA, double * R, double * C,
234  double * ROWCND, double * COLCND, double * AMAX, int * INFO) const {
235  DGEEQU_F77(&M, &N, A, &LDA, R, C, ROWCND, COLCND, AMAX, INFO);
236 }
237 //=============================================================================
238 void Epetra_LAPACK::GERFS(const char TRANS, const int N, const int NRHS, const float * A, const int LDA, const float * AF, const int LDAF,
239  const int * IPIV, const float * B, const int LDB, float * X, const int LDX,
240  float * FERR, float * BERR, float * WORK, int * IWORK, int * INFO) const {
241  SGERFS_F77(CHAR_MACRO(TRANS), &N, &NRHS, A, &LDA, AF, &LDAF, IPIV, B, &LDB, X, &LDX,
242  FERR, BERR, WORK, IWORK, INFO);
243 }
244 //=============================================================================
245 void Epetra_LAPACK::GERFS(const char TRANS, const int N, const int NRHS, const double * A, const int LDA, const double * AF, const int LDAF,
246  const int * IPIV, const double * B, const int LDB, double * X, const int LDX,
247  double * FERR, double * BERR, double * WORK, int * IWORK, int * INFO) const {
248  DGERFS_F77( CHAR_MACRO(TRANS), &N, &NRHS, A, &LDA, AF, &LDAF, IPIV, B, &LDB, X, &LDX,
249  FERR, BERR, WORK, IWORK, INFO);
250 }
251 //=============================================================================
252 void Epetra_LAPACK::GESVX(const char FACT, const char TRANS, const int N, const int NRHS, float * A, const int LDA, float * AF, const int LDAF,
253  int * IPIV, const char EQUED, float * R, float * C, float * B, const int LDB, float * X, const int LDX, float * RCOND,
254  float * FERR, float * BERR, float * WORK, int * IWORK, int * INFO) const {
255  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,
256  FERR, BERR, WORK, IWORK, INFO);
257 }
258 //=============================================================================
259 void Epetra_LAPACK::GESVX(const char FACT, const char TRANS, const int N, const int NRHS, double * A, const int LDA, double * AF, const int LDAF,
260  int * IPIV, const char EQUED, double * R, double * C, double * B, const int LDB, double * X, const int LDX, double * RCOND,
261  double * FERR, double * BERR, double * WORK, int * IWORK, int * INFO) const {
262  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,
263  FERR, BERR, WORK, IWORK, INFO);
264 }
265 
266 
267 
268 
269 //=============================================================================
270 void Epetra_LAPACK::GEHRD(const int N, const int ILO, const int IHI, float * A, const int LDA, float * TAU,
271  float * WORK, const int LWORK, int * INFO) const {
272  SGEHRD_F77(&N, &ILO, &IHI, A, &LDA, TAU, WORK, &LWORK, INFO);
273 }
274 //=============================================================================
275 void Epetra_LAPACK::GEHRD(const int N, const int ILO, const int IHI, double * A, const int LDA, double * TAU,
276  double * WORK, const int LWORK, int * INFO) const {
277  DGEHRD_F77(&N, &ILO, &IHI, A, &LDA, TAU, WORK, &LWORK, INFO);
278 }
279 //=============================================================================
280 void Epetra_LAPACK::HSEQR( const char JOB, const char COMPZ, const int N, const int ILO, const int IHI, float * H, const int LDH,
281  float * WR, float * WI, float * Z, const int LDZ, float * WORK, const int LWORK,
282  int * INFO) const {
283  SHSEQR_F77(CHAR_MACRO(JOB), CHAR_MACRO(COMPZ), &N, &ILO, &IHI, H, &LDH, WR, WI, Z, &LDZ, WORK, &LWORK, INFO);
284 }
285 //=============================================================================
286 void Epetra_LAPACK::HSEQR( const char JOB, const char COMPZ, const int N, const int ILO, const int IHI, double * H, const int LDH,
287  double * WR, double * WI, double * Z, const int LDZ, double * WORK, const int LWORK,
288  int * INFO) const {
289  DHSEQR_F77(CHAR_MACRO(JOB), CHAR_MACRO(COMPZ), &N, &ILO, &IHI, H, &LDH, WR, WI, Z, &LDZ, WORK, &LWORK, INFO);
290 }
291 //=============================================================================
292 void Epetra_LAPACK::ORGQR( const int M, const int N, const int K, float * A, const int LDA, float * TAU,
293  float * WORK, const int LWORK, int * INFO) const {
294  SORGQR_F77( &M, &N, &K, A, &LDA, TAU, WORK, &LWORK, INFO);
295 }
296 //=============================================================================
297 void Epetra_LAPACK::ORGQR( const int M, const int N, const int K, double * A, const int LDA, double * TAU,
298  double * WORK, const int LWORK, int * INFO) const {
299  DORGQR_F77( &M, &N, &K, A, &LDA, TAU, WORK, &LWORK, INFO);
300 }
301 //=============================================================================
302 void Epetra_LAPACK::ORGHR( const int N, const int ILO, const int IHI, float * A, const int LDA, float * TAU,
303  float * WORK, const int LWORK, int * INFO) const {
304  SORGHR_F77( &N, &ILO, &IHI, A, &LDA, TAU, WORK, &LWORK, INFO);
305 }
306 //=============================================================================
307 void Epetra_LAPACK::ORGHR( const int N, const int ILO, const int IHI, double * A, const int LDA, double * TAU,
308  double * WORK, const int LWORK, int * INFO) const {
309  DORGHR_F77( &N, &ILO, &IHI, A, &LDA, TAU, WORK, &LWORK, INFO);
310 }
311 //=============================================================================
312 void Epetra_LAPACK::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,
313  const float * TAU, float * C, const int LDC, float * WORK, const int LWORK, int * INFO) const {
314  SORMHR_F77(CHAR_MACRO(SIDE), CHAR_MACRO(TRANS), &M, &N, &ILO, &IHI, A, &LDA, TAU, C, &LDC, WORK, &LWORK, INFO);
315 }
316 //=============================================================================
317 void Epetra_LAPACK::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,
318  const double * TAU, double * C, const int LDC, double * WORK, const int LWORK, int * INFO) const {
319  DORMHR_F77(CHAR_MACRO(SIDE), CHAR_MACRO(TRANS), &M, &N, &ILO, &IHI, A, &LDA, TAU, C, &LDC, WORK, &LWORK, INFO);
320 }
321 //=============================================================================
322 void Epetra_LAPACK::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 {
323  SLARFT_F77(CHAR_MACRO(DIRECT), CHAR_MACRO(STOREV), &N, &K, V, &LDV, TAU, T, &LDT);
324 }
325 //=============================================================================
326 void Epetra_LAPACK::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 {
327  DLARFT_F77(CHAR_MACRO(DIRECT), CHAR_MACRO(STOREV), &N, &K, V, &LDV, TAU, T, &LDT);
328 }
329 //=============================================================================
330 void Epetra_LAPACK::TREVC( const char SIDE, const char HOWMNY, int * SELECT, const int N, const float * T, const int LDT, float *VL, const int LDVL,
331  float * VR, const int LDVR, const int MM, int * M, float * WORK, int * INFO) const {
332 
333  if (HOWMNY=='S') *INFO = -3; // We do not support 'S' since it requires a logical array (yuck!)
334 
335  else STREVC_F77(CHAR_MACRO(SIDE), CHAR_MACRO(HOWMNY), SELECT, &N, T, &LDT, VL, &LDVL, VR, &LDVR, &MM, M, WORK, INFO);
336 }
337 //=============================================================================
338 void Epetra_LAPACK::TREVC( const char SIDE, const char HOWMNY, int * SELECT, const int N, const double * T, const int LDT, double *VL, const int LDVL,
339  double * VR, const int LDVR, const int MM, int * M, double * WORK, int * INFO) const {
340 
341  if (HOWMNY=='S') *INFO = -3; // We do not support 'S' since it requires a logical array (yuck!)
342 
343  else DTREVC_F77(CHAR_MACRO(SIDE), CHAR_MACRO(HOWMNY), SELECT, &N, T, &LDT, VL, &LDVL, VR, &LDVR, &MM, M, WORK, INFO);
344 }
345 //=============================================================================
346 void Epetra_LAPACK::TREXC( const char COMPQ, const int N, float * T, const int LDT, float * Q, const int LDQ, int IFST, int ILST,
347  float * WORK, int * INFO) const {
348  STREXC_F77( CHAR_MACRO(COMPQ), &N, T, &LDT, Q, &LDQ, &IFST, &ILST, WORK, INFO);
349 }
350 //=============================================================================
351 void Epetra_LAPACK::TREXC( const char COMPQ, const int N, double * T, const int LDT, double * Q, const int LDQ, int IFST, int ILST,
352  double * WORK, int * INFO) const {
353  DTREXC_F77( CHAR_MACRO(COMPQ), &N, T, &LDT, Q, &LDQ, &IFST, &ILST, WORK, INFO);
354 }
355 //=============================================================================
356 void Epetra_LAPACK::LAMCH( const char CMACH, float & T) const {
357  T = SLAMCH_F77( CHAR_MACRO(CMACH));
358 }
359 //=============================================================================
360 void Epetra_LAPACK::LAMCH( const char CMACH, double & T) const {
361  T = DLAMCH_F77( CHAR_MACRO(CMACH));
362 }
363 //=============================================================================
364 void Epetra_LAPACK::GGSVD(const char JOBU, const char JOBV, const char JOBQ, const int M, const int N, const int P, int * K, int * L,
365  double* A, const int LDA, double* B, const int LDB,
366  double* ALPHA, double* BETA, double* U, const int LDU, double* V, const int LDV, double* Q, const int LDQ, double* WORK,
367  #ifdef HAVE_EPETRA_LAPACK_GSSVD3
368  const int LWORK,
369  #endif
370  int* IWORK, int* INFO) const {
371  DGGSVD_F77(CHAR_MACRO(JOBU), CHAR_MACRO(JOBV), CHAR_MACRO(JOBQ), &M, &N, &P, K, L, A, &LDA, B, &LDB,
372  ALPHA, BETA, U, &LDU, V, &LDV, Q, &LDQ, WORK,
373  #ifdef HAVE_EPETRA_LAPACK_GSSVD3
374  &LWORK,
375  #endif
376  IWORK, INFO);
377 }
378 
379 //=============================================================================
380 void Epetra_LAPACK::GGSVD(const char JOBU, const char JOBV, const char JOBQ, const int M, const int N, const int P, int * K, int * L,
381  float* A, const int LDA, float* B, const int LDB,
382  float* ALPHA, float* BETA, float* U, const int LDU, float* V, const int LDV, float* Q, const int LDQ, float* WORK,
383  #ifdef HAVE_EPETRA_LAPACK_GSSVD3
384  const int LWORK,
385  #endif
386  int* IWORK, int* INFO) const {
387  SGGSVD_F77(CHAR_MACRO(JOBU), CHAR_MACRO(JOBV), CHAR_MACRO(JOBQ), &M, &N, &P, K, L, A, &LDA, B, &LDB,
388  ALPHA, BETA, U, &LDU, V, &LDV, Q, &LDQ, WORK,
389  #ifdef HAVE_EPETRA_LAPACK_GSSVD3
390  &LWORK,
391  #endif
392  IWORK, INFO);
393 }
394 
395 //=============================================================================
396 void Epetra_LAPACK::GEEV(const char JOBVL, const char JOBVR, const int N, double* A, const int LDA, double* WR, double* WI,
397  double* VL, const int LDVL, double* VR, const int LDVR, double* WORK, const int LWORK, int* INFO) const {
398 
399  DGEEV_F77(CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), &N, A, &LDA, WR, WI, VL, &LDVL, VR, &LDVR,
400  WORK, &LWORK, INFO);
401 }
402 //=============================================================================
403 void Epetra_LAPACK::GEEV(const char JOBVL, const char JOBVR, const int N, float* A, const int LDA, float* WR, float* WI,
404  float* VL, const int LDVL, float* VR, const int LDVR, float* WORK, const int LWORK, int* INFO) const {
405 
406  SGEEV_F77(CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), &N, A, &LDA, WR, WI, VL, &LDVL, VR, &LDVR,
407  WORK, &LWORK, INFO);
408 }
409 //=============================================================================
410 void Epetra_LAPACK::SPEV(const char JOBZ, const char UPLO, const int N, double* AP, double* W, double* Z, int LDZ, double* WORK, int* INFO) const {
411 
412  DSPEV_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &N, AP, W, Z, &LDZ, WORK, INFO);
413 }
414 //=============================================================================
415 void Epetra_LAPACK::SPEV(const char JOBZ, const char UPLO, const int N, float* AP, float* W, float* Z, int LDZ, float* WORK, int* INFO) const {
416 
417  SSPEV_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &N, AP, W, Z, &LDZ, WORK, INFO);
418 }
419 //=============================================================================
420 void Epetra_LAPACK::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 {
421 
422  DSPGV_F77(&ITYPE, CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &N, AP, BP, W, Z, &LDZ, WORK, INFO);
423 }
424 //=============================================================================
425 void Epetra_LAPACK::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 {
426 
427  SSPGV_F77(&ITYPE, CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &N, AP, BP, W, Z, &LDZ, WORK, INFO);
428 }
429 //=============================================================================
430 void Epetra_LAPACK::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{
431 
432  DSYEV_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &N, A, &LDA, W, WORK, &LWORK, INFO);
433 
434 }
435 //=============================================================================
436 void Epetra_LAPACK::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{
437 
438  SSYEV_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &N, A, &LDA, W, WORK, &LWORK, INFO);
439 
440 }
441 //=============================================================================
442 void Epetra_LAPACK::SYEVD(const char JOBZ, const char UPLO, const int N, double* A, const int LDA, double* W, double* WORK,
443  const int LWORK, int* IWORK, const int LIWORK, int* INFO) const {
444 
445  DSYEVD_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &N, A, &LDA, W, WORK, &LWORK, IWORK, &LIWORK, INFO);
446 }
447 //=============================================================================
448 void Epetra_LAPACK::SYEVD(const char JOBZ, const char UPLO, const int N, float* A, const int LDA, float* W, float* WORK,
449  const int LWORK, int* IWORK, const int LIWORK, int* INFO) const {
450 
451  SSYEVD_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &N, A, &LDA, W, WORK, &LWORK, IWORK, &LIWORK, INFO);
452 }
453 //=============================================================================
454 void Epetra_LAPACK::SYEVX(const char JOBZ, const char RANGE, const char UPLO, const int N, double* A, const int LDA,
455  const double* VL, const double* VU, const int* IL, const int* IU,
456  const double ABSTOL, int * M, double* W, double* Z, const int LDZ, double* WORK,
457  const int LWORK, int* IWORK, int* IFAIL,
458  int* INFO) const {
459 
460  DSYEVX_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(RANGE), CHAR_MACRO(UPLO), &N, A, &LDA, VL, VU, IL, IU,
461  &ABSTOL, M, W, Z, &LDZ, WORK, &LWORK, IWORK, IFAIL, INFO);
462 }
463 //=============================================================================
464 void Epetra_LAPACK::SYEVX(const char JOBZ, const char RANGE, const char UPLO, const int N, float* A, const int LDA,
465  const float* VL, const float* VU, const int* IL, const int* IU,
466  const float ABSTOL, int * M, float* W, float* Z, const int LDZ, float* WORK,
467  const int LWORK, int* IWORK, int* IFAIL,
468  int* INFO) const {
469 
470  SSYEVX_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(RANGE), CHAR_MACRO(UPLO), &N, A, &LDA, VL, VU, IL, IU,
471  &ABSTOL, M, W, Z, &LDZ, WORK, &LWORK, IWORK, IFAIL, INFO);
472 }
473 //=============================================================================
474 void Epetra_LAPACK::SYGV(const int ITYPE, const char JOBZ, const char UPLO, const int N, double* A,
475  const int LDA, double* B, const int LDB, double* W, double* WORK, const int LWORK, int* INFO) const{
476 
477  DSYGV_F77(&ITYPE, CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &N, A, &LDA, B, &LDB, W, WORK, &LWORK, INFO);
478 }
479 //=============================================================================
480 void Epetra_LAPACK::SYGV(const int ITYPE, const char JOBZ, const char UPLO, const int N, float* A,
481  const int LDA, float* B, const int LDB, float* W, float* WORK, const int LWORK, int* INFO) const{
482 
483  SSYGV_F77(&ITYPE, CHAR_MACRO(JOBZ), CHAR_MACRO(UPLO), &N, A, &LDA, B, &LDB, W, WORK, &LWORK, INFO);
484 }
485 //=============================================================================
486 void Epetra_LAPACK::SYGVX(const int /* ITYPE */, const char /* JOBZ */, const char /* RANGE */, const char /* UPLO */, const int /* N */,
487  double* /* A */, const int /* LDA */, double* /* B */, const int /* LDB */, const double* /* VL */, const double* /* VU */,
488  const int* /* IL */, const int* /* IU */, const double /* ABSTOL */, int* /* M */, double* /* W */, double* /* Z */,
489  const int /* LDZ */, double* /* WORK */, const int /* LWORK */, int* /* IWORK */,
490  int* /* IFAIL */, int* /* INFO */) const {
491 
492 #ifdef EPETRA_LAPACK3
493  DSYGVX_F77(&ITYPE, CHAR_MACRO(JOBZ), CHAR_MACRO(RANGE), CHAR_MACRO(UPLO), &N, A, &LDA, B, &LDB, VL, VU,
494  IL, IU, &ABSTOL, M, W, Z, &LDZ, WORK, &LWORK, IWORK,
495  IFAIL, INFO);
496 #else
497 
498  Epetra_Object obj;
499  obj.ReportError("SYGVX requires LAPACK Version 3. Compile Epetra with -DEPETRA_LAPACK3 and link with LAPACK 3 library", -1);
500 #endif
501 }
502 //=============================================================================
503 void Epetra_LAPACK::SYGVX(const int /* ITYPE */, const char /* JOBZ */, const char /* RANGE */, const char /* UPLO */, const int /* N */,
504  float* /* A */, const int /* LDA */, float* /* B */, const int /* LDB */, const float* /* VL */, const float* /* VU */,
505  const int* /* IL */, const int* /* IU */, const float /* ABSTOL */, int* /* M */, float* /* W */, float* /* Z */,
506  const int /* LDZ */, float* /* WORK */, const int /* LWORK */, int* /* IWORK */,
507  int* /* IFAIL */, int* /* INFO */) const {
508 
509 #ifdef EPETRA_LAPACK3
510  SSYGVX_F77(&ITYPE, CHAR_MACRO(JOBZ), CHAR_MACRO(RANGE), CHAR_MACRO(UPLO), &N, A, &LDA, B, &LDB, VL, VU,
511  IL, IU, &ABSTOL, M, W, Z, &LDZ, WORK, &LWORK, IWORK,
512  IFAIL, INFO);
513 #else
514 
515  Epetra_Object obj;
516  obj.ReportError("SYGVX requires LAPACK Version 3. Compile Epetra with -DEPETRA_LAPACK3 and link with LAPACK 3 library", -1);
517 #endif
518 }
519 //=============================================================================
520 void Epetra_LAPACK::SYEVR(const char /* JOBZ */, const char /* RANGE */, const char /* UPLO */, const int /* N */, double* /* A */, const int /* LDA */,
521  const double* /* VL */, const double* /* VU */, const int */* IL */, const int */* IU */,
522  const double /* ABSTOL */, int* /* M */, double* /* W */, double* /* Z */, const int /* LDZ */, int* /* ISUPPZ */, double* /* WORK */, const int /* LWORK */, int* /* IWORK */,
523  const int /* LIWORK */, int* /* INFO */) const {
524 
525 #ifdef EPETRA_LAPACK3
526  DSYEVR_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(RANGE), CHAR_MACRO(UPLO), &N, A, &LDA, VL, VU, IL, IU,
527  &ABSTOL, M, W, Z, &LDZ, ISUPPZ, WORK, &LWORK, IWORK,
528  &LIWORK, INFO);
529 #else
530  Epetra_Object obj;
531  obj.ReportError("SYEVR requires LAPACK Version 3. Compile Epetra with -DEPETRA_LAPACK3 and link with LAPACK 3 library", -1);
532 #endif
533 }
534 //=============================================================================
535 void Epetra_LAPACK::SYEVR(const char /* JOBZ */, const char /* RANGE */, const char /* UPLO */, const int /* N */, float* /* A */, const int /* LDA */,
536  const float* /* VL */, const float* /* VU */, const int */* IL */, const int */* IU */,
537  const float /* ABSTOL */, int* /* M */, float* /* W */, float* /* Z */, const int /* LDZ */, int* /* ISUPPZ */, float* /* WORK */, const int /* LWORK */, int* /* IWORK */,
538  const int /* LIWORK */, int* /* INFO */) const {
539 
540 #ifdef EPETRA_LAPACK3
541  SSYEVR_F77(CHAR_MACRO(JOBZ), CHAR_MACRO(RANGE), CHAR_MACRO(UPLO), &N, A, &LDA, VL, VU, IL, IU,
542  &ABSTOL, M, W, Z, &LDZ, ISUPPZ, WORK, &LWORK, IWORK,
543  &LIWORK, INFO);
544 #else
545  Epetra_Object obj;
546  obj.ReportError("SYEVR requires LAPACK Version 3. Compile Epetra with -DEPETRA_LAPACK3 and link with LAPACK 3 library", -1);
547 #endif
548 }
549 //=============================================================================
550 void Epetra_LAPACK::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,
551  const int LDVL, double* VR, const int LDVR, int* ILO, int* IHI, double* SCALE, double* ABNRM, double* RCONDE,
552  double* RCONDV, double* WORK, const int LWORK, int* IWORK, int* INFO) const{
553 
554  DGEEVX_F77(CHAR_MACRO(BALANC), CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), CHAR_MACRO(SENSE), &N, A, &LDA, WR, WI, VL,
555  &LDVL, VR, &LDVR, ILO, IHI, SCALE, ABNRM, RCONDE,
556  RCONDV, WORK, &LWORK, IWORK, INFO);
557 }
558 //=============================================================================
559 void Epetra_LAPACK::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,
560  const int LDVL, float* VR, const int LDVR, int* ILO, int* IHI, float* SCALE, float* ABNRM, float* RCONDE,
561  float* RCONDV, float* WORK, const int LWORK, int* IWORK, int* INFO) const{
562 
563  SGEEVX_F77(CHAR_MACRO(BALANC), CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), CHAR_MACRO(SENSE), &N, A, &LDA, WR, WI, VL,
564  &LDVL, VR, &LDVR, ILO, IHI, SCALE, ABNRM, RCONDE,
565  RCONDV, WORK, &LWORK, IWORK, INFO);
566 }
567 //=============================================================================
568 void Epetra_LAPACK::GESDD(const char /* JOBZ */, const int /* M */, const int /* N */, double* /* A */, const int /* LDA */, double* /* S */,
569  double* /* U */, const int /* LDU */, double* /* VT */, const int /* LDVT */, double* /* WORK */,
570  const int /* LWORK */, int* /* IWORK */, int* /* INFO */) const{
571 
572 #ifdef EPETRA_LAPACK3
573  DGESDD_F77(CHAR_MACRO(JOBZ), &M, &N, A, &LDA, S, U, &LDU, VT, &LDVT, WORK,
574  &LWORK, IWORK, INFO);
575 #else
576  Epetra_Object obj;
577  obj.ReportError("GESDD requires LAPACK Version 3. Compile Epetra with -DEPETRA_LAPACK3 and link with LAPACK 3 library", -1);
578 #endif
579 }
580 //=============================================================================
581 void Epetra_LAPACK::GESDD(const char /* JOBZ */, const int /* M */, const int /* N */, float* /* A */, const int /* LDA */, float* /* S */,
582  float* /* U */, const int /* LDU */, float* /* VT */, const int /* LDVT */, float* /* WORK */,
583  const int /* LWORK */, int* /* IWORK */, int* /* INFO */) const{
584 
585 #ifdef EPETRA_LAPACK3
586  SGESDD_F77(CHAR_MACRO(JOBZ), &M, &N, A, &LDA, S, U, &LDU, VT, &LDVT, WORK,
587  &LWORK, IWORK, INFO);
588 #else
589  Epetra_Object obj;
590  obj.ReportError("GESDD requires LAPACK Version 3. Compile Epetra with -DEPETRA_LAPACK3 and link with LAPACK 3 library", -1);
591 #endif
592 }
593 //=============================================================================
594 void Epetra_LAPACK::GGEV(const char /* JOBVL */, const char /* JOBVR */, const int /* N */, double* /* A */,
595  const int /* LDA */, double* /* B */, const int /* LDB */, double* /* ALPHAR */, double* /* ALPHAI */,
596  double* /* BETA */, double* /* VL */, const int /* LDVL */, double* /* VR */, const int
597  /* LDVR */, double* /* WORK */, const int /* LWORK */, int* /* INFO */) const{
598 
599 #ifdef EPETRA_LAPACK3
600  DGGEV_F77(CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), &N, A, &LDA, B, &LDB, ALPHAR, ALPHAI,
601  BETA, VL, &LDVL, VR, &LDVR, WORK, &LWORK, INFO);
602 #else
603  Epetra_Object obj;
604  obj.ReportError("GGEV requires LAPACK Version 3. Compile Epetra with -DEPETRA_LAPACK3 and link with LAPACK 3 library", -1);
605 #endif
606 }
607 //=============================================================================
608 void Epetra_LAPACK::GGEV(const char /* JOBVL */, const char /* JOBVR */, const int /* N */, float* /* A */,
609  const int /* LDA */, float* /* B */, const int /* LDB */, float* /* ALPHAR */, float* /* ALPHAI */,
610  float* /* BETA */, float* /* VL */, const int /* LDVL */, float* /* VR */, const int
611  /* LDVR */, float* /* WORK */, const int /* LWORK */, int* /* INFO */) const {
612 
613 #ifdef EPETRA_LAPACK3
614  SGGEV_F77(CHAR_MACRO(JOBVL), CHAR_MACRO(JOBVR), &N, A, &LDA, B, &LDB, ALPHAR, ALPHAI,
615  BETA, VL, &LDVL, VR, &LDVR, WORK, &LWORK, INFO);
616 #else
617  Epetra_Object obj;
618  obj.ReportError("GGEV requires LAPACK Version 3. Compile Epetra with -DEPETRA_LAPACK3 and link with LAPACK 3 library", -1);
619 #endif
620 }
621 //=============================================================================
622 void Epetra_LAPACK::GGLSE(const int M, const int N, const int P, double* A, const int LDA, double* B, const int LDB,
623  double* C, double* D, double* X, double* WORK, const int LWORK, int* INFO) const{
624  DGGLSE_F77(&M, &N, &P, A, &LDA, B, &LDB, C, D, X, WORK, &LWORK, INFO);
625 }
626 //=============================================================================
627 void Epetra_LAPACK::GGLSE(const int M, const int N, const int P, float* A, const int LDA, float* B, const int LDB,
628  float* C, float* D, float* X, float* WORK, const int LWORK, int* INFO) const{
629  SGGLSE_F77(&M, &N, &P, A, &LDA, B, &LDB, C, D, X, WORK, &LWORK, INFO);
630 }
631 //=============================================================================
632 void Epetra_LAPACK::TRTRS(const char UPLO, const char TRANS, const char DIAG, const int N, const int NRHS, const float *A,
633  const int LDA, float *B, const int LDB, int *INFO) const{
634  STRTRS_F77(CHAR_MACRO(UPLO), CHAR_MACRO(TRANS), CHAR_MACRO(DIAG), &N, &NRHS, A, &LDA, B, &LDB, INFO);
635 }
636 //=============================================================================
637 void Epetra_LAPACK::TRTRS(const char UPLO, const char TRANS, const char DIAG, const int N, const int NRHS, const double *A,
638  const int LDA, double *B, const int LDB, int *INFO) const{
639  DTRTRS_F77(CHAR_MACRO(UPLO), CHAR_MACRO(TRANS), CHAR_MACRO(DIAG), &N, &NRHS, A, &LDA, B, &LDB, INFO);
640 }
641 
#define DGEHRD_F77
#define DGETRI_F77
void GGEV(const char JOBVL, const char JOBVR, const int N, double *A, const int LDA, double *B, const int LDB, double *ALPHAR, double *ALPHAI, double *BETA, double *VL, const int LDVL, double *VR, const int LDVR, double *WORK, const int LWORK, int *INFO) const
Epetra_LAPACK wrapper to compute for a pair of N-by-N real nonsymmetric matrices (A,B) the generalized eigenvalues, and optionally, the left and/or right generalized eigenvectors.
void PREFIX DSYGVX_F77(const int *itype, Epetra_fcd, Epetra_fcd, Epetra_fcd, const int *n, double *a, const int *lda, double *b, const int *ldb, const double *vl, const double *vu, const int *il, const int *iu, const double *abstol, int *m, double *w, double *z, const int *ldz, double *work, const int *lwork, int *iwork, int *ifail, int *info)
#define SGETRS_F77
#define SGETRF_F77
#define DGESV_F77
void SYEVX(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, const double ABSTOL, int *M, double *W, double *Z, const int LDZ, double *WORK, const int LWORK, int *IWORK, int *IFAIL, int *INFO) const
Epetra_LAPACK wrapper to compute selected eigenvalues and, optionally, eigenvectors of a real symmetr...
#define SSPEV_F77
#define SHSEQR_F77
void PREFIX SLARFT_F77(Epetra_fcd direct, Epetra_fcd storev, const int *n, const int *k, float *v, const int *ldv, float *tau, float *t, const int *ldt)
#define DORGQR_F77
#define SGESVX_F77
#define DPOTRI_F77
#define DGGEV_F77
void PREFIX SSPGV_F77(const int *itype, Epetra_fcd, Epetra_fcd, const int *n, float *ap, float *bp, float *w, float *z, const int *ldz, float *work, int *info)
#define SSYGV_F77
#define DGELS_F77
#define DSPEV_F77
void HSEQR(const char JOB, const char COMPZ, const int N, const int ILO, const int IHI, float *H, const int LDH, float *WR, float *WI, float *Z, const int LDZ, float *WORK, const int LWORK, int *INFO) const
Epetra_LAPACK wrapper for computing the eigenvalues of a real upper Hessenberg matrix (SHSEQR) ...
#define SGEHRD_F77
#define SORGQR_F77
void GESV(const int N, const int NRHS, float *A, const int LDA, int *IPIV, float *X, const int LDX, int *INFO) const
Epetra_LAPACK factor and solve for general matrix (SGESV)
#define DPOEQU_F77
#define SLAMCH_F77
#define DGETRF_F77
void GESVD(const char JOBU, const char JOBVT, const int M, const int N, float *A, const int LDA, float *S, float *U, const int LDU, float *VT, const int LDVT, float *WORK, const int *LWORK, int *INFO) const
Epetra_LAPACK wrapper for computing the singular value decomposition (SGESVD)
void SYGV(const int ITYPE, const char JOBZ, const char UPLO, const int N, double *A, const int LDA, double *B, const int LDB, double *W, double *WORK, const int LWORK, int *INFO) const
Epetra_LAPACK wrapper to compute all the eigenvalues, and optionally, the eigenvectors of a real gene...
#define DPOSVX_F77
#define SPOTRI_F77
#define SGGEV_F77
#define SPOTRS_F77
void GEQRF(const int M, const int N, float *A, const int LDA, float *TAU, float *WORK, const int lwork, int *INFO) const
Epetra_LAPACK QR factorization for general matrix (SGEQRF)
#define SGESVD_F77
void PREFIX DGGSVD_F77(Epetra_fcd, Epetra_fcd, Epetra_fcd, const int *m, const int *n, const int *p, int *k, int *l, double *a, const int *lda, double *b, const int *ldb, double *alpha, double *beta, double *u, const int *ldu, double *v, const int *ldv, double *q, const int *ldq, double *work, int *iwork, int *info)
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
Epetra_LAPACK wrapper to compute all the eigenvalues and, optionally, the eigenvectors of a real gene...
#define SPOSV_F77
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, const int LWORK, int *IWORK, int *INFO) const
Epetra_LAPACK wrapper to compute the singular value decomposition (SVD) of a real M-by-N matrix A...
#define DPOSV_F77
#define SGECON_F77
void POCON(const char UPLO, const int N, const float *A, const int LDA, const float ANORM, float *RCOND, float *WORK, int *IWORK, int *INFO) const
Epetra_LAPACK condition number estimator for positive definite matrix (SPOCON)
void ORMHR(const char SIDE, const char TRANS, const int M, const int N, const int ILO, const int IHI, const float *A, const int LDA, const float *TAU, float *C, const int LDC, float *WORK, const int LWORK, int *INFO) const
Epetra_LAPACK wrapper for applying an orthogonal matrix in-place (SORMHR)
void PREFIX SSYGVX_F77(const int *itype, Epetra_fcd, Epetra_fcd, Epetra_fcd, const int *n, float *a, const int *lda, float *b, const int *ldb, const float *vl, const float *vu, const int *il, const int *iu, const float *abstol, int *m, float *w, float *z, const int *ldz, float *work, const int *lwork, int *iwork, int *ifail, int *info)
void GETRF(const int M, const int N, float *A, const int LDA, int *IPIV, int *INFO) const
Epetra_LAPACK factorization for general matrix (SGETRF)
#define STREVC_F77
#define SORMHR_F77
void PREFIX DSYEVD_F77(Epetra_fcd, Epetra_fcd, const int *n, double *a, const int *lda, double *w, double *work, const int *lwork, int *iwork, const int *liwork, int *info)
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, const 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
Epetra_LAPACK solve driver for general matrix (SGESVX)
#define CHAR_MACRO(char_var)
void PREFIX DSYEVR_F77(Epetra_fcd, Epetra_fcd, Epetra_fcd, const int *n, double *a, const int *lda, const double *vl, const double *vu, const int *il, const int *iu, const double *abstol, int *m, double *w, double *z, const int *ldz, int *isuppz, double *work, const int *lwork, int *iwork, const int *liwork, int *info)
#define DGEQRF_F77
void GECON(const char NORM, const int N, const float *A, const int LDA, const float ANORM, float *RCOND, float *WORK, int *IWORK, int *INFO) const
Epetra_LAPACK condition number estimator for general matrix (SGECON)
#define DTREXC_F77
void GGLSE(const int M, const int N, const int P, double *A, const int LDA, double *B, const int LDB, double *C, double *D, double *X, double *WORK, const int LWORK, int *INFO) const
Epetra_LAPACK wrapper to solve the linear equality-constrained least squares (LSE) problem...
Epetra_Object: The base Epetra class.
Definition: Epetra_Object.h:65
void SPEV(const char JOBZ, const char UPLO, const int N, double *AP, double *W, double *Z, int LDZ, double *WORK, int *INFO) const
Epetra_LAPACK wrapper to compute all the eigenvalues and, optionally, eigenvectors of a real symmetri...
void PREFIX DSYEVX_F77(Epetra_fcd, Epetra_fcd, Epetra_fcd, const int *n, double *a, const int *lda, const double *vl, const double *vu, const int *il, const int *iu, const double *abstol, int *m, double *w, double *z, const int *ldz, double *work, const int *lwork, int *iwork, int *ifail, int *info)
#define STRTRS_F77
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, const double ABSTOL, int *M, double *W, double *Z, const int LDZ, int *ISUPPZ, double *WORK, const int LWORK, int *IWORK, const int LIWORK, int *INFO) const
Epetra_LAPACK wrapper to compute selected eigenvalues and, optionally, eigenvectors of a real symmetr...
#define SPOTRF_F77
void POEQU(const int N, const float *A, const int LDA, float *S, float *SCOND, float *AMAX, int *INFO) const
Epetra_LAPACK equilibration for positive definite matrix (SPOEQU)
#define DGEEV_F77
void PREFIX SGESDD_F77(Epetra_fcd, 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, const int *lwork, int *iwork, int *info)
#define DPOCON_F77
void SYGVX(const int ITYPE, const char JOBZ, const char RANGE, const char UPLO, const int N, double *A, const int LDA, double *B, const int LDB, const double *VL, const double *VU, const int *IL, const int *IU, const double ABSTOL, int *M, double *W, double *Z, const int LDZ, double *WORK, const int LWORK, int *IWORK, int *IFAIL, int *INFO) const
Epetra_LAPACK wrapper to compute selected eigenvalues, and optionally, eigenvectors of a real general...
#define SGEEVX_F77
#define SORGHR_F77
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
Epetra_LAPACK for forming the triangular factor of a product of elementary Householder reflectors (SL...
void PREFIX DLARFT_F77(Epetra_fcd direct, Epetra_fcd storev, const int *n, const int *k, double *v, const int *ldv, double *tau, double *t, const int *ldt)
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
Epetra_LAPACK wrapper for generating a real orthogonal matrix Q defined by elementary reflectors...
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
Epetra_LAPACK solve (after factorization) for general matrix (SGETRS)
void TREVC(const char SIDE, const char HOWMNY, int *SELECT, const int N, const float *T, const int LDT, float *VL, const int LDVL, float *VR, const int LDVR, const int MM, int *M, float *WORK, int *INFO) const
Epetra_LAPACK wrapper for computing eigenvectors of a quasi-triangular/triagnular matrix (STREVC) ...
void GETRI(const int N, float *A, const int LDA, int *IPIV, float *WORK, const int *LWORK, int *INFO) const
Epetra_LAPACK inversion for general matrix (SGETRI)
void TRTRS(const char UPLO, const char TRANS, const char DIAG, const int N, const int NRHS, const float *A, const int LDA, float *B, const int LDB, int *INFO) const
Epetra_LAPACK wrapper for TRTRS routine.
#define DGEEQU_F77
#define DGESVX_F77
#define SGETRI_F77
#define SPOCON_F77
#define SGGLSE_F77
void PREFIX SGGSVD_F77(Epetra_fcd, Epetra_fcd, Epetra_fcd, const int *m, const int *n, const int *p, int *k, int *l, float *a, const int *lda, float *b, const int *ldb, float *alpha, float *beta, float *u, const int *ldu, float *v, const int *ldv, float *q, const int *ldq, float *work, int *iwork, int *info)
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
Epetra_LAPACK equilibration for general matrix (SGEEQU)
void PREFIX DSPGV_F77(const int *itype, Epetra_fcd, Epetra_fcd, const int *n, double *ap, double *bp, double *w, double *z, const int *ldz, double *work, int *info)
#define DPOTRS_F77
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
Epetra_LAPACK wrapper for reduction to Hessenberg form (SGEHRD)
#define DSYEV_F77
void POSV(const char UPLO, const int N, const int NRHS, float *A, const int LDA, float *X, const int LDX, int *INFO) const
Epetra_LAPACK factor and solve for positive definite matrix (SPOSV)
#define DSYGV_F77
#define SPORFS_F77
#define DGESVD_F77
#define SGEQRF_F77
#define SSYEV_F77
void GELS(const char TRANS, const int M, const int N, const int NRHS, double *A, const int LDA, double *B, const int LDB, double *WORK, const int LWORK, int *INFO) const
Epetra_LAPACK simple driver to solve least-squares systems.
void SYEVD(const char JOBZ, const char UPLO, const int N, double *A, const int LDA, double *W, double *WORK, const int LWORK, int *IWORK, const int LIWORK, int *INFO) const
Epetra_LAPACK wrapper to compute all eigenvalues and, optionally, eigenvectors of a real symmetric ma...
void POTRF(const char UPLO, const int N, float *A, const int LDA, int *INFO) const
Epetra_LAPACK factorization for positive definite matrix (SPOTRF)
#define SGERFS_F77
#define DHSEQR_F77
void POSVX(const char FACT, const char UPLO, const int N, const int NRHS, float *A, const int LDA, float *AF, const int LDAF, const 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
Epetra_LAPACK solve driver for positive definite matrix (SPOSVX)
void PORFS(const char UPLO, const int N, const int NRHS, const float *A, const int LDA, const float *AF, const int LDAF, const float *B, const int LDB, float *X, const int LDX, float *FERR, float *BERR, float *WORK, int *IWORK, int *INFO) const
Epetra_LAPACK solve driver for positive definite matrix (SPOSVX)
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, double *ALPHA, double *BETA, double *U, const int LDU, double *V, const int LDV, double *Q, const int LDQ, double *WORK, int *IWORK, int *INFO) const
Epetra_LAPACK wrapper to compute the generalized singular value decomposition (GSVD) of an M-by-N rea...
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
Epetra_LAPACK solve (after factorization) for positive definite matrix (SPOTRS)
void PREFIX SSYEVR_F77(Epetra_fcd, Epetra_fcd, Epetra_fcd, const int *n, float *a, const int *lda, const float *vl, const float *vu, const int *il, const int *iu, const float *abstol, int *m, float *w, float *z, const int *ldz, int *isuppz, float *work, const int *lwork, int *iwork, const int *liwork, int *info)
void PREFIX SSYEVD_F77(Epetra_fcd, Epetra_fcd, const int *n, float *a, const int *lda, float *w, float *work, const int *lwork, int *iwork, const int *liwork, int *info)
#define SGEEQU_F77
#define DTREVC_F77
void GEEV(const char JOBVL, const char JOBVR, const int N, double *A, const int LDA, double *WR, double *WI, double *VL, const int LDVL, double *VR, const int LDVR, double *WORK, const int LWORK, int *INFO) const
Epetra_LAPACK wrapper to compute for an N-by-N real nonsymmetric matrix A, the eigenvalues and...
#define DPOTRF_F77
virtual int ReportError(const std::string Message, int ErrorCode) const
Error reporting method.
#define SGEEV_F77
#define DGECON_F77
#define DPORFS_F77
void GEEVX(const char BALANC, const char JOBVL, const char JOBVR, const char SENSE, const int N, double *A, const int LDA, double *WR, double *WI, double *VL, const int LDVL, double *VR, const int LDVR, int *ILO, int *IHI, double *SCALE, double *ABNRM, double *RCONDE, double *RCONDV, double *WORK, const int LWORK, int *IWORK, int *INFO) const
Epetra_LAPACK wrapper to compute for an N-by-N real nonsymmetric matrix A, the eigenvalues and...
void GERFS(const char TRANS, const int N, const int NRHS, const float *A, const int LDA, const float *AF, const int LDAF, const int *IPIV, const float *B, const int LDB, float *X, const int LDX, float *FERR, float *BERR, float *WORK, int *IWORK, int *INFO) const
Epetra_LAPACK Refine solution (GERFS)
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
Epetra_LAPACK wrapper to compute all eigenvalues and, optionally, eigenvectors of a real symmetric ma...
#define DGEEVX_F77
#define DORMHR_F77
#define SPOEQU_F77
void POTRI(const char UPLO, const int N, float *A, const int LDA, int *INFO) const
Epetra_LAPACK inversion for positive definite matrix (SPOTRI)
#define STREXC_F77
#define DORGHR_F77
void PREFIX SSYEVX_F77(Epetra_fcd, Epetra_fcd, Epetra_fcd, const int *n, float *a, const int *lda, const float *vl, const float *vu, const int *il, const int *iu, const float *abstol, int *m, float *w, float *z, const int *ldz, float *work, const int *lwork, int *iwork, int *ifail, int *info)
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
Epetra_LAPACK wrapper for generating a m x n real matrix Q with orthonormal columns, defined as the product of k elementary reflectors. (SORGQR)
#define SPOSVX_F77
void LAMCH(const char CMACH, float &T) const
Epetra_LAPACK wrapper for DLAMCH routine. On out, T holds machine double precision floating point cha...
#define SGESV_F77
#define DGETRS_F77
#define DGERFS_F77
#define DTRTRS_F77
#define DLAMCH_F77
void PREFIX DGESDD_F77(Epetra_fcd, 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, const int *lwork, int *iwork, int *info)
void TREXC(const char COMPQ, const int N, float *T, const int LDT, float *Q, const int LDQ, int IFST, int ILST, float *WORK, int *INFO) const
Epetra_LAPACK wrapper for reordering the real-Schur/Schur factorization of a matrix (STREXC) ...
#define DGGLSE_F77