Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Lapack.hpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File Lapack.hpp
4 //
5 // For more information, please see: http://www.nektar.info
6 //
7 // The MIT License
8 //
9 // Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10 // Department of Aeronautics, Imperial College London (UK), and Scientific
11 // Computing and Imaging Institute, University of Utah (USA).
12 //
13 // License for the specific language governing rights and limitations under
14 // Permission is hereby granted, free of charge, to any person obtaining a
15 // copy of this software and associated documentation files (the "Software"),
16 // to deal in the Software without restriction, including without limitation
17 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
18 // and/or sell copies of the Software, and to permit persons to whom the
19 // Software is furnished to do so, subject to the following conditions:
20 //
21 // The above copyright notice and this permission notice shall be included
22 // in all copies or substantial portions of the Software.
23 //
24 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
25 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
26 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
27 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
28 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
29 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
30 // DEALINGS IN THE SOFTWARE.
31 //
32 // Description: wrapper of functions around standard LAPACK routines
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 #ifndef NEKTAR_LIB_UTILITIES_LAPACK_HPP
36 #define NEKTAR_LIB_UTILITIES_LAPACK_HPP
37 
39 
40 namespace Lapack
41 {
42  extern "C"
43  {
44  // Matrix factorisation and solves
45  void F77NAME(dsptrf) (const char& uplo, const int& n,
46  double* ap, int *ipiv, int& info);
47  void F77NAME(dsptrs) (const char& uplo, const int& n,
48  const int& nrhs, const double* ap,
49  const int *ipiv, double* b,
50  const int& ldb, int& info);
51  void F77NAME(dsptri) (const char& uplo, const int& n,
52  const double* ap,
53  const int *ipiv, double* work,
54  int& info);
55  void F77NAME(dtrtrs) (const char& uplo, const char& trans, const char& diag,
56  const int& n, const int& nrhs, const double* a,
57  const int& lda, double* b, const int& ldb, int& info);
58  void F77NAME(dtptrs) (const char& uplo, const char& trans, const char& diag,
59  const int& n, const int& nrhs, const double* a,
60  double* b, const int& ldb, int& info);
61  void F77NAME(dpptrf) (const char& uplo, const int& n,
62  double* ap, int& info);
63  void F77NAME(dpptrs) (const char& uplo, const int& n,
64  const int& nrhs, const double* ap,
65  double* b, const int& ldb, int& info);
66  void F77NAME(dpbtrf) (const char& uplo, const int& n, const int& kd,
67  double* ab, const int& ldab, int& info);
68  void F77NAME(dpbtrs) (const char& uplo, const int& n,
69  const int& kd, const int& nrhs,
70  const double* ab, const int& ldab,
71  double* b, const int& ldb, int& info);
72  void F77NAME(dgbtrf) (const int& m, const int& n, const int& kl,
73  const int& ku, double* a, const int& lda,
74  int* ipiv, int& info);
75  void F77NAME(dgbtrs) (const char& trans, const int& n, const int& kl,
76  const int &ku, const int& nrhs, const double* a,
77  const int& lda, const int* ipiv, double* b,
78  const int& ldb, int& info);
79  void F77NAME(dgetrf) (const int& m, const int& n, double* a,
80  const int& lda, int* ipiv, int& info);
81  void F77NAME(dgetrs) (const char& trans, const int& n, const int& nrhs,
82  const double* a, const int& lda, int* ipiv,
83  double* b, const int& ldb, int& info);
84  void F77NAME(dgetri) (const int& n, double *a, const int& lda,
85  const int *ipiv, double *wk, const int& lwk,
86  int& info);
87  void F77NAME(dsterf) (const int& n, double *d, double *e, int& info);
88  void F77NAME(dgeev) (const char& uplo, const char& lrev, const int& n,
89  const double* a, const int& lda, double* wr, double* wi,
90  double* rev, const int& ldr,
91  double* lev, const int& ldv,
92  double* work, const int& lwork, int& info);
93 
94  void F77NAME(dspev) (const char& jobz, const char& uplo, const int& n,
95  double* ap, double* w, double* z, const int& ldz,
96  double* work, int& info);
97  void F77NAME(dsbev) (const char& jobz, const char& uplo, const int& kl,
98  const int& ku, double* ap, const int& lda,
99  double* w, double* z, const int& ldz,
100  double* work, int& info);
101  }
102 
103  // Non-standard versions.
104  void dgetrs(char trans, int matrixRows, int matrixColumns, const double* A, double* x);
105 
106 #ifdef NEKTAR_USING_LAPACK
107 
108  /// \brief factor a real packed-symmetric matrix using Bunch-Kaufman
109  /// pivoting.
110  static inline void Dsptrf (const char& uplo, const int& n,
111  double* ap, int *ipiv, int& info)
112  {
113  F77NAME(dsptrf) (uplo,n,ap,ipiv,info);
114  }
115 
116  /// \brief Solve a real symmetric matrix problem using Bunch-Kaufman
117  /// pivoting.
118  static inline void Dsptrs (const char& uplo, const int& n, const int& nrhs,
119  const double* ap, const int *ipiv, double* b,
120  const int& ldb, int& info)
121  {
122  F77NAME(dsptrs) (uplo,n,nrhs,ap,ipiv,b,ldb,info);
123  }
124 
125  /// \brief Invert a real symmetric matrix problem
126  static inline void Dsptri (const char& uplo, const int& n,
127  const double* ap, const int *ipiv, double* work,
128  int& info)
129  {
130  F77NAME(dsptri) (uplo,n,ap,ipiv,work,info);
131  }
132 
133  /// \brief Cholesky factor a real Positive Definite packed-symmetric matrix.
134  static inline void Dpptrf (const char& uplo, const int& n,
135  double *ap, int& info)
136  {
137  F77NAME(dpptrf) (uplo,n,ap,info);
138  }
139 
140  /// \brief Solve a real Positive defiinte symmetric matrix problem
141  /// using Cholesky factorization.
142  static inline void Dpptrs (const char& uplo, const int& n, const int& nrhs,
143  const double *ap, double *b, const int& ldb,
144  int& info)
145  {
146  F77NAME(dpptrs) (uplo,n,nrhs,ap,b,ldb,info);
147  }
148 
149  /// \brief Cholesky factorize a real positive-definite
150  /// banded-symmetric matrix
151  static inline void Dpbtrf (const char& uplo, const int& n, const int& kd,
152  double *ab, const int& ldab, int& info)
153  {
154  F77NAME(dpbtrf) (uplo,n,kd,ab,ldab,info);
155  }
156 
157 
158  /// \brief Solve a real, Positive definite banded symmetric matrix
159  /// problem using Cholesky factorization.
160  static inline void Dpbtrs (const char& uplo, const int& n,
161  const int& kd, const int& nrhs,
162  const double *ab, const int& ldab,
163  double *b, const int& ldb, int& info)
164  {
165  F77NAME(dpbtrs) (uplo,n,kd,nrhs,ab,ldab,b,ldb,info);
166  }
167 
168  /// \brief General banded matrix LU factorisation
169  static inline void Dgbtrf (const int& m, const int& n, const int& kl,
170  const int& ku, double* a, const int& lda,
171  int* ipiv, int& info)
172  {
173  F77NAME(dgbtrf)(m,n,kl,ku,a,lda,ipiv,info);
174  }
175 
176  /// \brief Solve general banded matrix using LU factorisation
177  static inline void Dgbtrs (const char& trans, const int& n, const int& kl,
178  const int &ku, const int& nrhs, const double* a,
179  const int& lda, const int* ipiv, double* b,
180  const int& ldb, int& info)
181  {
182  F77NAME(dgbtrs)(trans,n,kl,ku,nrhs,a,lda,ipiv,b,ldb,info);
183  }
184 
185  /// \brief General matrix LU factorisation
186  static inline void Dgetrf (const int& m, const int& n, double *a,
187  const int& lda, int *ipiv, int& info)
188  {
189  F77NAME(dgetrf) (m,n,a,lda,ipiv,info);
190  }
191 
192  /// \brief General matrix LU backsolve
193  static inline void Dgetrs (const char& trans, const int& n, const int& nrhs,
194  const double* a, const int& lda, int* ipiv,
195  double* b, const int& ldb, int& info)
196  {
197  F77NAME(dgetrs) (trans,n,nrhs,a,lda,ipiv,b,ldb,info);
198  }
199 
200  /// \brief Generate matrix inverse
201  static inline void Dgetri (const int& n, double *a, const int& lda,
202  const int *ipiv, double *wk, const int& lwk, int& info)
203  {
204  F77NAME(dgetri) (n, a, lda, ipiv, wk, lwk,info);
205  }
206 
207  /// \brief Find eigenvalues of symmetric tridiagonal matrix
208  static inline void Dsterf(const int& n, double *d, double *e, int& info)
209  {
210  F77NAME(dsterf)(n,d,e,info);
211  }
212 
213  /// \brief Solve general real matrix eigenproblem.
214  static inline void Dgeev (const char& uplo, const char& lrev, const int& n,
215  const double* a, const int& lda, double* wr, double* wi,
216  double* rev, const int& ldr,
217  double* lev, const int& ldv,
218  double* work, const int& lwork, int& info)
219  {
220  F77NAME(dgeev) (uplo, lrev, n, a, lda, wr, wi, rev,
221  ldr, lev, ldv, work, lwork, info);
222  }
223 
224  /// \brief Solve packed-symmetric real matrix eigenproblem.
225  static inline void Dspev (const char& jobz, const char& uplo, const int& n,
226  double* ap, double* w, double* z, const int& ldz,
227  double* work, int& info)
228  {
229  F77NAME(dspev) (jobz, uplo, n, ap, w, z, ldz, work, info);
230  }
231 
232  /// \brief Solve packed-banded real matrix eigenproblem.
233  static inline void Dsbev (const char& jobz, const char& uplo, const int& kl,
234  const int& ku, double* ap, const int& lda,
235  double* w, double* z, const int& ldz,
236  double* work, int& info)
237  {
238  F77NAME(dsbev) (jobz, uplo, kl, ku, ap, lda, w, z, ldz, work, info);
239  }
240 
241  static inline void Dtrtrs(const char& uplo, const char& trans, const char& diag,
242  const int& n, const int& nrhs, const double* a,
243  const int& lda, double* b, const int& ldb, int& info)
244  {
245  F77NAME(dtrtrs) (uplo, trans, diag, n, nrhs, a, lda, b, ldb, info);
246  }
247 
248  static inline void Dtptrs(const char& uplo, const char& trans, const char& diag,
249  const int& n, const int& nrhs, const double* a,
250  double* b, const int& ldb, int& info)
251  {
252  F77NAME(dtptrs) (uplo, trans, diag, n, nrhs, a, b, ldb, info);
253  }
254 
255 
256 #endif //NEKTAR_USING_LAPACK
257 }
258 #endif //NEKTAR_LIB_UTILITIES_LAPACK_HPP
259 
260 /***
261 $Log: Lapack.hpp,v $
262 Revision 1.5 2008/06/01 19:04:57 bnelson
263 Added triangular solvers.
264 
265 Revision 1.4 2008/04/30 02:57:15 bnelson
266 Fixed gcc compiler warning.
267 
268 Revision 1.3 2008/04/06 05:55:11 bnelson
269 Changed ConstArray to Array<const>
270 
271 Revision 1.2 2007/04/10 14:00:45 sherwin
272 Update to include SharedArray in all 2D element (including Nodal tris). Have also remvoed all new and double from 2D shapes in StdRegions
273 
274 Revision 1.1 2007/04/03 03:59:24 bnelson
275 Moved Lapack.hpp, Blas.hpp, Transf77.hpp to LinearAlgebra
276 
277 Revision 1.3 2007/02/04 00:15:40 bnelson
278 *** empty log message ***
279 
280 Revision 1.2 2006/06/01 13:44:28 kirby
281 *** empty log message ***
282 
283 Revision 1.1 2006/06/01 11:07:52 kirby
284 *** empty log message ***
285 
286 Revision 1.1 2006/05/04 18:57:43 kirby
287 *** empty log message ***
288 
289 Revision 1.4 2006/02/26 21:13:45 bnelson
290 Fixed a variety of compiler errors caused by updates to the coding standard.
291 
292 Revision 1.3 2006/02/15 08:07:15 sherwin
293 
294 Put codes into standard although have not yet been compiled
295 
296 Revision 1.2 2006/02/12 21:51:42 sherwin
297 
298 Added licence
299 
300 Revision 1.1 2006/02/12 15:06:12 sherwin
301 
302 Changed .h files to .hpp
303 
304 **/
void F77NAME() dsptrs(const char &uplo, const int &n, const int &nrhs, const double *ap, const int *ipiv, double *b, const int &ldb, int &info)
void F77NAME() dpbtrf(const char &uplo, const int &n, const int &kd, double *ab, const int &ldab, int &info)
void F77NAME() dgetrf(const int &m, const int &n, double *a, const int &lda, int *ipiv, int &info)
void F77NAME() dtrtrs(const char &uplo, const char &trans, const char &diag, const int &n, const int &nrhs, const double *a, const int &lda, double *b, const int &ldb, int &info)
#define F77NAME(x)
Fortran routines need an underscore.
Definition: TransF77.hpp:47
void F77NAME() dpbtrs(const char &uplo, const int &n, const int &kd, const int &nrhs, const double *ab, const int &ldab, double *b, const int &ldb, int &info)
void F77NAME() dsterf(const int &n, double *d, double *e, int &info)
void F77NAME() dtptrs(const char &uplo, const char &trans, const char &diag, const int &n, const int &nrhs, const double *a, double *b, const int &ldb, int &info)
void F77NAME() dpptrf(const char &uplo, const int &n, double *ap, int &info)
void F77NAME() dgeev(const char &uplo, const char &lrev, const int &n, const double *a, const int &lda, double *wr, double *wi, double *rev, const int &ldr, double *lev, const int &ldv, double *work, const int &lwork, int &info)
void F77NAME() dgetrs(const char &trans, const int &n, const int &nrhs, const double *a, const int &lda, int *ipiv, double *b, const int &ldb, int &info)
void F77NAME() dspev(const char &jobz, const char &uplo, const int &n, double *ap, double *w, double *z, const int &ldz, double *work, int &info)
void F77NAME() dsptri(const char &uplo, const int &n, const double *ap, const int *ipiv, double *work, int &info)
void F77NAME() dsbev(const char &jobz, const char &uplo, const int &kl, const int &ku, double *ap, const int &lda, double *w, double *z, const int &ldz, double *work, int &info)
void F77NAME() dgbtrf(const int &m, const int &n, const int &kl, const int &ku, double *a, const int &lda, int *ipiv, int &info)
void F77NAME() dpptrs(const char &uplo, const int &n, const int &nrhs, const double *ap, double *b, const int &ldb, int &info)
void F77NAME() dgetri(const int &n, double *a, const int &lda, const int *ipiv, double *wk, const int &lwk, int &info)
void F77NAME() dsptrf(const char &uplo, const int &n, double *ap, int *ipiv, int &info)
void F77NAME() dgbtrs(const char &trans, const int &n, const int &kl, const int &ku, const int &nrhs, const double *a, const int &lda, const int *ipiv, double *b, const int &ldb, int &info)