KIDS  ver-0.0.1
KIDS : Kernel Integrated Dynamics Simulator
Loading...
Searching...
No Matches
array_macro.h
Go to the documentation of this file.
1
12#ifndef Array_Macro_H
13#define Array_Macro_H
14#include <type_traits>
15
16#define ARRAY_0_OP1(_A, _op1, _C, _n) \
17 ({ \
18 for (int _i = 0; _i < (_n); ++_i) (_A)[_i] _op1(_C); \
19 })
20
21#define ARRAY_1_OP1(_A, _op1, _C, _n) \
22 ({ \
23 for (int _i = 0; _i < (_n); ++_i) (_A)[_i] _op1(_C)[_i]; \
24 })
25
26#define ARRAY_0_OP2(_A, _B, _op2, _C, _n) \
27 ({ \
28 for (int _i = 0; _i < (_n); ++_i) (_A)[_i] = (_B)[_i] _op2(_C); \
29 })
30
31#define ARRAY_1_OP2(_A, _B, _op2, _C, _n) \
32 ({ \
33 for (int _i = 0; _i < (_n); ++_i) (_A)[_i] = (_B)[_i] _op2(_C)[_i]; \
34 })
35
36#define ARRAY_CLEAR(_A, _n) (memset((_A), 0, (_n) * sizeof(*(_A))))
37#define ARRAY_0_COPY(_A, _C, _n) (ARRAY_0_OP1(_A, =, _C, _n))
38#define ARRAY_1_COPY(_A, _C, _n) (ARRAY_1_OP1(_A, =, _C, _n))
39#define ARRAY_1_ADD(_A, _B, _C, _n) (ARRAY_1_OP2(_A, _B, +, _C, _n))
40#define ARRAY_1_MIN(_A, _B, _C, _n) (ARRAY_1_OP2(_A, _B, -, _C, _n))
41#define ARRAY_1_MUL(_A, _B, _C, _n) (ARRAY_1_OP2(_A, _B, *, _C, _n))
42#define ARRAY_1_DIV(_A, _B, _C, _n) (ARRAY_1_OP2(_A, _B, /, _C, _n))
43
44// @noused
45#define ARRAY_SUM(_A, _n) \
46 ({ \
47 auto _res = (std::remove_reference<decltype(*_A)>::type)(0); \
48 for (int _i = 0; _i < (_n); ++_i) _res += (_A)[_i]; \
49 _res; \
50 })
51
52#define ARRAY_0_FUN(_A, _fun, _C, _n) \
53 ({ \
54 for (int _i = 0; _i < (_n); ++_i) (_A)[_i] = _fun((_C)); \
55 })
56
57#define ARRAY_1_FUN(_A, _fun, _C, _n) \
58 ({ \
59 for (int _i = 0; _i < (_n); ++_i) (_A)[_i] = _fun((_C)[_i]); \
60 })
61
62
63// specicals for complex array // @noused
64#define ARRAY_SIN(_A, _C, _n) (ARRAY_1_FUN(_A, std::sin, _C, _n))
65#define ARRAY_COS(_A, _C, _n) (ARRAY_1_FUN(_A, std::cos, _C, _n))
66#define ARRAY_EXP(_A, _C, _n) (ARRAY_1_FUN(_A, std::exp, _C, _n))
67#define ARRAY_LOG(_A, _C, _n) (ARRAY_1_FUN(_A, std::log, _C, _n))
68#define ARRAY_CONJ(_A, _C, _n) (ARRAY_1_FUN(_A, CONJ_OF, _C, _n))
69#define ARRAY_REAL(_A, _C, _n) (ARRAY_1_FUN(_A, REAL_OF, _C, _n))
70#define ARRAY_IMAG(_A, _C, _n) (ARRAY_1_FUN(_A, IMAG_OF, _C, _n))
71#define ARRAY_NORM(_A, _C, _n) (ARRAY_1_FUN(_A, NORM_OF, _C, _n))
72#define ARRAY_PHAS(_A, _C, _n) (ARRAY_1_FUN(_A, PHAS_OF, _C, _n))
73#define ARRAY_ABS(_A, _C, _n) (ARRAY_1_FUN(_A, ABS_OF, _C, _n))
74
75#define ARRAY_EYE(_A, _n) \
76 ({ \
77 for (int _i = 0, _idx = 0; _i < (_n); ++_i) { \
78 for (int _j = 0; _j < (_n); ++_j, ++_idx) { (_A)[_idx] = ((_i == _j) ? 1.0f : 0.0f); } \
79 } \
80 })
81
82// @noused
83#define ARRAY_ONES(_A, _n1, _n2) \
84 ({ \
85 for (int _i = 0, _idx = 0; _i < (_n1); ++_i) { \
86 for (int _j = 0; _j < (_n2); ++_j, ++_idx) { (_A)[_idx] = 1.0f; } \
87 } \
88 })
89
90// @noused
91#define ARRAY_ZEROS(_A, _n1, _n2) \
92 ({ \
93 for (int _i = 0, _idx = 0; _i < (_n1); ++_i) { \
94 for (int _j = 0; _j < (_n2); ++_j, ++_idx) { (_A)[_idx] = 0.0f; } \
95 } \
96 })
97
98// @noused
99// array A is set to diagonal of matrix B
100#define ARRAY_VDIAGM(_A, _B, _n) \
101 ({ \
102 int _idx = 0, _add = (_n) + 1; \
103 for (int _i = 0; _i < (_n); ++_i) { \
104 (_A)[_i] = (_B)[_idx]; \
105 _idx += _add; \
106 } \
107 })
108
109// @noused
110// diagonal of array A is set to array B
111#define ARRAY_MDIAGV(_A, _B, _n) \
112 ({ \
113 int _idx = 0, _add = (_n) + 1; \
114 for (int _i = 0; _i < (_n); ++_i) { \
115 (_A)[_idx] = (_B)[i]; \
116 _idx += _add; \
117 } \
118 })
119
120// return Tr[B*C]: B[n1,n2], C[n2, n1]
121#define ARRAY_TRACE2(_B, _C, _n1, _n2) \
122 ({ \
123 auto _res = (std::remove_reference<decltype(*_B)>::type)(0); \
124 int _idxB = 0; \
125 for (int _i = 0; _i < (_n1); ++_i) { \
126 int _idxC = _i; \
127 for (int _j = 0; _j < (_n2); ++_j) { \
128 _res += (_B)[_idxB++] * (_C)[_idxC]; \
129 _idxC += (_n1); \
130 } \
131 } \
132 _res; \
133 })
134
135// A = B*C^, B[n1], c[n2]
136#define ARRAY_OUTER_CONJ2(_A, _B, _C, _n1, _n2) \
137 ({ \
138 int _idxA = 0; \
139 for (int _i = 0; _i < (_n1); ++_i) { \
140 for (int _j = 0; _j < (_n2); ++_j) { (_A)[_idxA++] = (_B)[_i] * CONJ_OF((_C)[_j]); } \
141 } \
142 })
143
144// @NOTE: it is well optimized.
145// @brief: A[n1,n3] = B[n1,n2] * C[n2,n3]
146#define ARRAY_MATMUL(_A, _B, _C, _n1, _n2, _n3) \
147 ({ \
148 int _size_A = (_n1) * (_n3); \
149 int _idxB = 0; \
150 ARRAY_CLEAR(_A, _size_A); \
151 for (int _i = 0; _i < (_n1); ++_i) { \
152 int _idxC = 0; \
153 for (int _j = 0; _j < (_n2); ++_j) { \
154 auto _tmp = (_B)[_idxB++]; \
155 int _idxA = _i * (_n3); \
156 for (int _k = 0; _k < (_n3); ++_k) { (_A)[_idxA++] += _tmp * (_C)[_idxC++]; } \
157 } \
158 } \
159 })
160
161// @brief: A[n1,n3] = B^[n1,n2] * C[n2,n3], so B is B[n2,n1]
162#define ARRAY_MATMUL_TRANS1(_A, _B, _C, _n1, _n2, _n3) \
163 ({ \
164 int _size_A = (_n1) * (_n3); \
165 ARRAY_CLEAR(_A, _size_A); \
166 for (int _i = 0; _i < (_n1); ++_i) { \
167 int _idxB = _i; \
168 int _idxC = 0; \
169 for (int _j = 0; _j < (_n2); ++_j) { \
170 auto _tmp = CONJ_OF((_B)[_idxB]); \
171 _idxB += (_n1); \
172 int _idxA = _i * (_n3); \
173 for (int _k = 0; _k < (_n3); ++_k) { (_A)[_idxA++] += _tmp * (_C)[_idxC++]; } \
174 } \
175 } \
176 })
177
178// @brief: A[n1,n3] = B[n1,n2] * C^[n2,n3], so C is C[n3,n2]
179#define ARRAY_MATMUL_TRANS2(_A, _B, _C, _n1, _n2, _n3) \
180 ({ \
181 int _size_A = (_n1) * (_n3); \
182 int _idxB = 0; \
183 ARRAY_CLEAR(_A, _size_A); \
184 for (int _i = 0; _i < (_n1); ++_i) { \
185 for (int _j = 0; _j < (_n2); ++_j) { \
186 auto _tmp = (_B)[_idxB++]; \
187 int _idxC = _j; \
188 int _idxA = _i * (_n3); \
189 for (int _k = 0; _k < (_n3); ++_k) { \
190 (_A)[_idxA++] += _tmp * CONJ_OF((_C)[_idxC]); \
191 _idxC += (_n2); \
192 } \
193 } \
194 } \
195 })
196
197// @brief: A = B^*C*D, here C is diagonal, B[n2,n1], C[n2,n2], D[n2,n3]
198// if C[n2,n2], please set n0 = n2, if C[n2] please set n0 = 0.
199#define ARRAY_MATMUL3_TRANS1(_A, _B, _C, _D, _n1, _n2, _n0, _n3) \
200 ({ \
201 int _size_A = (_n1) * (_n3); \
202 ARRAY_CLEAR(_A, _size_A); \
203 int _iaddC = (_n0) + 1; \
204 for (int _i = 0; _i < (_n1); ++_i) { \
205 int _idxB = _i; \
206 int _idxC = 0, _idxD = 0; \
207 for (int _j = 0; _j < (_n2); ++_j) { \
208 auto _tmp = CONJ_OF((_B)[_idxB]) * (_C)[_idxC]; \
209 _idxB += (_n1); \
210 _idxC += _iaddC; \
211 int _idxA = _i * (_n3); \
212 for (int _k = 0; _k < (_n3); ++_k) { (_A)[_idxA++] += _tmp * (_D)[_idxD++]; } \
213 } \
214 } \
215 })
216
217// @brief: A = B*C*D^, here C is diagonal, B[n1,n2], C[n2, n2], D[n3,n2]
218// if C[n2,n2], please set n0 = n2, if C[n2] please set n0 = 0.
219#define ARRAY_MATMUL3_TRANS2(_A, _B, _C, _D, _n1, _n2, _n0, _n3) \
220 ({ \
221 int _size_A = (_n1) * (_n3); \
222 ARRAY_CLEAR(_A, _size_A); \
223 int _idxB = 0; \
224 int _iaddC = (_n0) + 1; \
225 for (int _i = 0; _i < (_n1); ++_i) { \
226 int _idxC = 0; \
227 for (int _j = 0; _j < (_n2); ++_j) { \
228 auto _tmp = (_B)[_idxB++] * (_C)[_idxC]; \
229 _idxC += _iaddC; \
230 int _idxD = _j; \
231 int _idxA = _i * (_n3); \
232 for (int _k = 0; _k < (_n3); ++_k) { \
233 (_A)[_idxA++] += _tmp * CONJ_OF((_D)[_idxD]); \
234 _idxD += (_n2); \
235 } \
236 } \
237 } \
238 })
239
240#endif // Array_Macro_H