GetFEM  5.5
gmm_lapack_interface.h
Go to the documentation of this file.
1 /* -*- c++ -*- (enables emacs c++ mode) */
2 /*===========================================================================
3 
4  Copyright (C) 2003-2026 Yves Renard
5 
6  This file is a part of GetFEM
7 
8  GetFEM is free software; you can redistribute it and/or modify it
9  under the terms of the GNU Lesser General Public License as published
10  by the Free Software Foundation; either version 3 of the License, or
11  (at your option) any later version along with the GCC Runtime Library
12  Exception either version 3.1 or (at your option) any later version.
13  This program is distributed in the hope that it will be useful, but
14  WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
15  or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
16  License and GCC Runtime Library Exception for more details.
17  You should have received a copy of the GNU Lesser General Public License
18  along with this program. If not, see https://www.gnu.org/licenses/.
19 
20  As a special exception, you may use this file as it is a part of a free
21  software library without restriction. Specifically, if other files
22  instantiate templates or use macros or inline functions from this file,
23  or you compile this file and link it with other files to produce an
24  executable, this file does not by itself cause the resulting executable
25  to be covered by the GNU Lesser General Public License. This exception
26  does not however invalidate any other reasons why the executable file
27  might be covered by the GNU Lesser General Public License.
28 
29 ===========================================================================*/
30 
31 /**@file gmm_lapack_interface.h
32  @author Yves Renard <[email protected]>
33  @date October 7, 2003.
34  @brief gmm interface for LAPACK
35 */
36 
37 #ifndef GMM_LAPACK_INTERFACE_H
38 #define GMM_LAPACK_INTERFACE_H
39 
40 #include "gmm_blas_interface.h"
41 #include "gmm_dense_lu.h"
42 #include "gmm_dense_qr.h"
43 
44 #if defined(GMM_USES_LAPACK) && !defined(GMM_MATLAB_INTERFACE)
45 
46 namespace gmm {
47 
48  /* ********************************************************************** */
49  /* Operations interfaced for T = float, double, std::complex<float> */
50  /* or std::complex<double> : */
51  /* */
52  /* lu_factor(dense_matrix<T>, std::vector<long>) */
53  /* lu_solve(dense_matrix<T>, std::vector<T>, std::vector<T>) */
54  /* lu_solve(dense_matrix<T>, std::vector<long>, std::vector<T>, */
55  /* std::vector<T>) */
56  /* lu_solve_transposed(dense_matrix<T>, std::vector<long>, std::vector<T>,*/
57  /* std::vector<T>) */
58  /* lu_inverse(dense_matrix<T>) */
59  /* lu_inverse(dense_matrix<T>, std::vector<long>, dense_matrix<T>) */
60  /* */
61  /* qr_factor(dense_matrix<T>, dense_matrix<T>, dense_matrix<T>) */
62  /* */
63  /* implicit_qr_algorithm(dense_matrix<T>, std::vector<T>) */
64  /* implicit_qr_algorithm(dense_matrix<T>, std::vector<T>, */
65  /* dense_matrix<T>) */
66  /* implicit_qr_algorithm(dense_matrix<T>, std::vector<std::complex<T> >) */
67  /* implicit_qr_algorithm(dense_matrix<T>, std::vector<std::complex<T> >, */
68  /* dense_matrix<T>) */
69  /* */
70  /* geev_interface_right */
71  /* geev_interface_left */
72  /* */
73  /* schur(dense_matrix<T>, dense_matrix<T>, dense_matrix<T>) */
74  /* */
75  /* svd(dense_matrix<T>, dense_matrix<T>, dense_matrix<T>, std::vector<T>) */
76  /* svd(dense_matrix<T>, dense_matrix<T>, dense_matrix<T>, */
77  /* std::vector<std::complex<T> >) */
78  /* */
79  /* ********************************************************************** */
80 
81  /* ********************************************************************** */
82  /* LAPACK functions used. */
83  /* ********************************************************************** */
84 
85  extern "C" {
86  void sgetrf_(...); void dgetrf_(...); void cgetrf_(...); void zgetrf_(...);
87  void sgetrs_(...); void dgetrs_(...); void cgetrs_(...); void zgetrs_(...);
88  void sgetri_(...); void dgetri_(...); void cgetri_(...); void zgetri_(...);
89  void sgeqrf_(...); void dgeqrf_(...); void cgeqrf_(...); void zgeqrf_(...);
90  void sorgqr_(...); void dorgqr_(...); void cungqr_(...); void zungqr_(...);
91  void sormqr_(...); void dormqr_(...); void cunmqr_(...); void zunmqr_(...);
92  void sgees_ (...); void dgees_ (...); void cgees_ (...); void zgees_ (...);
93  void sgeev_ (...); void dgeev_ (...); void cgeev_ (...); void zgeev_ (...);
94  void sgeesx_(...); void dgeesx_(...); void cgeesx_(...); void zgeesx_(...);
95  void sgesvd_(...); void dgesvd_(...); void cgesvd_(...); void zgesvd_(...);
96  }
97 
98  /* ********************************************************************** */
99  /* LU decomposition. */
100  /* ********************************************************************** */
101 
102 # define getrf_interface(lapack_name, base_type) \
103  inline size_type \
104  lu_factor(dense_matrix<base_type> &A, lapack_ipvt &ipvt) { \
105  GMMLAPACK_TRACE("getrf_interface"); \
106  const BLAS_INT m=BLAS_INT(mat_nrows(A)), n=BLAS_INT(mat_ncols(A)), lda(m);\
107  BLAS_INT info(-1); \
108  if (m && n) lapack_name(&m, &n, &A(0,0), &lda, &ipvt[0], &info); \
109  return size_type(abs(info)); \
110  }
111 
112  getrf_interface(sgetrf_, BLAS_S)
113  getrf_interface(dgetrf_, BLAS_D)
114  getrf_interface(cgetrf_, BLAS_C)
115  getrf_interface(zgetrf_, BLAS_Z)
116 
117  /* ********************************************************************* */
118  /* LU solve. */
119  /* ********************************************************************* */
120 
121 # define getrs_interface(f_name, NorT, lapack_name, base_type) \
122  inline void \
123  f_name(const dense_matrix<base_type> &A, const lapack_ipvt &ipvt, \
124  std::vector<base_type> &x, const std::vector<base_type> &b) { \
125  GMMLAPACK_TRACE("getrs_interface"); \
126  const char t=NorT; const BLAS_INT n=BLAS_INT(mat_nrows(A)), nrhs(1); \
127  BLAS_INT info(0); gmm::copy(b, x); \
128  if (n) lapack_name(&t, &n, &nrhs, &A(0,0), &n, &ipvt[0], &x[0], &n, \
129  &info); \
130  }
131 
132  getrs_interface(lu_solve, 'N', sgetrs_, BLAS_S)
133  getrs_interface(lu_solve, 'N', dgetrs_, BLAS_D)
134  getrs_interface(lu_solve, 'N', cgetrs_, BLAS_C)
135  getrs_interface(lu_solve, 'N', zgetrs_, BLAS_Z)
136  getrs_interface(lu_solve_transposed, 'T', sgetrs_, BLAS_S)
137  getrs_interface(lu_solve_transposed, 'T', dgetrs_, BLAS_D)
138  getrs_interface(lu_solve_transposed, 'T', cgetrs_, BLAS_C)
139  getrs_interface(lu_solve_transposed, 'T', zgetrs_, BLAS_Z)
140 
141  /* ********************************************************************* */
142  /* LU inverse. */
143  /* ********************************************************************* */
144 
145 # define getri_interface(lapack_name, base_type) \
146  inline void lu_inverse(const dense_matrix<base_type> &LU, \
147  const lapack_ipvt &ipvt, \
148  dense_matrix<base_type> &A) { \
149  GMMLAPACK_TRACE("getri_interface"); \
150  const BLAS_INT n=BLAS_INT(mat_nrows(A)); \
151  BLAS_INT info(0), lwork(-1); base_type work1; \
152  if (n) { \
153  gmm::copy(LU, A); \
154  lapack_name(&n, &A(0,0), &n, &ipvt[0], &work1, &lwork, &info); \
155  const size_type worksize=size_type(gmm::real(work1)); \
156  std::vector<base_type> work(worksize); lwork=BLAS_INT(worksize); \
157  lapack_name(&n, &A(0,0), &n, &ipvt[0], &work[0], &lwork, &info); \
158  } \
159  }
160 
161  getri_interface(sgetri_, BLAS_S)
162  getri_interface(dgetri_, BLAS_D)
163  getri_interface(cgetri_, BLAS_C)
164  getri_interface(zgetri_, BLAS_Z)
165 
166  /* ********************************************************************** */
167  /* QR factorization. */
168  /* ********************************************************************** */
169 
170 # define geqrf_interface(lapack_name, base_type) \
171  inline void qr_factor(dense_matrix<base_type> &A) { \
172  GMMLAPACK_TRACE("geqrf_interface"); \
173  const size_type mm(mat_nrows(A)), nn(mat_ncols(A)); \
174  if (mm && nn) { \
175  const BLAS_INT m=BLAS_INT(mm), n=BLAS_INT(nn); \
176  BLAS_INT info(0), lwork(-1); \
177  std::vector<base_type> tau(nn); base_type work1; \
178  lapack_name(&m, &n, &A(0,0), &m, &tau[0], &work1, &lwork, &info); \
179  const size_type worksize=size_type(gmm::real(work1)); \
180  std::vector<base_type> work(worksize); lwork=BLAS_INT(worksize); \
181  lapack_name(&m, &n, &A(0,0), &m, &tau[0], &work[0], &lwork, &info); \
182  GMM_ASSERT1(!info, "QR factorization failed"); \
183  } \
184  }
185 
186  geqrf_interface(sgeqrf_, BLAS_S)
187  geqrf_interface(dgeqrf_, BLAS_D)
188  // For complex values, housholder vectors are not the same as in
189  // gmm::lu_factor. Impossible to interface for the moment.
190  // geqrf_interface(cgeqrf_, BLAS_C)
191  // geqrf_interface(zgeqrf_, BLAS_Z)
192 
193 # define geqrf_interface2(lapack_name1, lapack_name2, base_type) inline \
194  void qr_factor(const dense_matrix<base_type> &A, \
195  dense_matrix<base_type> &Q, dense_matrix<base_type> &R) { \
196  GMMLAPACK_TRACE("geqrf_interface2"); \
197  const size_type mm(mat_nrows(A)), nn(mat_ncols(A)); \
198  if (mm && nn) { \
199  const BLAS_INT m=BLAS_INT(mm), n=BLAS_INT(nn); \
200  BLAS_INT info(0), lwork(-1); \
201  std::copy(A.begin(), A.end(), Q.begin()); \
202  std::vector<base_type> tau(nn); base_type work1; \
203  lapack_name1(&m, &n, &Q(0,0), &m, &tau[0], &work1, &lwork, &info); \
204  const size_type worksize=size_type(gmm::real(work1)); \
205  std::vector<base_type> work(worksize); lwork=BLAS_INT(worksize); \
206  lapack_name1(&m, &n, &Q(0,0), &m, &tau[0], &work[0], &lwork, &info); \
207  GMM_ASSERT1(!info, "QR factorization failed"); \
208  base_type *p = &R(0,0), *q = &Q(0,0); \
209  for (BLAS_INT j = 0; j < n; ++j, q += m-n) \
210  for (BLAS_INT i = 0; i < n; ++i, ++p, ++q) \
211  *p = (j < i) ? base_type(0) : *q; \
212  lapack_name2(&m, &n, &n, &Q(0,0), &m,&tau[0],&work[0],&lwork,&info); \
213  } \
214  else gmm::clear(Q); \
215  }
216 
217  geqrf_interface2(sgeqrf_, sorgqr_, BLAS_S)
218  geqrf_interface2(dgeqrf_, dorgqr_, BLAS_D)
219  geqrf_interface2(cgeqrf_, cungqr_, BLAS_C)
220  geqrf_interface2(zgeqrf_, zungqr_, BLAS_Z)
221 
222  /* ********************************************************************** */
223  /* QR algorithm for eigenvalues search. */
224  /* ********************************************************************** */
225 
226 # define gees_interface(lapack_name, base_type) \
227  template <typename VECT> inline void implicit_qr_algorithm( \
228  const dense_matrix<base_type> &A, VECT &eigval_, \
229  dense_matrix<base_type> &Q, \
230  double tol=gmm::default_tol(base_type()), bool compvect = true) { \
231  GMMLAPACK_TRACE("gees_interface"); \
232  typedef bool (*L_fp)(...); L_fp p = 0; \
233  const size_type nn(mat_nrows(A)); \
234  if (!nn) return; \
235  dense_matrix<base_type> H(nn,nn); gmm::copy(A, H); \
236  const char jobvs=(compvect ? 'V' : 'N'), sort='N'; \
237  const BLAS_INT n=BLAS_INT(nn); \
238  std::vector<double> rwork(nn), eigv1(nn), eigv2(nn); \
239  BLAS_INT info(0), lwork(-1), sdim; base_type work1; \
240  lapack_name(&jobvs, &sort, p, &n, &H(0,0), &n, &sdim, &eigv1[0], \
241  &eigv2[0], &Q(0,0), &n, &work1, &lwork, &rwork[0], &info); \
242  const size_type worksize=size_type(gmm::real(work1)); \
243  std::vector<base_type> work(worksize); lwork=BLAS_INT(worksize); \
244  lapack_name(&jobvs, &sort, p, &n, &H(0,0), &n, &sdim, &eigv1[0], \
245  &eigv2[0], &Q(0,0), &n, &work[0], &lwork, &rwork[0],&info);\
246  GMM_ASSERT1(!info, "QR algorithm failed"); \
247  extract_eig(H, eigval_, tol); \
248  }
249 
250 # define gees_interface_cplx(lapack_name, base_type) \
251  template <typename VECT> inline void implicit_qr_algorithm( \
252  const dense_matrix<base_type> &A, VECT &eigval_, \
253  dense_matrix<base_type> &Q, \
254  double tol=gmm::default_tol(base_type()), bool compvect = true) { \
255  GMMLAPACK_TRACE("gees_interface2"); \
256  typedef bool (*L_fp)(...); L_fp p = 0; \
257  const size_type nn(mat_nrows(A)); \
258  if (!nn) return; \
259  dense_matrix<base_type> H(nn,nn); gmm::copy(A, H); \
260  const char jobvs=(compvect ? 'V' : 'N'), sort='N'; \
261  const BLAS_INT n=BLAS_INT(nn); \
262  std::vector<double> rwork(nn), eigvv(nn*2); \
263  BLAS_INT info(0), lwork(-1), sdim; base_type work1; \
264  lapack_name(&jobvs, &sort, p, &n, &H(0,0), &n, &sdim, &eigvv[0], \
265  &Q(0,0), &n, &work1, &lwork, &rwork[0], &rwork[0], &info); \
266  const size_type worksize=size_type(gmm::real(work1)); \
267  std::vector<base_type> work(worksize); lwork=BLAS_INT(worksize); \
268  lapack_name(&jobvs, &sort, p, &n, &H(0,0), &n, &sdim, &eigvv[0], \
269  &Q(0,0), &n, &work[0], &lwork, &rwork[0], &rwork[0],&info);\
270  GMM_ASSERT1(!info, "QR algorithm failed"); \
271  extract_eig(H, eigval_, tol); \
272  }
273 
274  gees_interface(sgees_, BLAS_S)
275  gees_interface(dgees_, BLAS_D)
276  gees_interface_cplx(cgees_, BLAS_C)
277  gees_interface_cplx(zgees_, BLAS_Z)
278 
279 
280 # define geev_interface(f_name, lapack_name, lapack_name_cplx, \
281  base_type, cplx_type, _jobvl, _jobvr) \
282  template <typename VECT> inline void \
283  f_name(const dense_matrix<base_type> &A, VECT &eigval_, \
284  dense_matrix<base_type> &Q) { \
285  GMMLAPACK_TRACE("geev_interface"); \
286  const size_type nn(mat_nrows(A)); \
287  if (!nn) return; \
288  dense_matrix<base_type> H(nn,nn); gmm::copy(A, H); \
289  const char jobvl=_jobvl, jobvr=_jobvl; \
290  const BLAS_INT n=BLAS_INT(nn); BLAS_INT info(0), lwork(-1); \
291  std::vector<base_type> eigvr(nn), eigvi(nn); \
292  base_type work1; \
293  lapack_name(&jobvl, &jobvr, &n, &H(0,0), &n, &eigvr[0], &eigvi[0], \
294  &Q(0,0), &n, &Q(0,0), &n, &work1, &lwork, &info); \
295  const size_type worksize=size_type(gmm::real(work1)); \
296  std::vector<base_type> work(worksize); lwork=BLAS_INT(worksize); \
297  lapack_name(&jobvl, &jobvr, &n, &H(0,0), &n, &eigvr[0], &eigvi[0], \
298  &Q(0,0), &n, &Q(0,0), &n, &work[0], &lwork, &info); \
299  GMM_ASSERT1(!info, "QR algorithm failed"); \
300  gmm::copy(eigvr, gmm::real_part(eigval_)); \
301  gmm::copy(eigvi, gmm::imag_part(eigval_)); \
302  } \
303  template <typename VECT> inline void \
304  f_name(const dense_matrix<cplx_type> &A, VECT &eigval_, \
305  dense_matrix<cplx_type> &Q) { \
306  GMMLAPACK_TRACE("geev_interface"); \
307  const size_type nn(mat_nrows(A)); \
308  if (!nn) return; \
309  dense_matrix<cplx_type> H(nn,nn); gmm::copy(A, H); \
310  const char jobvl=_jobvl, jobvr=_jobvl; \
311  const BLAS_INT n=BLAS_INT(nn); BLAS_INT info(0), lwork(-1); \
312  std::vector<cplx_type> eigv(nn); std::vector<base_type> rwork(2*nn); \
313  cplx_type work1; \
314  lapack_name(&jobvl, &jobvr, &n, &H(0,0), &n, &eigv[0], &Q(0,0), &n, \
315  &Q(0,0), &n, &work1, &lwork, &rwork[0], &info); \
316  const size_type worksize=size_type(gmm::real(work1)); \
317  std::vector<cplx_type> work(worksize); lwork=BLAS_INT(worksize); \
318  lapack_name_cplx(&jobvl, &jobvr, &n, &H(0,0), &n, &eigv[0], &Q(0,0), \
319  &n, &Q(0,0), &n, &work[0], &lwork, &rwork[0], &info);\
320  GMM_ASSERT1(!info, "QR algorithm failed"); \
321  gmm::copy(eigv, eigval_); \
322  }
323 
324  geev_interface(geev_interface_right, sgeev_, cgeev_, BLAS_S, BLAS_C, 'N', 'V')
325  geev_interface(geev_interface_right, dgeev_, zgeev_, BLAS_D, BLAS_Z, 'N', 'V')
326  geev_interface(geev_interface_left, sgeev_, cgeev_, BLAS_S, BLAS_C, 'V', 'N')
327  geev_interface(geev_interface_left, dgeev_, zgeev_, BLAS_D, BLAS_Z, 'V', 'N')
328 
329  /* ********************************************************************** */
330  /* SCHUR algorithm: */
331  /* A = Q*S*(Q^T), with Q orthogonal and S upper quasi-triangula */
332  /* ********************************************************************** */
333 
334 # define geesx_interface(lapack_name, lapack_name_cplx, \
335  base_type, cplx_type) \
336  inline void schur(const dense_matrix<base_type> &A, \
337  dense_matrix<base_type> &S, \
338  dense_matrix<base_type> &Q) { \
339  GMMLAPACK_TRACE("geesx_interface"); \
340  const size_type mm(mat_nrows(A)), nn(mat_ncols(A)), worksize(8*nn); \
341  GMM_ASSERT1(mm==nn, "Schur decomposition requires square matrix"); \
342  resize(S, nn, nn); copy(A, S); resize(Q, nn, nn); \
343  const char jobvs='V', sort='N', sense='N'; const bool select=false; \
344  const BLAS_INT n=BLAS_INT(nn), lwork=BLAS_INT(worksize), liwork(1); \
345  BLAS_INT info(0), sdim(0); base_type rconde(0), rcondv(0); \
346  std::vector<base_type> work(worksize), wr(nn), wi(nn); \
347  std::vector<BLAS_INT> iwork(1), bwork(1); \
348  lapack_name(&jobvs, &sort, &select, &sense, &n, &S(0,0), &n, \
349  &sdim, &wr[0], &wi[0], &Q(0,0), &n, &rconde, &rcondv, \
350  &work[0], &lwork, &iwork[0], &liwork, &bwork[0], &info);\
351  GMM_ASSERT1(!info, "SCHUR algorithm failed"); \
352  } \
353  inline void schur(const dense_matrix<cplx_type> &A, \
354  dense_matrix<cplx_type> &S, \
355  dense_matrix<cplx_type> &Q) { \
356  GMMLAPACK_TRACE("geesx_interface"); \
357  const size_type mm(mat_nrows(A)), nn(mat_ncols(A)), worksize(8*nn); \
358  GMM_ASSERT1(mm==nn, "Schur decomposition requires square matrix"); \
359  resize(S, nn, nn); copy(A, S); resize(Q, nn, nn); \
360  const char jobvs='V', sort='N', sense='N'; const bool select=false; \
361  const BLAS_INT n=BLAS_INT(nn), lwork=BLAS_INT(worksize); \
362  BLAS_INT info(0), sdim(0); base_type rconde(0), rcondv(0); \
363  std::vector<cplx_type> work(worksize), w(nn); \
364  std::vector<base_type> rwork(worksize); \
365  std::vector<BLAS_INT> bwork(1); \
366  lapack_name_cplx(&jobvs, &sort, &select, &sense, &n, &S(0,0), &n, \
367  &sdim, &w[0], &Q(0,0), &n, &rconde, &rcondv, \
368  &work[0], &lwork, &rwork[0], &bwork[0], &info); \
369  GMM_ASSERT1(!info, "SCHUR algorithm failed"); \
370  }
371 
372  geesx_interface(sgeesx_, cgeesx_, BLAS_S, BLAS_C)
373  geesx_interface(dgeesx_, zgeesx_, BLAS_D, BLAS_Z)
374 
375 
376  /* ********************************************************************** */
377  /* Interface to SVD. Does not correspond to a Gmm++ functionality. */
378  /* Author : Sebastian Nowozin <[email protected]> */
379  /* ********************************************************************** */
380 
381 # define gesvd_interface(lapack_name, lapack_name_cplx, \
382  base_type, cplx_type) \
383  inline void svd(const dense_matrix<base_type> &X, \
384  dense_matrix<base_type> &U, \
385  dense_matrix<base_type> &Vtransposed, \
386  std::vector<base_type> &sigma) { \
387  GMMLAPACK_TRACE("gesvd_interface"); \
388  const size_type mm(mat_nrows(X)), nn(mat_ncols(X)), \
389  mn_min(mm < nn ? mm : nn), worksize(15*mn_min); \
390  sigma.resize(mn_min); resize(U, mm, mm); resize(Vtransposed, nn, nn); \
391  const char job='A'; const BLAS_INT m=BLAS_INT(mm), n=BLAS_INT(nn), \
392  lwork=BLAS_INT(worksize); \
393  std::vector<base_type> work(worksize); BLAS_INT info(0); \
394  lapack_name(&job, &job, &m, &n, &X(0,0), &m, &sigma[0], &U(0,0), \
395  &m, &Vtransposed(0,0), &n, &work[0], &lwork, &info); \
396  } \
397  inline void svd(const dense_matrix<cplx_type> &X, \
398  dense_matrix<cplx_type> &U, \
399  dense_matrix<cplx_type> &Vtransposed, \
400  std::vector<base_type> &sigma) { \
401  GMMLAPACK_TRACE("gesvd_interface"); \
402  const size_type mm(mat_nrows(X)), nn(mat_ncols(X)), \
403  mn_min(mm < nn ? mm : nn), worksize(15*mn_min); \
404  sigma.resize(mn_min); resize(U, mm, mm); resize(Vtransposed, nn, nn); \
405  const char job='A'; const BLAS_INT m=BLAS_INT(mm), n=BLAS_INT(nn), \
406  lwork=BLAS_INT(worksize); \
407  std::vector<cplx_type> work(worksize); \
408  std::vector<base_type> rwork(5*mn_min); BLAS_INT info(0); \
409  lapack_name_cplx(&job, &job, &m, &n, &X(0,0), &m, &sigma[0], &U(0,0), \
410  &m, &Vtransposed(0,0), &n, &work[0], &lwork, \
411  &rwork[0], &info); \
412  }
413 
414  gesvd_interface(sgesvd_, cgesvd_, BLAS_S, BLAS_C)
415  gesvd_interface(dgesvd_, zgesvd_, BLAS_D, BLAS_Z)
416 
417 }
418 
419 #else
420 
421 namespace gmm
422 {
423 template <typename MAT>
424 void schur(const MAT &, MAT &, MAT &)
425 {
426  GMM_ASSERT1(false, "Use of function schur(A,S,Q) requires GetFEM "
427  "to be built with Lapack");
428 }
429 
430 template <typename BLAS_TYPE>
431 inline void svd(dense_matrix<BLAS_TYPE> &, dense_matrix<BLAS_TYPE> &,
432  dense_matrix<BLAS_TYPE> &, std::vector<BLAS_TYPE> &)
433 {
434  GMM_ASSERT1(false, "Use of function svd(X,U,Vtransposed,sigma) requires GetFEM "
435  "to be built with Lapack");
436 }
437 
438 }// namespace gmm
439 
440 #endif // GMM_USES_LAPACK
441 
442 #endif // GMM_LAPACK_INTERFACE_H
gmm interface for fortran BLAS.
LU factorizations and determinant computation for dense matrices.
Dense QR factorization.