5 #ifndef FML_CPU_LINALG_INTERNALS_LAPACK_H
6 #define FML_CPU_LINALG_INTERNALS_LAPACK_H
10 #include "lapack_prototypes.h"
17 inline void getrf(
const int m,
const int n,
float *a,
const int lda,
20 sgetrf_(&m, &n, a, &lda, ipiv, info);
23 inline void getrf(
const int m,
const int n,
double *a,
const int lda,
26 dgetrf_(&m, &n, a, &lda, ipiv, info);
31 inline void gesdd(
const char jobz,
const int m,
const int n,
float *a,
32 const int lda,
float *s,
float *u,
const int ldu,
float *vt,
33 const int ldvt,
float *work,
const int lwork,
int *iwork,
int *info)
35 sgesdd_(&jobz, &m, &n, a, &lda, s, u, &ldu, vt, &ldvt, work, &lwork, iwork,
39 inline void gesdd(
const char jobz,
const int m,
const int n,
double *a,
40 const int lda,
double *s,
double *u,
const int ldu,
double *vt,
41 const int ldvt,
double *work,
const int lwork,
int *iwork,
int *info)
43 dgesdd_(&jobz, &m, &n, a, &lda, s, u, &ldu, vt, &ldvt, work, &lwork, iwork,
49 inline void syevr(
const char jobz,
const char range,
const char uplo,
50 const int n,
float *a,
const int lda,
const float vl,
const float vu,
51 const int il,
const int iu,
const float abstol,
int *m,
float *w,
52 float *z,
const int ldz,
int *isuppz,
float *work,
const int lwork,
53 int *iwork,
const int liwork,
int *info)
55 ssyevr_(&jobz, &range, &uplo, &n, a, &lda, &vl, &vu, &il, &iu, &abstol, m,
56 w, z, &ldz, isuppz, work, &lwork, iwork, &liwork, info);
59 inline void syevr(
const char jobz,
const char range,
const char uplo,
60 const int n,
double *a,
const int lda,
const double vl,
const double vu,
61 const int il,
const int iu,
const double abstol,
int *m,
double *w,
62 double *z,
const int ldz,
int *isuppz,
double *work,
const int lwork,
63 int *iwork,
const int liwork,
int *info)
65 dsyevr_(&jobz, &range, &uplo, &n, a, &lda, &vl, &vu, &il, &iu, &abstol, m,
66 w, z, &ldz, isuppz, work, &lwork, iwork, &liwork, info);
71 inline void geev(
const char jobvl,
const char jobvr,
const int n,
float *a,
72 const int lda,
float *wr,
float *wi,
float *vl,
const int ldvl,
float *vr,
73 const int ldvr,
float *work,
const int lwork,
int *info)
75 sgeev_(&jobvl, &jobvr, &n, a, &lda, wr, wi, vl, &ldvl, vr, &ldvr, work,
79 inline void geev(
const char jobvl,
const char jobvr,
const int n,
double *a,
80 const int lda,
double *wr,
double *wi,
double *vl,
const int ldvl,
double *vr,
81 const int ldvr,
double *work,
const int lwork,
int *info)
83 dgeev_(&jobvl, &jobvr, &n, a, &lda, wr, wi, vl, &ldvl, vr, &ldvr, work,
89 inline void getri(
const int n,
float *a,
const int lda,
int *ipiv,
90 float *work,
const int lwork,
int *info)
92 sgetri_(&n, a, &lda, ipiv, work, &lwork, info);
95 inline void getri(
const int n,
double *a,
const int lda,
int *ipiv,
96 double *work,
const int lwork,
int *info)
98 dgetri_(&n, a, &lda, ipiv, work, &lwork, info);
103 inline void gesv(
const int n,
const int nrhs,
float *a,
const int lda,
104 int *ipiv,
float *b,
const int ldb,
int *info)
106 sgesv_(&n, &nrhs, a, &lda, ipiv, b, &ldb, info);
109 inline void gesv(
const int n,
const int nrhs,
double *a,
const int lda,
110 int *ipiv,
double *b,
const int ldb,
int *info)
112 dgesv_(&n, &nrhs, a, &lda, ipiv, b, &ldb, info);
117 inline void lacpy(
const char uplo,
const int m,
const int n,
const float *a,
118 const int lda,
float *b,
const int ldb)
120 slacpy_(&uplo, &m, &n, a, &lda, b, &ldb);
123 inline void lacpy(
const char uplo,
const int m,
const int n,
const double *a,
124 const int lda,
double *b,
const int ldb)
126 dlacpy_(&uplo, &m, &n, a, &lda, b, &ldb);
131 inline void geqp3(
const int m,
const int n,
float *x,
const int lda,
132 int *pivot,
float *qraux,
float *work,
const int lwork,
int *info)
134 sgeqp3_(&m, &n, x, &lda, pivot, qraux, work, &lwork, info);
137 inline void geqp3(
const int m,
const int n,
double *x,
const int lda,
138 int *pivot,
double *qraux,
double *work,
const int lwork,
int *info)
140 dgeqp3_(&m, &n, x, &lda, pivot, qraux, work, &lwork, info);
145 inline void geqrf(
const int m,
const int n,
float *x,
const int lda,
146 float *tau,
float *work,
const int lwork,
int *info)
148 sgeqrf_(&m, &n, x, &lda, tau, work, &lwork, info);
151 inline void geqrf(
const int m,
const int n,
double *x,
const int lda,
152 double *tau,
double *work,
const int lwork,
int *info)
154 dgeqrf_(&m, &n, x, &lda, tau, work, &lwork, info);
159 inline void ormqr(
const char side,
const char trans,
const int m,
160 const int n,
const int k,
const float *x,
const int lda,
161 const float *tau,
float *c,
const int ldc,
float *work,
162 const int lwork,
int *info)
164 sormqr_(&side, &trans, &m, &n, &k, x, &lda, tau, c, &ldc, work, &lwork, info);
167 inline void ormqr(
const char side,
const char trans,
const int m,
168 const int n,
const int k,
const double *x,
const int lda,
169 const double *tau,
double *c,
const int ldc,
double *work,
170 const int lwork,
int *info)
172 dormqr_(&side, &trans, &m, &n, &k, x, &lda, tau, c, &ldc, work, &lwork, info);
177 inline void orgqr(
const int m,
const int n,
const int k,
float *A,
178 const int lda,
const float *tau,
float *work,
const int ldwork,
int *info)
180 sorgqr_(&m, &n, &k, A, &lda, tau, work, &ldwork, info);
183 inline void orgqr(
const int m,
const int n,
const int k,
double *A,
184 const int lda,
const double *tau,
double *work,
const int ldwork,
int *info)
186 dorgqr_(&m, &n, &k, A, &lda, tau, work, &ldwork, info);
191 inline void orglq(
const int m,
const int n,
const int k,
float *A,
192 const int lda,
const float *tau,
float *work,
const int ldwork,
int *info)
194 sorglq_(&m, &n, &k, A, &lda, tau, work, &ldwork, info);
197 inline void orglq(
const int m,
const int n,
const int k,
double *A,
198 const int lda,
const double *tau,
double *work,
const int ldwork,
int *info)
200 dorglq_(&m, &n, &k, A, &lda, tau, work, &ldwork, info);
205 inline void gelqf(
const int m,
const int n,
float *x,
const int lda,
206 float *tau,
float *work,
const int lwork,
int *info)
208 sgelqf_(&m, &n, x, &lda, tau, work, &lwork, info);
211 inline void gelqf(
const int m,
const int n,
double *x,
const int lda,
212 double *tau,
double *work,
const int lwork,
int *info)
214 dgelqf_(&m, &n, x, &lda, tau, work, &lwork, info);
219 inline void ormlq(
const char side,
const char trans,
const int m,
220 const int n,
const int k,
const float *x,
const int lda,
221 const float *tau,
float *c,
const int ldc,
float *work,
222 const int lwork,
int *info)
224 sormlq_(&side, &trans, &m, &n, &k, x, &lda, tau, c, &ldc, work, &lwork, info);
227 inline void ormlq(
const char side,
const char trans,
const int m,
228 const int n,
const int k,
const double *x,
const int lda,
229 const double *tau,
double *c,
const int ldc,
double *work,
230 const int lwork,
int *info)
232 dormlq_(&side, &trans, &m, &n, &k, x, &lda, tau, c, &ldc, work, &lwork, info);
237 inline void potrf(
const char uplo,
const int n,
float *A,
const int lda,
240 spotrf_(&uplo, &n, A, &lda, info);
243 inline void potrf(
const char uplo,
const int n,
double *A,
const int lda,
246 dpotrf_(&uplo, &n, A, &lda, info);
251 inline void lassq(
const int n,
const float *x,
const int incx,
252 float *scale,
float *sumsq)
254 slassq_(&n, x, &incx, scale, sumsq);
257 inline void lassq(
const int n,
const double *x,
const int incx,
258 double *scale,
double *sumsq)
260 dlassq_(&n, x, &incx, scale, sumsq);
265 inline void gecon(
const char norm,
const int n,
const float *A,
266 const int lda,
const float anorm,
float *rcond,
float *work,
267 float *work2,
int *info)
269 sgecon_(&norm, &n, A, &lda, &anorm, rcond, work, work2, info);
272 inline void gecon(
const char norm,
const int n,
const double *A,
273 const int lda,
const double anorm,
double *rcond,
double *work,
274 double *work2,
int *info)
276 dgecon_(&norm, &n, A, &lda, &anorm, rcond, work, work2, info);
281 inline void trcon(
const char norm,
const char uplo,
const char diag,
282 const int n,
const float *A,
const int lda,
float *rcond,
float *work,
283 float *work2,
int *info)
285 strcon_(&norm, &uplo, &diag, &n, A, &lda, rcond, work, work2, info);
288 inline void trcon(
const char norm,
const char uplo,
const char diag,
289 const int n,
const double *A,
const int lda,
double *rcond,
double *work,
290 double *work2,
int *info)
292 dtrcon_(&norm, &uplo, &diag, &n, A, &lda, rcond, work, work2, info);
297 inline void trtri(
const char uplo,
const char diag,
const int n,
298 float *A,
const int lda,
int *info)
300 strtri_(&uplo, &diag, &n, A, &lda, info);
303 inline void trtri(
const char uplo,
const char diag,
const int n,
304 double *A,
const int lda,
int *info)
306 dtrtri_(&uplo, &diag, &n, A, &lda, info);