fml  0.1-0
Fused Matrix Library
scalapack.hh
1 // This file is part of fml which is released under the Boost Software
2 // License, Version 1.0. See accompanying file LICENSE or copy at
3 // https://www.boost.org/LICENSE_1_0.txt
4 
5 #ifndef FML_MPI_LINALG_INTERNALS_SCALAPACK_H
6 #define FML_MPI_LINALG_INTERNALS_SCALAPACK_H
7 #pragma once
8 
9 
10 #include "scalapack_prototypes.h"
11 
12 
13 namespace fml
14 {
15  namespace scalapack
16  {
17  inline void getrf(const int m, const int n, float *a, int *desca, int *ipiv,
18  int *info)
19  {
20  int ij = 1;
21  psgetrf_(&m, &n, a, &ij, &ij, desca, ipiv, info);
22  }
23 
24  inline void getrf(const int m, const int n, double *a, int *desca, int *ipiv,
25  int *info)
26  {
27  int ij = 1;
28  pdgetrf_(&m, &n, a, &ij, &ij, desca, ipiv, info);
29  }
30 
31 
32 
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)
36  {
37  int ij = 1;
38  psgesvd_(&jobu, &jobvt, &m, &n, a, &ij, &ij, desca, s, u, &ij, &ij, descu,
39  vt, &ij, &ij, descvt, work, &lwork, info);
40  }
41 
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)
45  {
46  int ij = 1;
47  pdgesvd_(&jobu, &jobvt, &m, &n, a, &ij, &ij, desca, s, u, &ij, &ij, descu,
48  vt, &ij, &ij, descvt, work, &lwork, info);
49  }
50 
51 
52 
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)
58  {
59  int ij = 1;
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);
62  }
63 
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)
69  {
70  int ij = 1;
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);
73  }
74 
75 
76 
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)
79  {
80  int ij = 1;
81  psgetri_(&n, a, &ij, &ij, desca, ipiv, work, &lwork, iwork, &liwork, info);
82  }
83 
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)
86  {
87  int ij = 1;
88  pdgetri_(&n, a, &ij, &ij, desca, ipiv, work, &lwork, iwork, &liwork, info);
89  }
90 
91 
92 
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)
95  {
96  int ij = 1;
97  psgesv_(&n, &nrhs, a, &ij, &ij, desca, ipvt, b, &ij, &ij, descb, info);
98  }
99 
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)
102  {
103  int ij = 1;
104  pdgesv_(&n, &nrhs, a, &ij, &ij, desca, ipvt, b, &ij, &ij, descb, info);
105  }
106 
107 
108 
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)
111  {
112  int ij = 1;
113  pslacpy_(&uplo, &m, &n, a, &ij, &ij, desca, b, &ij, &ij, descb);
114  }
115 
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)
118  {
119  int ij = 1;
120  pdlacpy_(&uplo, &m, &n, a, &ij, &ij, desca, b, &ij, &ij, descb);
121  }
122 
123 
124 
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)
127  {
128  int ij = 1;
129  psgeqpf_(&m, &n, a, &ij, &ij, desca, ipiv, tau, work, &lwork, info);
130  }
131 
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)
134  {
135  int ij = 1;
136  pdgeqpf_(&m, &n, a, &ij, &ij, desca, ipiv, tau, work, &lwork, info);
137  }
138 
139 
140 
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)
143  {
144  int ij = 1;
145  psgeqrf_(&m, &n, a, &ij, &ij, desca, tau, work, &lwork, info);
146  }
147 
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)
150  {
151  int ij = 1;
152  pdgeqrf_(&m, &n, a, &ij, &ij, desca, tau, work, &lwork, info);
153  }
154 
155 
156 
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)
160  {
161  int ij = 1;
162  psormqr_(&side, &trans, &m, &n, &k, a, &ij, &ij, desca, tau, c, &ij, &ij,
163  descc, work, &lwork, info);
164  }
165 
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)
169  {
170  int ij = 1;
171  pdormqr_(&side, &trans, &m, &n, &k, a, &ij, &ij, desca, tau, c, &ij, &ij,
172  descc, work, &lwork, info);
173  }
174 
175 
176 
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)
179  {
180  int ij = 1;
181  psorgqr_(&m, &n, &k, a, &ij, &ij, desca, tau, work, &lwork, info);
182  }
183 
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)
186  {
187  int ij = 1;
188  pdorgqr_(&m, &n, &k, a, &ij, &ij, desca, tau, work, &lwork, info);
189  }
190 
191 
192 
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)
195  {
196  int ij = 1;
197  psorglq_(&m, &n, &k, a, &ij, &ij, desca, tau, work, &lwork, info);
198  }
199 
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)
202  {
203  int ij = 1;
204  pdorglq_(&m, &n, &k, a, &ij, &ij, desca, tau, work, &lwork, info);
205  }
206 
207 
208 
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)
211  {
212  int ij = 1;
213  psgelqf_(&m, &n, a, &ij, &ij, desca, tau, work, &lwork, info);
214  }
215 
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)
218  {
219  int ij = 1;
220  pdgelqf_(&m, &n, a, &ij, &ij, desca, tau, work, &lwork, info);
221  }
222 
223 
224 
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)
228  {
229  int ij = 1;
230  psormlq_(&side, &trans, &m, &n, &k, a, &ij, &ij, desca, tau, c, &ij, &ij,
231  descc, work, &lwork, info);
232  }
233 
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)
237  {
238  int ij = 1;
239  pdormlq_(&side, &trans, &m, &n, &k, a, &ij, &ij, desca, tau, c, &ij, &ij,
240  descc, work, &lwork, info);
241  }
242 
243 
244 
245  inline void potrf(const char uplo, const int n, float *a, const int *desca,
246  int *info)
247  {
248  int ij = 1;
249  pspotrf_(&uplo, &n, a, &ij, &ij, desca, info);
250  }
251 
252  inline void potrf(const char uplo, const int n, double *a, const int *desca,
253  int *info)
254  {
255  int ij = 1;
256  pdpotrf_(&uplo, &n, a, &ij, &ij, desca, info);
257  }
258 
259 
260 
261  inline void lassq(const int n, const float *x, const int *descx,
262  const int incx, float *scale, float *sumsq)
263  {
264  int ij = 1;
265  pslassq_(&n, x, &ij, &ij, descx, &incx, scale, sumsq);
266  }
267 
268  inline void lassq(const int n, const double *x, const int *descx,
269  const int incx, double *scale, double *sumsq)
270  {
271  int ij = 1;
272  pdlassq_(&n, x, &ij, &ij, descx, &incx, scale, sumsq);
273  }
274 
275 
276 
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,
280  int *info)
281  {
282  int ij = 1;
283  psgecon_(&norm, &n, A, &ij, &ij, desca, &anorm, rcond, work, &lwork,
284  iwork, &liwork, info);
285  }
286 
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,
290  int *info)
291  {
292  int ij = 1;
293  pdgecon_(&norm, &n, A, &ij, &ij, desca, &anorm, rcond, work, &lwork,
294  iwork, &liwork, info);
295  }
296 
297 
298 
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)
303  {
304  int ij = 1;
305  pstrcon_(&norm, &uplo, &diag, &n, A, &ij, &ij, desca, rcond, work, &lwork,
306  iwork, &liwork, info);
307  }
308 
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)
313  {
314  int ij = 1;
315  pdtrcon_(&norm, &uplo, &diag, &n, A, &ij, &ij, desca, rcond, work, &lwork,
316  iwork, &liwork, info);
317  }
318 
319 
320 
321  inline void trtri(const char uplo, const char diag, const int n, float *A,
322  const int *desca, int *info)
323  {
324  int ij = 1;
325  pstrtri_(&uplo, &diag, &n, A, &ij, &ij, desca, info);
326  }
327 
328  inline void trtri(const char uplo, const char diag, const int n, double *A,
329  const int *desca, int *info)
330  {
331  int ij = 1;
332  pdtrtri_(&uplo, &diag, &n, A, &ij, &ij, desca, info);
333  }
334  }
335 }
336 
337 
338 #endif
fml
Core namespace.
Definition: dimops.hh:10