5 #ifndef FML_MPI_LINALG_INTERNALS_SCALAPACK_H
6 #define FML_MPI_LINALG_INTERNALS_SCALAPACK_H
10 #include "scalapack_prototypes.h"
17 inline void getrf(
const int m,
const int n,
float *a,
int *desca,
int *ipiv,
21 psgetrf_(&m, &n, a, &ij, &ij, desca, ipiv, info);
24 inline void getrf(
const int m,
const int n,
double *a,
int *desca,
int *ipiv,
28 pdgetrf_(&m, &n, a, &ij, &ij, desca, ipiv, info);
33 inline void gesvd(
const char jobu,
const char jobvt,
const int m,
const int n,
34 float *a,
int *desca,
float *s,
float *u,
int *descu,
float *vt,
35 int *descvt,
float *work,
int lwork,
int *info)
38 psgesvd_(&jobu, &jobvt, &m, &n, a, &ij, &ij, desca, s, u, &ij, &ij, descu,
39 vt, &ij, &ij, descvt, work, &lwork, info);
42 inline void gesvd(
const char jobu,
const char jobvt,
const int m,
const int n,
43 double *a,
int *desca,
double *s,
double *u,
int *descu,
double *vt,
44 int *descvt,
double *work,
int lwork,
int *info)
47 pdgesvd_(&jobu, &jobvt, &m, &n, a, &ij, &ij, desca, s, u, &ij, &ij, descu,
48 vt, &ij, &ij, descvt, work, &lwork, info);
53 inline void syevr(
const char jobz,
const char range,
const char uplo,
54 const int n,
float *a,
const int *desca,
const float vl,
const float vu,
55 const int il,
const int iu,
int *m,
int *nz,
float *w,
float *z,
56 const int *descz,
float *work,
const int lwork,
int *iwork,
57 const int liwork,
int *info)
60 pssyevr_(&jobz, &range, &uplo, &n, a, &ij, &ij, desca, &vl, &vu, &il, &iu,
61 m, nz, w, z, &ij, &ij, descz, work, &lwork, iwork, &liwork, info);
64 inline void syevr(
const char jobz,
const char range,
const char uplo,
65 const int n,
double *a,
const int *desca,
const double vl,
const double vu,
66 const int il,
const int iu,
int *m,
int *nz,
double *w,
double *z,
67 const int *descz,
double *work,
const int lwork,
int *iwork,
68 const int liwork,
int *info)
71 pdsyevr_(&jobz, &range, &uplo, &n, a, &ij, &ij, desca, &vl, &vu, &il, &iu,
72 m, nz, w, z, &ij, &ij, descz, work, &lwork, iwork, &liwork, info);
77 inline void getri(
const int n,
float *a,
const int *desca,
int *ipiv,
78 float *work,
const int lwork,
int *iwork,
const int liwork,
int *info)
81 psgetri_(&n, a, &ij, &ij, desca, ipiv, work, &lwork, iwork, &liwork, info);
84 inline void getri(
const int n,
double *a,
const int *desca,
int *ipiv,
85 double *work,
const int lwork,
int *iwork,
const int liwork,
int *info)
88 pdgetri_(&n, a, &ij, &ij, desca, ipiv, work, &lwork, iwork, &liwork, info);
93 inline void gesv(
const int n,
const int nrhs,
float *a,
const int *desca,
94 int *ipvt,
float *b,
const int *descb,
int *info)
97 psgesv_(&n, &nrhs, a, &ij, &ij, desca, ipvt, b, &ij, &ij, descb, info);
100 inline void gesv(
const int n,
const int nrhs,
double *a,
const int *desca,
101 int *ipvt,
double *b,
const int *descb,
int *info)
104 pdgesv_(&n, &nrhs, a, &ij, &ij, desca, ipvt, b, &ij, &ij, descb, info);
109 inline void lacpy(
const char uplo,
const int m,
const int n,
110 const float *a,
const int *desca,
float *b,
const int *descb)
113 pslacpy_(&uplo, &m, &n, a, &ij, &ij, desca, b, &ij, &ij, descb);
116 inline void lacpy(
const char uplo,
const int m,
const int n,
117 const double *a,
const int *desca,
double *b,
const int *descb)
120 pdlacpy_(&uplo, &m, &n, a, &ij, &ij, desca, b, &ij, &ij, descb);
125 inline void geqpf(
const int m,
const int n,
float *a,
const int *desca,
126 int *ipiv,
float *tau,
float *work,
const int lwork,
int *info)
129 psgeqpf_(&m, &n, a, &ij, &ij, desca, ipiv, tau, work, &lwork, info);
132 inline void geqpf(
const int m,
const int n,
double *a,
const int *desca,
133 int *ipiv,
double *tau,
double *work,
const int lwork,
int *info)
136 pdgeqpf_(&m, &n, a, &ij, &ij, desca, ipiv, tau, work, &lwork, info);
141 inline void geqrf(
const int m,
const int n,
float *a,
const int *desca,
142 float *tau,
float *work,
const int lwork,
int *info)
145 psgeqrf_(&m, &n, a, &ij, &ij, desca, tau, work, &lwork, info);
148 inline void geqrf(
const int m,
const int n,
double *a,
const int *desca,
149 double *tau,
double *work,
const int lwork,
int *info)
152 pdgeqrf_(&m, &n, a, &ij, &ij, desca, tau, work, &lwork, info);
157 inline void ormqr(
const char side,
const char trans,
const int m,
158 const int n,
const int k,
const float *a,
const int *desca,
float *tau,
159 float *c,
const int *descc,
float *work,
const int lwork,
int *info)
162 psormqr_(&side, &trans, &m, &n, &k, a, &ij, &ij, desca, tau, c, &ij, &ij,
163 descc, work, &lwork, info);
166 inline void ormqr(
const char side,
const char trans,
const int m,
167 const int n,
const int k,
const double *a,
const int *desca,
double *tau,
168 double *c,
const int *descc,
double *work,
const int lwork,
int *info)
171 pdormqr_(&side, &trans, &m, &n, &k, a, &ij, &ij, desca, tau, c, &ij, &ij,
172 descc, work, &lwork, info);
177 inline void orgqr(
const int m,
const int n,
const int k,
float *a,
178 const int *desca,
float *tau,
float *work,
const int lwork,
int *info)
181 psorgqr_(&m, &n, &k, a, &ij, &ij, desca, tau, work, &lwork, info);
184 inline void orgqr(
const int m,
const int n,
const int k,
double *a,
185 const int *desca,
double *tau,
double *work,
const int lwork,
int *info)
188 pdorgqr_(&m, &n, &k, a, &ij, &ij, desca, tau, work, &lwork, info);
193 inline void orglq(
const int m,
const int n,
const int k,
float *a,
194 const int *desca,
float *tau,
float *work,
const int lwork,
int *info)
197 psorglq_(&m, &n, &k, a, &ij, &ij, desca, tau, work, &lwork, info);
200 inline void orglq(
const int m,
const int n,
const int k,
double *a,
201 const int *desca,
double *tau,
double *work,
const int lwork,
int *info)
204 pdorglq_(&m, &n, &k, a, &ij, &ij, desca, tau, work, &lwork, info);
209 inline void gelqf(
const int m,
const int n,
float *a,
const int *desca,
210 float *tau,
float *work,
const int lwork,
int *info)
213 psgelqf_(&m, &n, a, &ij, &ij, desca, tau, work, &lwork, info);
216 inline void gelqf(
const int m,
const int n,
double *a,
const int *desca,
217 double *tau,
double *work,
const int lwork,
int *info)
220 pdgelqf_(&m, &n, a, &ij, &ij, desca, tau, work, &lwork, info);
225 inline void ormlq(
const char side,
const char trans,
const int m,
226 const int n,
const int k,
const float *a,
const int *desca,
float *tau,
227 float *c,
const int *descc,
float *work,
const int lwork,
int *info)
230 psormlq_(&side, &trans, &m, &n, &k, a, &ij, &ij, desca, tau, c, &ij, &ij,
231 descc, work, &lwork, info);
234 inline void ormlq(
const char side,
const char trans,
const int m,
235 const int n,
const int k,
const double *a,
const int *desca,
double *tau,
236 double *c,
const int *descc,
double *work,
const int lwork,
int *info)
239 pdormlq_(&side, &trans, &m, &n, &k, a, &ij, &ij, desca, tau, c, &ij, &ij,
240 descc, work, &lwork, info);
245 inline void potrf(
const char uplo,
const int n,
float *a,
const int *desca,
249 pspotrf_(&uplo, &n, a, &ij, &ij, desca, info);
252 inline void potrf(
const char uplo,
const int n,
double *a,
const int *desca,
256 pdpotrf_(&uplo, &n, a, &ij, &ij, desca, info);
261 inline void lassq(
const int n,
const float *x,
const int *descx,
262 const int incx,
float *scale,
float *sumsq)
265 pslassq_(&n, x, &ij, &ij, descx, &incx, scale, sumsq);
268 inline void lassq(
const int n,
const double *x,
const int *descx,
269 const int incx,
double *scale,
double *sumsq)
272 pdlassq_(&n, x, &ij, &ij, descx, &incx, scale, sumsq);
277 inline void gecon(
const char norm,
const int n,
278 const float *
const restrict A,
const int *desca,
const float anorm,
279 float *rcond,
float *work,
const int lwork,
int *iwork,
const int liwork,
283 psgecon_(&norm, &n, A, &ij, &ij, desca, &anorm, rcond, work, &lwork,
284 iwork, &liwork, info);
287 inline void gecon(
const char norm,
const int n,
288 const double *
const restrict A,
const int *desca,
const double anorm,
289 double *rcond,
double *work,
const int lwork,
int *iwork,
const int liwork,
293 pdgecon_(&norm, &n, A, &ij, &ij, desca, &anorm, rcond, work, &lwork,
294 iwork, &liwork, info);
299 inline void trcon(
const char norm,
const char uplo,
const char diag,
300 const int n,
const float *
const restrict A,
const int *desca,
301 float *rcond,
float *work,
const int lwork,
int *iwork,
302 const int liwork,
int *info)
305 pstrcon_(&norm, &uplo, &diag, &n, A, &ij, &ij, desca, rcond, work, &lwork,
306 iwork, &liwork, info);
309 inline void trcon(
const char norm,
const char uplo,
const char diag,
310 const int n,
const double *
const restrict A,
const int *desca,
311 double *rcond,
double *work,
const int lwork,
int *iwork,
312 const int liwork,
int *info)
315 pdtrcon_(&norm, &uplo, &diag, &n, A, &ij, &ij, desca, rcond, work, &lwork,
316 iwork, &liwork, info);
321 inline void trtri(
const char uplo,
const char diag,
const int n,
float *A,
322 const int *desca,
int *info)
325 pstrtri_(&uplo, &diag, &n, A, &ij, &ij, desca, info);
328 inline void trtri(
const char uplo,
const char diag,
const int n,
double *A,
329 const int *desca,
int *info)
332 pdtrtri_(&uplo, &diag, &n, A, &ij, &ij, desca, info);