QDP++
sse_blas_vaxmbyz4_double.cc
Go to the documentation of this file.
1
6
8#include <xmmintrin.h>
10
11namespace QDP {
12
13
14#ifndef L2BY2
15#define L2BY2 1365 /* L2 / 2 in SPINORS */
16#endif
17
18
19 void vaxmbyz4(REAL64 *z,REAL64 *a,REAL64 *x, REAL64 *b, REAL64 *y, int n_4vec)
20{
21 __m128d a_sse;
22 __m128d b_sse;
23
24 __m128d tmp1;
25 __m128d tmp2;
26 __m128d tmp3;
27 __m128d tmp4;
28
29 __m128d x1;
30 __m128d y1;
31 __m128d z1;
32
33 __m128d x2;
34 __m128d y2;
35 __m128d z2;
36
37 __m128d x3;
38 __m128d y3;
39 __m128d z3;
40
41 // Load the scalar into low bytes of scalar
42 a_sse = _mm_load_sd(a);
43 b_sse = _mm_load_sd(b);
44
45 // cross components into tmp
46 // Zero tmp
47 tmp1 = _mm_setzero_pd();
48
49 tmp1 = _mm_shuffle_pd(a_sse, a_sse, 0x1);
50 a_sse = _mm_add_pd(a_sse, tmp1);
51
52 tmp2 = _mm_setzero_pd();
53 tmp2 = _mm_shuffle_pd(b_sse, b_sse, 0x1);
54 b_sse = _mm_add_pd(b_sse, tmp2);
55
56
57 // Do n_3vec 3vectors.
58 double *x_p=x;
59 double *y_p=y;
60 double *z_p=z;
61
62 if( n_4vec < L2BY2 ) {
63 PREFETCHW(((const char *)z_p)+8);
64 PREFETCHW(((const char *)z_p)+16);
65 for(int i=0; i < 3*n_4vec; i++) {
66
67 y1 = _mm_load_pd(y_p);
68 x1 = _mm_load_pd(x_p);
69 z1 = _mm_mul_pd(a_sse,x1);
70 tmp1 = _mm_mul_pd(b_sse,y1);
71 z1 = _mm_sub_pd(z1,tmp1);
72 _mm_store_pd(z_p, z1);
73
74 y2 = _mm_load_pd(y_p+2);
75 x2 = _mm_load_pd(x_p+2);
76 z2 = _mm_mul_pd(a_sse,x2);
77 tmp2 = _mm_mul_pd(b_sse,y2);
78 z2 = _mm_sub_pd(z2,tmp2);
79 _mm_store_pd(z_p+2, z2);
80
81
82 y3 = _mm_load_pd(y_p+4);
83 x3 = _mm_load_pd(x_p+4);
84 z3 = _mm_mul_pd(a_sse,x3);
85 tmp3 = _mm_mul_pd(b_sse,y3);
86 z3 = _mm_sub_pd(z3,tmp3);
87 _mm_store_pd(z_p+4, z3);
88
89
90 y1 = _mm_load_pd(y_p+6);
91 x1 = _mm_load_pd(x_p+6);
92 z1 = _mm_mul_pd(a_sse,x1);
93 tmp1 = _mm_mul_pd(b_sse,y1);
94 z1 = _mm_sub_pd(z1,tmp1);
95 _mm_store_pd(z_p+6, z1);
96
97 x_p+=8; y_p+=8; z_p+=8;
98
99 }
100
101 }
102 else {
103
104 for(int i=0; i < 3*n_4vec; i++) {
105
106 PREFETCHNTA(((const char *)x_p)+56);
107 PREFETCHNTA(((const char *)y_p)+56);
108
109 y1 = _mm_load_pd(y_p);
110 x1 = _mm_load_pd(x_p);
111 z1 = _mm_mul_pd(a_sse,x1);
112 tmp1 = _mm_mul_pd(b_sse,y1);
113 z1 = _mm_sub_pd(z1,tmp1);
114 _mm_stream_pd(z_p, z1);
115
116 y2 = _mm_load_pd(y_p+2);
117 x2 = _mm_load_pd(x_p+2);
118 z2 = _mm_mul_pd(a_sse,x2);
119 tmp2 = _mm_mul_pd(b_sse,y2);
120 z2 = _mm_sub_pd(z2,tmp2);
121 _mm_stream_pd(z_p+2, z2);
122
123
124 y3 = _mm_load_pd(y_p+4);
125 x3 = _mm_load_pd(x_p+4);
126 z3 = _mm_mul_pd(a_sse,x3);
127 tmp3 = _mm_mul_pd(b_sse,y3);
128 z3 = _mm_sub_pd(z3,tmp3);
129 _mm_stream_pd(z_p+4, z3);
130
131
132 y1 = _mm_load_pd(y_p+6);
133 x1 = _mm_load_pd(x_p+6);
134 z1 = _mm_mul_pd(a_sse,x1);
135 tmp1 = _mm_mul_pd(b_sse,y1);
136 z1 = _mm_sub_pd(z1,tmp1);
137 _mm_stream_pd(z_p+6, z1);
138
139 x_p+=8; y_p+=8; z_p+=8;
140
141 }
142
143
144 }
145}
146
147 void vaxmby4(REAL64 *y,REAL64 *a,REAL64 *x, REAL64 *b, int n_4vec)
148{
149 __m128d a_sse;
150 __m128d b_sse;
151
152 __m128d tmp1;
153 __m128d tmp2;
154 __m128d tmp3;
155 __m128d tmp4;
156
157 __m128d x1;
158 __m128d y1;
159 __m128d z1;
160
161 __m128d x2;
162 __m128d y2;
163 __m128d z2;
164
165 __m128d x3;
166 __m128d y3;
167 __m128d z3;
168
169 // Load the scalar into low bytes of scalar
170 a_sse = _mm_load_sd(a);
171 b_sse = _mm_load_sd(b);
172
173 // cross components into tmp
174 // Zero tmp
175 tmp1 = _mm_setzero_pd();
176 tmp1 = _mm_shuffle_pd(a_sse, a_sse, 0x1);
177 a_sse = _mm_add_pd(a_sse, tmp1);
178
179 tmp2 = _mm_setzero_pd();
180 tmp2 = _mm_shuffle_pd(b_sse, b_sse, 0x1);
181 b_sse = _mm_add_pd(b_sse, tmp2);
182
183
184 // Do n_3vec 3vectors.
185 double *x_p=x;
186 double *y_p=y;
187
188
189 if( n_4vec < L2BY2 ) {
190 for(int i=0; i < 3*n_4vec; i++) {
191 PREFETCHW(((const char *)y_p)+16);
192
193 y1 = _mm_load_pd(y_p);
194 x1 = _mm_load_pd(x_p);
195 z1 = _mm_mul_pd(a_sse,x1);
196 tmp1 = _mm_mul_pd(b_sse,y1);
197 z1 = _mm_sub_pd(z1,tmp1);
198 _mm_store_pd(y_p, z1);
199
200 y2 = _mm_load_pd(y_p+2);
201 x2 = _mm_load_pd(x_p+2);
202 z2 = _mm_mul_pd(a_sse,x2);
203 tmp2 = _mm_mul_pd(b_sse,y2);
204 z2 = _mm_sub_pd(z2,tmp2);
205 _mm_store_pd(y_p+2, z2);
206
207
208 y3 = _mm_load_pd(y_p+4);
209 x3 = _mm_load_pd(x_p+4);
210 z3 = _mm_mul_pd(a_sse,x3);
211 tmp3 = _mm_mul_pd(b_sse,y3);
212 z3 = _mm_sub_pd(z3,tmp3);
213 _mm_store_pd(y_p+4, z3);
214
215
216 y1 = _mm_load_pd(y_p+6);
217 x1 = _mm_load_pd(x_p+6);
218 z1 = _mm_mul_pd(a_sse,x1);
219 tmp1 = _mm_mul_pd(b_sse,y1);
220 z1 = _mm_sub_pd(z1,tmp1);
221 _mm_store_pd(y_p+6, z1);
222
223 x_p +=8; y_p+=8;
224 }
225
226 }
227 else {
228 for(int i=0; i < 3*n_4vec; i++) {
229
230 PREFETCHNTA(((const char *)y_p)+56);
231 PREFETCHNTA(((const char *)x_p)+56);
232
233 y1 = _mm_load_pd(y_p);
234 x1 = _mm_load_pd(x_p);
235 z1 = _mm_mul_pd(a_sse,x1);
236 tmp1 = _mm_mul_pd(b_sse,y1);
237 z1 = _mm_sub_pd(z1,tmp1);
238 _mm_store_pd(y_p, z1);
239
240 y2 = _mm_load_pd(y_p+2);
241 x2 = _mm_load_pd(x_p+2);
242 z2 = _mm_mul_pd(a_sse,x2);
243 tmp2 = _mm_mul_pd(b_sse,y2);
244 z2 = _mm_sub_pd(z2,tmp2);
245 _mm_store_pd(y_p+2, z2);
246
247
248 y3 = _mm_load_pd(y_p+4);
249 x3 = _mm_load_pd(x_p+4);
250 z3 = _mm_mul_pd(a_sse,x3);
251 tmp3 = _mm_mul_pd(b_sse,y3);
252 z3 = _mm_sub_pd(z3,tmp3);
253 _mm_store_pd(y_p+4, z3);
254
255
256 y1 = _mm_load_pd(y_p+6);
257 x1 = _mm_load_pd(x_p+6);
258 z1 = _mm_mul_pd(a_sse,x1);
259 tmp1 = _mm_mul_pd(b_sse,y1);
260 z1 = _mm_sub_pd(z1,tmp1);
261 _mm_store_pd(y_p+6, z1);
262
263 x_p +=8; y_p+=8;
264 }
265
266 }
267
268}
269
270}
271
272
273
double REAL64
Yet another random number generator.
void vaxmby4(REAL64 *y, REAL64 *a, REAL64 *x, REAL64 *b, int n_4vec)
void vaxmbyz4(REAL64 *z, REAL64 *a, REAL64 *x, REAL64 *b, REAL64 *y, int n_4vec)
#define L2BY2
Generic Scalar VAXPY routine.
#define PREFETCHW(a)
#define PREFETCHNTA(a)
Definition sse_prefetch.h:9