MADNESS  0.10.1
clapack_fortran.h
Go to the documentation of this file.
1 /*
2  This file is part of MADNESS.
3 
4  Copyright (C) 2007,2010 Oak Ridge National Laboratory
5 
6  This program is free software; you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation; either version 2 of the License, or
9  (at your option) any later version.
10 
11  This program is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with this program; if not, write to the Free Software
18  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19 
20  For more information please contact:
21 
22  Robert J. Harrison
23  Oak Ridge National Laboratory
24  One Bethel Valley Road
25  P.O. Box 2008, MS-6367
26 
27  email: harrisonrj@ornl.gov
28  tel: 865-241-3937
29  fax: 865-572-0680
30 
31 
32  $Id$
33 */
34 
35 
36 #ifndef MADNESS_LINALG_CLAPACK_FORTRAN_H__INCLUDED
37 #define MADNESS_LINALG_CLAPACK_FORTRAN_H__INCLUDED
38 
39 /// \file clapack_fortran.h
40 /// \brief Legacy C++ prototypes for Fortran LAPACK with associated typedefs and macos
41 
42 #include <madness/fortran_ctypes.h>
43 
44 #ifdef FORTRAN_LINKAGE_LC
45 # define sgesvd_ sgesvd
46 # define dgesvd_ dgesvd
47 # define cgesvd_ cgesvd
48 # define zgesvd_ zgesvd
49 
50 # define sgesv_ sgesv
51 # define dgesv_ dgesv
52 # define cgesv_ cgesv
53 # define zgesv_ zgesv
54 
55 # define sgelss_ sgelss
56 # define dgelss_ dgelss
57 # define cgelss_ cgelss
58 # define zgelss_ zgelss
59 
60 # define sgels_ sgels
61 # define dgels_ dgels
62 # define cgels_ cgels
63 # define zgels_ zgels
64 
65 # define ssyev_ ssyev
66 # define dsyev_ dsyev
67 # define cheev_ cheev
68 # define zheev_ zheev
69 
70 # define sggev_ sggev
71 # define dggev_ dggev
72 # define cggev_ cggev
73 # define zggev_ zggev
74 
75 #ifndef MADNESS_HAS_ELEMENTAL
76 # define ssygv_ ssygv
77 # define dsygv_ dsygv
78 # define chegv_ chegv
79 # define zhegv_ zhegv
80 #endif
81 
82 # define spotrf_ spotrf
83 # define cpotrf_ cpotrf
84 # define dpotrf_ dpotrf
85 # define zpotrf_ zpotrf
86 
87 # define sgetrf_ sgetrf
88 # define cgetrf_ cgetrf
89 # define dgetrf_ dgetrf
90 # define zgetrf_ zgetrf
91 
92 # define sgetri_ sgetri
93 # define cgetri_ cgetri
94 # define dgetri_ dgetri
95 # define zgetri_ zgetri
96 
97 # define strsm_ strsm
98 # define ctrsm_ ctrsm
99 # define dtrsm_ dtrsm
100 # define ztrsm_ ztrsm
101 
102 # define dlamch_ dlamch
103 # define slamch_ slamch
104 
105 # define sgeev_ sgeev
106 # define cgeev_ cgeev
107 # define dgeev_ dgeev
108 # define zgeev_ zgeev
109 
110 #else
111  // only lowercase with zero and one underscores are handled -- if detected another convention complain loudly
112 # ifndef FORTRAN_LINKAGE_LCU
113 # error "clapack.h does not support the current Fortran symbol convention -- please, edit and check in the changes."
114 # endif
115 #endif
116 
117 // SUBROUTINE DLAMCH( CMACH, RESULT )
118 
119 // PURPOSE
120 // DLAMCH determines double precision machine parameters.
121 
122 extern "C"
123  float slamch_(const char* mode, int modelen);
124 extern "C"
125  double dlamch_(const char* mode, int modelen);
126 
127 // SUBROUTINE DGESVD( JOBU, JOBVT, M, N, A, LDA, S, U, LDU, VT, LDVT, WORK, LWORK, INFO )
128 
129 // PURPOSE
130 // DGESVD computes the singular value decomposition (SVD) of a real
131 // M-by-N matrix A, optionally computing the left and/or right singular
132 // vectors.
133 
134 extern "C"
135  void sgesvd_(const char *jobu, const char *jobvt, integer *m, integer *n,
136  real4 *a, integer *lda, real4 *s, real4 *u, integer *ldu,
137  real4 *vt, integer *ldvt, real4 *work, integer *lwork,
138  integer *info, char_len jobulen, char_len jobvtlen);
139 extern "C"
140  void dgesvd_(const char *jobu, const char *jobvt, integer *m, integer *n,
141  real8 *a, integer *lda, real8 *s, real8 *u, integer *ldu,
142  real8 *vt, integer *ldvt, real8 *work, integer *lwork,
143  integer *info, char_len jobulen, char_len jobvtlen);
144 extern "C"
145  void cgesvd_(const char *jobu, const char *jobvt, integer *m, integer *n,
147  integer *ldu, complex_real4 *vt, integer *ldvt, complex_real4 *work,
148  integer *lwork, real4 *rwork,
149  integer *info, char_len jobulen, char_len jobvtlen);
150 extern "C"
151  void zgesvd_(const char *jobu, const char *jobvt, integer *m, integer *n,
153  integer *ldu, complex_real8 *vt, integer *ldvt, complex_real8 *work,
154  integer *lwork, real8 *rwork,
155  integer *info, char_len jobulen, char_len jobvtlen);
156 
157 // SUBROUTINE DGESV( N, NRHS, A, LDA, IPIV, B, LDB, INFO )
158 
159 // PURPOSE
160 // DGESV computes the solution to a real system of linear equations
161 // A * X = B,
162 // where A is an N-by-N matrix and X and B are N-by-NRHS matrices.
163 
164 extern "C"
165  void sgesv_(integer* n, integer* nrhs, real4* AT, integer* lda,
166  integer* piv, real4* x, integer* ldx, integer* info);
167 extern "C"
168  void dgesv_(integer* n, integer* nrhs, real8* AT, integer* lda,
169  integer* piv, real8* x, integer* ldx, integer* info);
170 extern "C"
171  void cgesv_(integer* n, integer* nrhs, complex_real4* AT, integer* lda,
172  integer* piv, complex_real4* x, integer* ldx, integer* info);
173 extern "C"
174  void zgesv_(integer* n, integer* nrhs, complex_real8* AT, integer* lda,
175  integer* piv, complex_real8* x, integer* ldx, integer* info);
176 
177 // SUBROUTINE DGELSS( M, N, NRHS, A, LDA, B, LDB, S, RCOND, RANK, WORK, LWORK, INFO )
178 
179 // PURPOSE
180 // DGELSS computes the minimum norm solution to a real linear least
181 // squares problem:
182 // Minimize 2-norm(| b - A*x |)
183 // using the singular value decomposition (SVD) of A. A is an M-by-N
184 // matrix which may be rank-deficient.
185 
186 extern "C"
187  void sgelss_(integer *m, integer *n, integer *nrhs,
188  real4 *a, integer *lda, real4 *b, integer *ldb, real4 *sOUT,
189  real4 *rcondIN, integer *rankOUT, real4 *work,
190  integer *lwork, integer *infoOUT);
191 extern "C"
192  void dgelss_(integer *m, integer *n, integer *nrhs,
193  real8 *a, integer *lda, real8 *b, integer *ldb, real8 *sOUT,
194  real8 *rcondIN, integer *rankOUT, real8 *work,
195  integer *lwork, integer *infoOUT);
196 extern "C"
197  void cgelss_(integer *m, integer *n, integer *nrhs,
198  complex_real4 *a, integer *lda, complex_real4 *b, integer *ldb,
199  real4 *sOUT,
200  real4 *rcondIN, integer *rankOUT, complex_real4 *work,
201  integer *lwork, real4 *rwork, integer *infoOUT);
202 extern "C"
203  void zgelss_(integer *m, integer *n, integer *nrhs,
204  complex_real8 *a, integer *lda, complex_real8 *b, integer *ldb,
205  real8 *sOUT,
206  real8 *rcondIN, integer *rankOUT, complex_real8 *work,
207  integer *lwork, real8 *rwork, integer *infoOUT);
208 
209 // SUBROUTINE DGELS( TRANS, M, N, NRHS, A, LDA, B, LDB, WORK, LWORK, INFO )
210 
211 // PURPOSE
212 // DGELS solves overdetermined or underdetermined real linear systems
213 // involving an M-by-N matrix A, or its transpose, using a QR or LQ
214 // factorization of A. It is assumed that A has full rank.
215 
216 extern "C"
217  void sgels_(const char *trans, integer *m, integer *n, integer *nrhs,
218  real4 *a, integer *lda, real4 *b, integer *ldb, real4 *work,
219  integer *lwork, integer *infoOUT, char_len translen);
220 extern "C"
221  void dgels_(const char *trans, integer *m, integer *n, integer *nrhs,
222  real8 *a, integer *lda, real8 *b, integer *ldb, real8 *work,
223  integer *lwork, integer *infoOUT, char_len translen);
224 extern "C"
225  void cgels_(const char *trans, integer *m, integer *n, integer *nrhs,
226  complex_real4 *a, integer *lda, complex_real4 *b, integer *ldb,
227  complex_real4 *work,
228  integer *lwork, real4 *rwork, integer *infoOUT, char_len translen);
229 extern "C"
230  void zgels_(const char *trans, integer *m, integer *n, integer *nrhs,
231  complex_real8 *a, integer *lda, complex_real8 *b, integer *ldb,
232  complex_real8 *work,
233  integer *lwork, real8 *rwork, integer *infoOUT, char_len translen);
234 
235 // SUBROUTINE DGGEV( JOBZ, UPLO, N, A, LDA, B, LDB, ALPHAR, ALPHAI, BETA, VL, LDVL, VR, LDVR, WORK, LWORK, INFO )
236 
237 // PURPOSE
238 // DGGEV computes for a pair of N-by-N real nonsymmetric matrices (A,B)
239 // the generalized eigenvalues, and optionally, the left and/or right
240 // generalized eigenvectors.
241 
242 extern "C"
243  void sggev_(const char* jobz, const char* uplo, integer *n,
244  real4* a, integer* lda, real4* b, integer* ldb,
245  real4* alphar, real4* alphai, real4* beta,
246  real4* vl, integer* ldvl, real4* vr, integer* ldvr,
247  real4* work, integer* lwork, integer* info,
248  char_len jobzlen, char_len uplo_len);
249 extern "C"
250  void dggev_(const char* jobl, const char* jobr, integer *n,
251  real8 *a, integer *lda, real8 *b, integer *ldb,
252  real8 *w_real, real8 *w_imag, real8 *beta,
253  real8 *vl, integer *ldvl, real8 *vr, integer *ldvr,
254  real8 *work, integer *lwork, integer *info,
255  char_len jobzlen, char_len uplo_len);
256 extern "C"
257  void cggev_(const char* jobz, const char* uplo, integer *n,
258  complex_real4* a, integer* lda, complex_real4* b, integer* ldb,
260  complex_real4* vl, integer* ldvl, complex_real4* vr, integer* ldvr,
261  complex_real4* work, integer* lwork, real4* rwork, integer* info,
262  char_len jobzlen, char_len uplo_len);
263 extern "C"
264  void zggev_(const char* jobz, const char* uplo, integer *n,
265  complex_real8* a, integer* lda, complex_real8* b, integer* ldb,
267  complex_real8* vl, integer* ldvl, complex_real8* vr, integer* ldvr,
268  complex_real8* work, integer* lwork, real8* rwork, integer* info,
269  char_len jobzlen, char_len uplo_len);
270 
271 // SUBROUTINE DGEEV( JOBZ, UPLO, N, A, LDA, W, WORK, LWORK, INFO )
272 
273 // PURPOSE
274 // DGEEV computes for an N-by-N real nonsymmetric matrix A, the
275 // eigenvalues and, optionally, the left and/or right eigenvectors.
276 
277 extern "C"
278  void sgeev_(const char* jobz, const char* uplo, integer *n, real4* a, integer* lda,
279  real4* w_real, real4* w_imag, real4* v, integer* ldv, real4* vr, integer* ldvr,
280  real4* work, integer* lwork, integer* info,
281  char_len jobzlen, char_len uplo_len );
282 extern "C"
283  void dgeev_(const char* jobz, const char* uplo, integer *n,
284  real8* a, integer* lda, real8* w_real, real8* w_imag, real8* v, integer* ldv,
285  real8* vr, integer* ldvr, real8* work, integer* lwork, integer* info,
286  char_len jobzlen, char_len uplo_len );
287 extern "C"
288  void cgeev_(const char* jobz, const char* uplo, integer *n, complex_real4* a, integer* lda,
289  complex_real4* w, complex_real4* vl, integer* ldvl, complex_real4* vr, integer* ldvr,
290  complex_real4* work, integer* lwork, real4* rwork, integer* info,
291  char_len jobzlen, char_len uplo_len );
292 extern "C"
293  void zgeev_(const char* jobz, const char* uplo, integer *n, complex_real8* a, integer* lda,
294  complex_real8* w, complex_real8* vl, integer* ldvl, complex_real8* vr, integer* ldvr,
295  complex_real8* work, integer* lwork, real8* rwork, integer* info,
296  char_len jobzlen, char_len uplo_len );
297 
298 // SUBROUTINE DSYEV( JOBZ, UPLO, N, A, LDA, W, WORK, LWORK, INFO )
299 
300 // PURPOSE
301 // DSYEV computes all eigenvalues and, optionally, eigenvectors of a
302 // real symmetric matrix A.
303 
304 extern "C"
305  void ssyev_(const char* jobz, const char* uplo, integer *n,
306  real4 *a, integer *lda, real4 *w, real4 *work, integer *lwork,
307  integer *info, char_len jobzlen, char_len uplo_len );
308 extern "C"
309  void dsyev_(const char* jobz, const char* uplo, integer *n,
310  real8 *a, integer *lda, real8 *w, real8 *work, integer *lwork,
311  integer *info, char_len jobzlen, char_len uplo_len );
312 extern "C"
313  void cheev_(const char* jobz, const char* uplo, integer *n,
314  complex_real4 *a, integer *lda, real4 *w, complex_real4 *work,
315  integer *lwork, real4 *rwork,
316  integer *info, char_len jobzlen, char_len uplo_len );
317 extern "C"
318  void zheev_(const char* jobz, const char* uplo, integer *n,
319  complex_real8 *a, integer *lda, real8 *w, complex_real8 *work,
320  integer *lwork, real8 *rwork,
321  integer *info, char_len jobzlen, char_len uplo_len );
322 
323 // SUBROUTINE DSYGV( ITYPE, JOBZ, UPLO, N, A, LDA, B, LDB, W, WORK, LWORK, INFO )
324 // SUBROUTINE ZHEGV( ITYPE, JOBZ, UPLO, N, A, LDA, B, LDB, W, WORK, LWORK, INFO )
325 
326 // PURPOSE
327 // DSYGV computes all the eigenvalues, and optionally, the eigenvectors
328 // of a real generalized symmetric-definite eigenproblem, of the form
329 // A*x=(lambda)*B*x, A*Bx=(lambda)*x, or B*A*x=(lambda)*x.
330 
331 extern "C"
332  void ssygv_(integer *itype, const char* jobz, const char* uplo, integer *n,
333  real4 *a, integer *lda, real4 *b, integer *ldb,
334  real4 *w, real4 *work, integer *lwork,
335  integer *info, char_len jobzlen, char_len uplo_len );
336 extern "C"
337  void dsygv_(integer *itype, const char* jobz, const char* uplo, integer *n,
338  real8 *a, integer *lda, real8 *b, integer *ldb,
339  real8 *w, real8 *work, integer *lwork,
340  integer *info, char_len jobzlen, char_len uplo_len );
341 extern "C"
342  void chegv_(integer *itype, const char* jobz, const char* uplo, integer *n,
343  complex_real4 *a, integer *lda, complex_real4 *b, integer *ldb,
344  real4 *w, complex_real4 *work, integer *lwork, real4 *rwork,
345  integer *info, char_len jobzlen, char_len uplo_len );
346 extern "C"
347  void zhegv_(integer *itype, const char* jobz, const char* uplo, integer *n,
348  complex_real8 *a, integer *lda, complex_real8 *b, integer *ldb,
349  real8 *w, complex_real8 *work, integer *lwork, real8 *rwork,
350  integer *info, char_len jobzlen, char_len uplo_len );
351 
352 // SUBROUTINE DGEQRF (M, N, A, LDA, TAU, WORK, LWORK, INFO)
353 //
354 // PURPOSE
355 // DGEQRF computes a QR factorization of a real M-by-N matrix A:
356 // A = Q * R.
357 
358 extern "C"
360  real4 *a, integer *lda, real4 *tau,
361  real4 *work, integer *lwork, integer *infoOUT);
362 
363 extern "C"
365  real8 *a, integer *lda, real8 *tau,
366  real8 *work, integer *lwork, integer *infoOUT);
367 
368 extern "C"
370  complex_real4 *a, integer *lda, complex_real4 *tau,
371  complex_real4 *work, integer *lwork, integer *infoOUT);
372 
373 extern "C"
375  complex_real8 *a, integer *lda, complex_real8 *tau,
376  complex_real8 *work, integer *lwork, integer *infoOUT);
377 
378 // SUBROUTINE DGEQP3(M, N, A, LDA, JPVT, TAU, WORK, LWORK, INFO );
379 
380 // PURPOSE
381 // DGEQP3 computes a QR factorization with column pivoting of a
382 // matrix A: A*P = Q*R using Level 3 BLAS.
383 
384 extern "C"
386  real4 *a, integer *lda, integer *jpvt, real4 *tau,
387  real4 *work, integer *lwork, integer *infoOUT);
388 
389 extern "C"
391  real8 *a, integer *lda, integer *jpvt, real8 *tau,
392  real8 *work, integer *lwork, integer *infoOUT);
393 
394 extern "C"
396  integer *lda, integer *jpvt, complex_real4 *tau,
397  complex_real4 *work, integer *lwork, real4 *rwork,
398  integer *infoOUT);
399 
400 extern "C"
402  integer *lda, integer *jpvt, complex_real8 *tau,
403  complex_real8 *work, integer *lwork, real8 *rwork,
404  integer *infoOUT);
405 
406 // SUBROUTINE DORGQR( M, N, K, A, LDA, TAU, WORK, LWORK, INFO )
407 // SUBROUTINE ZUNGQR( M, N, K, A, LDA, TAU, WORK, LWORK, INFO )
408 
409 // PURPOSE
410 // DORGQR generates an M-by-N real matrix Q with orthonormal
411 // columns, which is defined as the first N columns of a pro-
412 // duct of K elementary reflectors of order M
413 
414 extern "C"
416  real4 *a, integer *lda, real4 *tau,
417  real4 *work, integer *lwork, integer *info);
418 extern "C"
420  real8 *a, integer *lda, real8 *tau,
421  real8 *work, integer *lwork, integer *info);
422 extern "C"
424  complex_real4 *a, integer *lda, complex_real4 *tau,
425  complex_real4 *work, integer *lwork, integer *info);
426 extern "C"
428  complex_real8 *a, integer *lda, complex_real8 *tau,
429  complex_real8 *work, integer *lwork, integer *info);
430 
431 // SUBROUTINE DPOTRF( UPLO, N, A, LDA, INFO )
432 
433 // PURPOSE
434 // computes the Cholesky factorization of a real symmetric
435 // positive definite matrix A.
436 
437 extern "C"
438 void spotrf_(const char *uplo, const integer* n, real4 *a, const integer *lda, integer *info, char_len uplo_len);
439 extern "C"
440 void dpotrf_(const char *uplo, const integer* n, real8 *a, const integer *lda, integer *info, char_len uplo_len);
441 extern "C"
442 void cpotrf_(const char *uplo, const integer* n, complex_real4 *a, const integer *lda, integer *info, char_len uplo_len);
443 extern "C"
444 void zpotrf_(const char *uplo, const integer* n, complex_real8 *a, const integer *lda, integer *info, char_len uplo_len);
445 
446 // SUBROUTINE DPSTRF( UPLO, N, A, LDA, IPIV, RANK, TOL, WORK, INFO )
447 
448 // PURPOSE
449 // DPSTRF computes the Cholesky factorization with complete
450 // pivoting of a real symmetric positive semidefinite matrix A.
451 
452 extern "C"
453 void spstrf_(const char *uplo, const integer* n, real4 *a, const integer *lda, integer* ipiv, integer* rank, real4* tol,
454  real4* work, integer *info);
455 extern "C"
456 void dpstrf_(const char *uplo, const integer* n, real8 *a, const integer *lda, integer* ipiv, integer* rank, real8* tol,
457  real8* work, integer *info);
458 extern "C"
459 void cpstrf_(const char *uplo, const integer* n, complex_real4 *a, const integer *lda, integer* ipiv, integer* rank, real4* tol,
460  complex_real4* work, integer *info);
461 extern "C"
462 void zpstrf_(const char *uplo, const integer* n, complex_real8 *a, const integer *lda, integer* ipiv, integer* rank, real8* tol,
463  complex_real8* work, integer *info);
464 
465 // SUBROUTINE DGETRF( M, N, A, LDA, IPIV, INFO )
466 
467 // PURPOSE
468 // DGETRF computes an LU factorization of a general M-by-N matrix A
469 // using partial pivoting with row interchanges.
470 
471 extern "C"
472 void sgetrf_(const integer* m, const integer* n, real4 *a, const integer *lda,
473  integer* ipiv, integer *info);
474 extern "C"
475 void dgetrf_(const integer* m, const integer* n, real8 *a, const integer *lda,
476  integer* ipiv, integer *info);
477 extern "C"
478 void cgetrf_(const integer* m, const integer* n, complex_real4 *a, const integer *lda,
479  integer* ipiv, integer *info);
480 extern "C"
481 void zgetrf_(const integer* m, const integer* n, complex_real8 *a, const integer *lda,
482  integer* ipiv, integer *info);
483 
484 // SUBROUTINE DGETRI( N, A, LDA, IPIV, WORK, LWORK, INFO )
485 
486 // PURPOSE
487 // DGETRI computes the inverse of a matrix using the LU factorization
488 // computed by DGETRF.
489 
490 extern "C"
491 void sgetri_(const integer* n, real4 *a, const integer *lda, const integer* ipiv,
492  real4 *work, const integer *lwork, integer *info);
493 extern "C"
494 void dgetri_(const integer* n, real8 *a, const integer *lda, const integer* ipiv,
495  real8 *work, const integer *lwork, integer *info);
496 extern "C"
497 void cgetri_(const integer* n, complex_real4 *a, const integer *lda, const integer* ipiv,
498  complex_real4 *work, const integer *lwork, integer *info);
499 extern "C"
500 void zgetri_(const integer* n, complex_real8 *a, const integer *lda, const integer* ipiv,
501  complex_real8 *work, const integer *lwork, integer *info);
502 
503 // SUBROUTINE DTRSM( SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB )
504 
505 // PURPOSE
506 // DTRSM solves one of the matrix equations
507 // op( A )*X = alpha*B, or X*op( A ) = alpha*B,
508 // where alpha is a scalar, X and B are m by n matrices, A is a unit, or
509 // non-unit, upper or lower triangular matrix and op( A ) is one of
510 // op( A ) = A or op( A ) = A**T.
511 
512 extern "C"
513 void strsm_(const char* side, const char* uplo, const char* transa, const char* diag,
514  const integer* m, const integer* n, const real4* alpha,
515  const real4* a, const integer* lda, real4* b, const integer* ldb,
516  char_len sidelen, char_len uplolen, char_len transalen, char_len diaglen);
517 extern "C"
518 void dtrsm_(const char* side, const char* uplo, const char* transa, const char* diag,
519  const integer* m, const integer* n, const real8* alpha,
520  const real8* a, const integer* lda, real8* b, const integer* ldb,
521  char_len sidelen, char_len uplolen, char_len transalen, char_len diaglen);
522 extern "C"
523 void ctrsm_(const char* side, const char* uplo, const char* transa, const char* diag,
524  const integer* m, const integer* n, const complex_real4* alpha,
525  const complex_real4* a, const integer* lda, complex_real4* b, const integer* ldb,
526  char_len sidelen, char_len uplolen, char_len transalen, char_len diaglen);
527 extern "C"
528 void ztrsm_(const char* side, const char* uplo, const char* transa, const char* diag,
529  const integer* m, const integer* n, const complex_real8* alpha,
530  const complex_real8* a, const integer* lda, complex_real8* b, const integer* ldb,
531  char_len sidelen, char_len uplolen, char_len transalen, char_len diaglen);
532 
533 // SUBROUTINE DTRTRI( UPLO, DIAG, N, A, LDA, INFO )
534 
535 // PURPOSE
536 // DTRTRI computes the inverse of a real upper or lower triangular
537 // matrix A.
538 
539 extern "C"
540 void strtri_(const char* uplo, const char* diag, const integer* n, const real4* a,
541  const integer* lda, integer *info);
542 extern "C"
543 void dtrtri_(const char* uplo, const char* diag, const integer* n, const real8* a,
544  const integer* lda, integer *info);
545 extern "C"
546 void ctrtri_(const char* uplo, const char* diag, const integer* n, const complex_real4* a,
547  const integer* lda, integer *info);
548 extern "C"
549 void ztrtri_(const char* uplo, const char* diag, const integer* n, const complex_real8* a,
550  const integer* lda, integer *info);
551 
552 #endif // MADNESS_LINALG_CLAPACK_FORTRAN_H__INCLUDED
double w(double t, double eps)
Definition: DKops.h:22
void sgeqrf_(integer *m, integer *n, real4 *a, integer *lda, real4 *tau, real4 *work, integer *lwork, integer *infoOUT)
void sgeev_(const char *jobz, const char *uplo, integer *n, real4 *a, integer *lda, real4 *w_real, real4 *w_imag, real4 *v, integer *ldv, real4 *vr, integer *ldvr, real4 *work, integer *lwork, integer *info, char_len jobzlen, char_len uplo_len)
void cgeqrf_(integer *m, integer *n, complex_real4 *a, integer *lda, complex_real4 *tau, complex_real4 *work, integer *lwork, integer *infoOUT)
void cgesvd_(const char *jobu, const char *jobvt, integer *m, integer *n, complex_real4 *a, integer *lda, real4 *s, complex_real4 *u, integer *ldu, complex_real4 *vt, integer *ldvt, complex_real4 *work, integer *lwork, real4 *rwork, integer *info, char_len jobulen, char_len jobvtlen)
void dgelss_(integer *m, integer *n, integer *nrhs, real8 *a, integer *lda, real8 *b, integer *ldb, real8 *sOUT, real8 *rcondIN, integer *rankOUT, real8 *work, integer *lwork, integer *infoOUT)
void sgeqp3_(integer *m, integer *n, real4 *a, integer *lda, integer *jpvt, real4 *tau, real4 *work, integer *lwork, integer *infoOUT)
void dgels_(const char *trans, integer *m, integer *n, integer *nrhs, real8 *a, integer *lda, real8 *b, integer *ldb, real8 *work, integer *lwork, integer *infoOUT, char_len translen)
void zgeev_(const char *jobz, const char *uplo, integer *n, complex_real8 *a, integer *lda, complex_real8 *w, complex_real8 *vl, integer *ldvl, complex_real8 *vr, integer *ldvr, complex_real8 *work, integer *lwork, real8 *rwork, integer *info, char_len jobzlen, char_len uplo_len)
void dsyev_(const char *jobz, const char *uplo, integer *n, real8 *a, integer *lda, real8 *w, real8 *work, integer *lwork, integer *info, char_len jobzlen, char_len uplo_len)
void dsygv_(integer *itype, const char *jobz, const char *uplo, integer *n, real8 *a, integer *lda, real8 *b, integer *ldb, real8 *w, real8 *work, integer *lwork, integer *info, char_len jobzlen, char_len uplo_len)
void zggev_(const char *jobz, const char *uplo, integer *n, complex_real8 *a, integer *lda, complex_real8 *b, integer *ldb, complex_real8 *alpha, complex_real8 *beta, complex_real8 *vl, integer *ldvl, complex_real8 *vr, integer *ldvr, complex_real8 *work, integer *lwork, real8 *rwork, integer *info, char_len jobzlen, char_len uplo_len)
void ctrtri_(const char *uplo, const char *diag, const integer *n, const complex_real4 *a, const integer *lda, integer *info)
void spstrf_(const char *uplo, const integer *n, real4 *a, const integer *lda, integer *ipiv, integer *rank, real4 *tol, real4 *work, integer *info)
void cungqr_(integer *m, integer *n, integer *k, complex_real4 *a, integer *lda, complex_real4 *tau, complex_real4 *work, integer *lwork, integer *info)
void cpotrf_(const char *uplo, const integer *n, complex_real4 *a, const integer *lda, integer *info, char_len uplo_len)
void sgelss_(integer *m, integer *n, integer *nrhs, real4 *a, integer *lda, real4 *b, integer *ldb, real4 *sOUT, real4 *rcondIN, integer *rankOUT, real4 *work, integer *lwork, integer *infoOUT)
void dgesv_(integer *n, integer *nrhs, real8 *AT, integer *lda, integer *piv, real8 *x, integer *ldx, integer *info)
double dlamch_(const char *mode, int modelen)
void ctrsm_(const char *side, const char *uplo, const char *transa, const char *diag, const integer *m, const integer *n, const complex_real4 *alpha, const complex_real4 *a, const integer *lda, complex_real4 *b, const integer *ldb, char_len sidelen, char_len uplolen, char_len transalen, char_len diaglen)
void sggev_(const char *jobz, const char *uplo, integer *n, real4 *a, integer *lda, real4 *b, integer *ldb, real4 *alphar, real4 *alphai, real4 *beta, real4 *vl, integer *ldvl, real4 *vr, integer *ldvr, real4 *work, integer *lwork, integer *info, char_len jobzlen, char_len uplo_len)
void zgesv_(integer *n, integer *nrhs, complex_real8 *AT, integer *lda, integer *piv, complex_real8 *x, integer *ldx, integer *info)
void sgesvd_(const char *jobu, const char *jobvt, integer *m, integer *n, real4 *a, integer *lda, real4 *s, real4 *u, integer *ldu, real4 *vt, integer *ldvt, real4 *work, integer *lwork, integer *info, char_len jobulen, char_len jobvtlen)
void zpotrf_(const char *uplo, const integer *n, complex_real8 *a, const integer *lda, integer *info, char_len uplo_len)
void cgetrf_(const integer *m, const integer *n, complex_real4 *a, const integer *lda, integer *ipiv, integer *info)
void zungqr_(integer *m, integer *n, integer *k, complex_real8 *a, integer *lda, complex_real8 *tau, complex_real8 *work, integer *lwork, integer *info)
void cgelss_(integer *m, integer *n, integer *nrhs, complex_real4 *a, integer *lda, complex_real4 *b, integer *ldb, real4 *sOUT, real4 *rcondIN, integer *rankOUT, complex_real4 *work, integer *lwork, real4 *rwork, integer *infoOUT)
void zgetri_(const integer *n, complex_real8 *a, const integer *lda, const integer *ipiv, complex_real8 *work, const integer *lwork, integer *info)
void zhegv_(integer *itype, const char *jobz, const char *uplo, integer *n, complex_real8 *a, integer *lda, complex_real8 *b, integer *ldb, real8 *w, complex_real8 *work, integer *lwork, real8 *rwork, integer *info, char_len jobzlen, char_len uplo_len)
void ssygv_(integer *itype, const char *jobz, const char *uplo, integer *n, real4 *a, integer *lda, real4 *b, integer *ldb, real4 *w, real4 *work, integer *lwork, integer *info, char_len jobzlen, char_len uplo_len)
float slamch_(const char *mode, int modelen)
void cgetri_(const integer *n, complex_real4 *a, const integer *lda, const integer *ipiv, complex_real4 *work, const integer *lwork, integer *info)
void sgetrf_(const integer *m, const integer *n, real4 *a, const integer *lda, integer *ipiv, integer *info)
void dtrsm_(const char *side, const char *uplo, const char *transa, const char *diag, const integer *m, const integer *n, const real8 *alpha, const real8 *a, const integer *lda, real8 *b, const integer *ldb, char_len sidelen, char_len uplolen, char_len transalen, char_len diaglen)
void cheev_(const char *jobz, const char *uplo, integer *n, complex_real4 *a, integer *lda, real4 *w, complex_real4 *work, integer *lwork, real4 *rwork, integer *info, char_len jobzlen, char_len uplo_len)
void dgeqrf_(integer *m, integer *n, real8 *a, integer *lda, real8 *tau, real8 *work, integer *lwork, integer *infoOUT)
void cggev_(const char *jobz, const char *uplo, integer *n, complex_real4 *a, integer *lda, complex_real4 *b, integer *ldb, complex_real4 *alpha, complex_real4 *beta, complex_real4 *vl, integer *ldvl, complex_real4 *vr, integer *ldvr, complex_real4 *work, integer *lwork, real4 *rwork, integer *info, char_len jobzlen, char_len uplo_len)
void zgeqrf_(integer *m, integer *n, complex_real8 *a, integer *lda, complex_real8 *tau, complex_real8 *work, integer *lwork, integer *infoOUT)
void cgeev_(const char *jobz, const char *uplo, integer *n, complex_real4 *a, integer *lda, complex_real4 *w, complex_real4 *vl, integer *ldvl, complex_real4 *vr, integer *ldvr, complex_real4 *work, integer *lwork, real4 *rwork, integer *info, char_len jobzlen, char_len uplo_len)
void dgetri_(const integer *n, real8 *a, const integer *lda, const integer *ipiv, real8 *work, const integer *lwork, integer *info)
void zgesvd_(const char *jobu, const char *jobvt, integer *m, integer *n, complex_real8 *a, integer *lda, real8 *s, complex_real8 *u, integer *ldu, complex_real8 *vt, integer *ldvt, complex_real8 *work, integer *lwork, real8 *rwork, integer *info, char_len jobulen, char_len jobvtlen)
void dorgqr_(integer *m, integer *n, integer *k, real8 *a, integer *lda, real8 *tau, real8 *work, integer *lwork, integer *info)
void dpstrf_(const char *uplo, const integer *n, real8 *a, const integer *lda, integer *ipiv, integer *rank, real8 *tol, real8 *work, integer *info)
void dgeqp3_(integer *m, integer *n, real8 *a, integer *lda, integer *jpvt, real8 *tau, real8 *work, integer *lwork, integer *infoOUT)
void dtrtri_(const char *uplo, const char *diag, const integer *n, const real8 *a, const integer *lda, integer *info)
void zheev_(const char *jobz, const char *uplo, integer *n, complex_real8 *a, integer *lda, real8 *w, complex_real8 *work, integer *lwork, real8 *rwork, integer *info, char_len jobzlen, char_len uplo_len)
void zpstrf_(const char *uplo, const integer *n, complex_real8 *a, const integer *lda, integer *ipiv, integer *rank, real8 *tol, complex_real8 *work, integer *info)
void dgesvd_(const char *jobu, const char *jobvt, integer *m, integer *n, real8 *a, integer *lda, real8 *s, real8 *u, integer *ldu, real8 *vt, integer *ldvt, real8 *work, integer *lwork, integer *info, char_len jobulen, char_len jobvtlen)
void cpstrf_(const char *uplo, const integer *n, complex_real4 *a, const integer *lda, integer *ipiv, integer *rank, real4 *tol, complex_real4 *work, integer *info)
void ztrsm_(const char *side, const char *uplo, const char *transa, const char *diag, const integer *m, const integer *n, const complex_real8 *alpha, const complex_real8 *a, const integer *lda, complex_real8 *b, const integer *ldb, char_len sidelen, char_len uplolen, char_len transalen, char_len diaglen)
void sgetri_(const integer *n, real4 *a, const integer *lda, const integer *ipiv, real4 *work, const integer *lwork, integer *info)
void spotrf_(const char *uplo, const integer *n, real4 *a, const integer *lda, integer *info, char_len uplo_len)
void zgels_(const char *trans, integer *m, integer *n, integer *nrhs, complex_real8 *a, integer *lda, complex_real8 *b, integer *ldb, complex_real8 *work, integer *lwork, real8 *rwork, integer *infoOUT, char_len translen)
void strtri_(const char *uplo, const char *diag, const integer *n, const real4 *a, const integer *lda, integer *info)
void cgels_(const char *trans, integer *m, integer *n, integer *nrhs, complex_real4 *a, integer *lda, complex_real4 *b, integer *ldb, complex_real4 *work, integer *lwork, real4 *rwork, integer *infoOUT, char_len translen)
void zgelss_(integer *m, integer *n, integer *nrhs, complex_real8 *a, integer *lda, complex_real8 *b, integer *ldb, real8 *sOUT, real8 *rcondIN, integer *rankOUT, complex_real8 *work, integer *lwork, real8 *rwork, integer *infoOUT)
void dpotrf_(const char *uplo, const integer *n, real8 *a, const integer *lda, integer *info, char_len uplo_len)
void sgesv_(integer *n, integer *nrhs, real4 *AT, integer *lda, integer *piv, real4 *x, integer *ldx, integer *info)
void cgeqp3_(integer *m, integer *n, complex_real4 *a, integer *lda, integer *jpvt, complex_real4 *tau, complex_real4 *work, integer *lwork, real4 *rwork, integer *infoOUT)
void cgesv_(integer *n, integer *nrhs, complex_real4 *AT, integer *lda, integer *piv, complex_real4 *x, integer *ldx, integer *info)
void dgeev_(const char *jobz, const char *uplo, integer *n, real8 *a, integer *lda, real8 *w_real, real8 *w_imag, real8 *v, integer *ldv, real8 *vr, integer *ldvr, real8 *work, integer *lwork, integer *info, char_len jobzlen, char_len uplo_len)
void zgeqp3_(integer *m, integer *n, complex_real8 *a, integer *lda, integer *jpvt, complex_real8 *tau, complex_real8 *work, integer *lwork, real8 *rwork, integer *infoOUT)
void ssyev_(const char *jobz, const char *uplo, integer *n, real4 *a, integer *lda, real4 *w, real4 *work, integer *lwork, integer *info, char_len jobzlen, char_len uplo_len)
void sorgqr_(integer *m, integer *n, integer *k, real4 *a, integer *lda, real4 *tau, real4 *work, integer *lwork, integer *info)
void zgetrf_(const integer *m, const integer *n, complex_real8 *a, const integer *lda, integer *ipiv, integer *info)
void strsm_(const char *side, const char *uplo, const char *transa, const char *diag, const integer *m, const integer *n, const real4 *alpha, const real4 *a, const integer *lda, real4 *b, const integer *ldb, char_len sidelen, char_len uplolen, char_len transalen, char_len diaglen)
void chegv_(integer *itype, const char *jobz, const char *uplo, integer *n, complex_real4 *a, integer *lda, complex_real4 *b, integer *ldb, real4 *w, complex_real4 *work, integer *lwork, real4 *rwork, integer *info, char_len jobzlen, char_len uplo_len)
void dgetrf_(const integer *m, const integer *n, real8 *a, const integer *lda, integer *ipiv, integer *info)
void sgels_(const char *trans, integer *m, integer *n, integer *nrhs, real4 *a, integer *lda, real4 *b, integer *ldb, real4 *work, integer *lwork, integer *infoOUT, char_len translen)
void ztrtri_(const char *uplo, const char *diag, const integer *n, const complex_real8 *a, const integer *lda, integer *info)
void dggev_(const char *jobl, const char *jobr, integer *n, real8 *a, integer *lda, real8 *b, integer *ldb, real8 *w_real, real8 *w_imag, real8 *beta, real8 *vl, integer *ldvl, real8 *vr, integer *ldvr, real8 *work, integer *lwork, integer *info, char_len jobzlen, char_len uplo_len)
int integer
Definition: crayio.c:25
Correspondence between C++ and Fortran types.
std::complex< double > complex_real8
Fortran double complex.
Definition: fortran_ctypes.h:83
std::complex< float > complex_real4
Fortran single complex.
Definition: fortran_ctypes.h:88
int char_len
Type of variable appended to argument list for length of fortran character strings.
Definition: fortran_ctypes.h:93
double real8
Fortran double precision.
Definition: fortran_ctypes.h:73
float real4
Fortran single precision.
Definition: fortran_ctypes.h:78
const double m
Definition: gfit.cc:199
const double beta
Definition: gygi_soltion.cc:62
static const double v
Definition: hatom_sf_dirac.cc:20
static const double b
Definition: nonlinschro.cc:119
static const double a
Definition: nonlinschro.cc:118
static const long k
Definition: rk.cc:44
const double alpha
Definition: test_coulomb.cc:54
double u(const double x, const double expnt)
Definition: testperiodic.cc:56