Eigenvalue SoLvers for Petaflop-Applications (ELPA)  2018.05.001
elpa_generic.h
Go to the documentation of this file.
1 #pragma once
2 
12 #define elpa_set(e, name, value, error) _Generic((value), \
13  int: \
14  elpa_set_integer, \
15  \
16  double: \
17  elpa_set_double \
18  )(e, name, value, error)
19 
20 
30 #define elpa_get(e, name, value, error) _Generic((value), \
31  int*: \
32  elpa_get_integer, \
33  \
34  double*: \
35  elpa_get_double \
36  )(e, name, value, error)
37 
38 
49 #define elpa_eigenvectors(handle, a, ev, q, error) _Generic((a), \
50  double*: \
51  elpa_eigenvectors_d, \
52  \
53  float*: \
54  elpa_eigenvectors_f, \
55  \
56  double complex*: \
57  elpa_eigenvectors_dc, \
58  \
59  float complex*: \
60  elpa_eigenvectors_fc \
61  )(handle, a, ev, q, error)
62 
63 
76 #define elpa_generalized_eigenvectors(handle, a, b, ev, q, is_already_decomposed, error) _Generic((a), \
77  double*: \
78  elpa_generalized_eigenvectors_d, \
79  \
80  float*: \
81  elpa_generalized_eigenvectors_f, \
82  \
83  double complex*: \
84  elpa_generalized_eigenvectors_dc, \
85  \
86  float complex*: \
87  elpa_generalized_eigenvectors_fc \
88  )(handle, a, b, ev, q, is_already_decomposed, error)
89 
90 
100 #define elpa_eigenvalues(handle, a, ev, error) _Generic((a), \
101  double*: \
102  elpa_eigenvalues_d, \
103  \
104  float*: \
105  elpa_eigenvalues_f, \
106  \
107  double complex*: \
108  elpa_eigenvalues_dc, \
109  \
110  float complex*: \
111  elpa_eigenvalues_fc \
112  )(handle, a, ev, error)
113 
114 /* \brief generic C method for elpa_cholesky
115  *
116  * \details
117  * \param handle handle of the ELPA object, which defines the problem
118  * \param a float/double float complex/double complex pointer to matrix a, for which
119  * the cholesky factorizaion will be computed
120  * \param error on return the error code, which can be queried with elpa_strerr()
121  * \result void
122  */
123 #define elpa_cholesky(handle, a, error) _Generic((a), \
124  double*: \
125  elpa_cholesky_d, \
126  \
127  float*: \
128  elpa_cholesky_f, \
129  \
130  double complex*: \
131  elpa_cholesky_dc, \
132  \
133  float complex*: \
134  elpa_cholesky_fc \
135  )(handle, a, error)
136 
137 
155 #define elpa_hermitian_multiply(handle, uplo_a, uplo_c, ncb, a, b, nrows_b, ncols_b, c, nrows_c, ncols_c, error) _Generic((a), \
156  double*: \
157  elpa_hermitian_multiply_d, \
158  \
159  float*: \
160  elpa_hermitian_multiply_f, \
161  \
162  double complex*: \
163  elpa_hermitian_multiply_dc, \
164  \
165  float complex*: \
166  elpa_hermitian_multiply_fc \
167  )(handle, a, error)
168 
169 
179 #define elpa_invert_triangular(handle, a, error) _Generic((a), \
180  double*: \
181  elpa_invert_trm_d, \
182  \
183  float*: \
184  elpa_invert_trm_f, \
185  \
186  double complex*: \
187  elpa_invert_trm_dc, \
188  \
189  float complex*: \
190  elpa_invert_trm_fc \
191  )(handle, a, error)