~ubuntu-branches/ubuntu/vivid/atlas/vivid

« back to all changes in this revision

Viewing changes to src/blas/ger/ATL_gerk_Mlt16.c

  • Committer: Package Import Robot
  • Author(s): Sébastien Villemot, Sylvestre Ledru, Sébastien Villemot
  • Date: 2013-06-11 15:58:16 UTC
  • mfrom: (1.1.4) (25 sid)
  • mto: This revision was merged to the branch mainline in revision 26.
  • Revision ID: package-import@ubuntu.com-20130611155816-8xeeiziu1iml040c
Tags: 3.10.1-1
[ Sylvestre Ledru ]
* New upstream release (Closes: #609287)

[ Sébastien Villemot ]
* Provide architectural defaults (i.e. precomputed timings) for all
  release archs (except armel and mips for the time being, due to slow
  porterboxes). This will make the package build much faster and should
  eliminate transient build failures due to excessive variance in the
  timings.
* Move symlinks for lib{cblas,f77blas,atlas,lapack_atlas} out of the
  libblas.so.3 alternative and make them always present, so that
  software relying on these libs do not break when another alternative
  is selected for BLAS
* ATLAS now has improved ARM support with native asm constructs. This required
  the following tunes:
  + armel-is-v4t.diff: new patch, prevents FTBFS on armel; otherwise,
    ATLAS uses asm constructs too recent for the platform (armel is only v4t)
  + debian/rules: on armhf, define the ATL_ARM_HARDFP flag; otherwise the asm
    constructs use the soft-float ABI for passing floating points
  + on armhf, ensure that -mfloat-abi=softfp and -mcpu=vfpv3 flags are never
    used; this is implemented via a patch (armhf.diff) and by the use of fixed
    archdefs
* The generic package is now built without multi-threading, because otherwise
  the package fails to build on some single-processor machines (this required
  the introduction of a patch: fix-non-threaded-build.diff). As a side effect,
  the build of the custom package gracefully handles non-threaded
  builds. (Closes: #602524)
* Add libblas.a as slave in the libblas.so alternative (Closes: #701921)
* Add symlinks for lib{f77blas,atlas}.a in /usr/lib (Closes: #666203)
* Modify shlibs file of libatlas3-base, such that packages using
  libblas/liblapack depend on any BLAS/LAPACK alternative, while packages
  depending on ATLAS-specific libraries (e.g. libatlas.so) depend specifically
  on libatlas3-base.
* corei1.diff: remove patch, applied upstream
* Use my @debian.org email address
* Remove obsolete DM-Upload-Allowed flag
* Switch VCS to git
* Remove Conflicts/Replaces against pre-squeeze packages
* libatlas-base-dev now provides libblas.so, as libblas-dev
* No longer use -Wa,--noexecstack in CFLAGS, it makes the package FTBFS
* Do not use POWER3 arch for powerpcspe port (Closes: #701068)
* Bump to debhelper compat level 9
* README.Debian: mention that devscripts is needed to compile the custom
  package (Closes: #697431)
* Bump Standards-Version to 3.9.4. As a consequence, add Built-Using
  fields because the package embeds stuff from liblapack-pic

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
#include "atlas_misc.h"
 
2
#include "atlas_lvl2.h"
 
3
#ifdef TCPLX
 
4
   #include "atlas_reflevel2.h"
 
5
#endif
 
6
 
 
7
typedef void (*gerk_t)
 
8
   (const int M, const int N, const SCALAR alpha, const TYPE *X, const int incX,
 
9
    const TYPE *Y, const int incY, TYPE *A, const int lda);
 
10
 
 
11
#ifdef TREAL
 
12
 
 
13
static void ATL_gerk_Meq1
 
14
   (const int M, const int N, const SCALAR alpha, const TYPE *X, const int incX,
 
15
    const TYPE *Y, const int incY, TYPE *A, const int lda)
 
16
{
 
17
   register int j;
 
18
   register TYPE y0, x0;
 
19
 
 
20
   if (alpha == ATL_rone) goto ALPHA1;
 
21
   else if (alpha != ATL_rnone) goto ALPHAX;
 
22
   x0 = -X[0*incX];
 
23
X_IS_LOADED:
 
24
 
 
25
   for (j=0; j < N; j++, A += lda)
 
26
   {
 
27
      y0 = Y[incY*j];
 
28
      A[0] += x0 * y0;
 
29
   }
 
30
   return;
 
31
 
 
32
ALPHAX:
 
33
   y0 = alpha;
 
34
   x0 = X[0*incX] * y0;
 
35
   goto X_IS_LOADED;
 
36
ALPHA1:
 
37
   x0 = X[0*incX];
 
38
   goto X_IS_LOADED;
 
39
}
 
40
static void ATL_gerk_Meq2
 
41
   (const int M, const int N, const SCALAR alpha, const TYPE *X, const int incX,
 
42
    const TYPE *Y, const int incY, TYPE *A, const int lda)
 
43
{
 
44
   register int j;
 
45
   register TYPE y0, x0, x1;
 
46
 
 
47
   if (alpha == ATL_rone) goto ALPHA1;
 
48
   else if (alpha != ATL_rnone) goto ALPHAX;
 
49
   x0 = -X[0*incX];
 
50
   x1 = -X[1*incX];
 
51
X_IS_LOADED:
 
52
 
 
53
   for (j=0; j < N; j++, A += lda)
 
54
   {
 
55
      y0 = Y[incY*j];
 
56
      A[0] += x0 * y0;
 
57
      A[1] += x1 * y0;
 
58
   }
 
59
   return;
 
60
 
 
61
ALPHAX:
 
62
   y0 = alpha;
 
63
   x0 = X[0*incX] * y0;
 
64
   y0 = alpha;
 
65
   x1 = X[1*incX] * y0;
 
66
   goto X_IS_LOADED;
 
67
ALPHA1:
 
68
   x0 = X[0*incX];
 
69
   x1 = X[1*incX];
 
70
   goto X_IS_LOADED;
 
71
}
 
72
static void ATL_gerk_Meq3
 
73
   (const int M, const int N, const SCALAR alpha, const TYPE *X, const int incX,
 
74
    const TYPE *Y, const int incY, TYPE *A, const int lda)
 
75
{
 
76
   register int j;
 
77
   register TYPE y0, x0, x1, x2;
 
78
 
 
79
   if (alpha == ATL_rone) goto ALPHA1;
 
80
   else if (alpha != ATL_rnone) goto ALPHAX;
 
81
   x0 = -X[0*incX];
 
82
   x1 = -X[1*incX];
 
83
   x2 = -X[2*incX];
 
84
X_IS_LOADED:
 
85
 
 
86
   for (j=0; j < N; j++, A += lda)
 
87
   {
 
88
      y0 = Y[incY*j];
 
89
      A[0] += x0 * y0;
 
90
      A[1] += x1 * y0;
 
91
      A[2] += x2 * y0;
 
92
   }
 
93
   return;
 
94
 
 
95
ALPHAX:
 
96
   y0 = alpha;
 
97
   x0 = X[0*incX] * y0;
 
98
   y0 = alpha;
 
99
   x1 = X[1*incX] * y0;
 
100
   y0 = alpha;
 
101
   x2 = X[2*incX] * y0;
 
102
   goto X_IS_LOADED;
 
103
ALPHA1:
 
104
   x0 = X[0*incX];
 
105
   x1 = X[1*incX];
 
106
   x2 = X[2*incX];
 
107
   goto X_IS_LOADED;
 
108
}
 
109
static void ATL_gerk_Meq4
 
110
   (const int M, const int N, const SCALAR alpha, const TYPE *X, const int incX,
 
111
    const TYPE *Y, const int incY, TYPE *A, const int lda)
 
112
{
 
113
   register int j;
 
114
   register TYPE y0, x0, x1, x2, x3;
 
115
 
 
116
   if (alpha == ATL_rone) goto ALPHA1;
 
117
   else if (alpha != ATL_rnone) goto ALPHAX;
 
118
   x0 = -X[0*incX];
 
119
   x1 = -X[1*incX];
 
120
   x2 = -X[2*incX];
 
121
   x3 = -X[3*incX];
 
122
X_IS_LOADED:
 
123
 
 
124
   for (j=0; j < N; j++, A += lda)
 
125
   {
 
126
      y0 = Y[incY*j];
 
127
      A[0] += x0 * y0;
 
128
      A[1] += x1 * y0;
 
129
      A[2] += x2 * y0;
 
130
      A[3] += x3 * y0;
 
131
   }
 
132
   return;
 
133
 
 
134
ALPHAX:
 
135
   y0 = alpha;
 
136
   x0 = X[0*incX] * y0;
 
137
   y0 = alpha;
 
138
   x1 = X[1*incX] * y0;
 
139
   y0 = alpha;
 
140
   x2 = X[2*incX] * y0;
 
141
   y0 = alpha;
 
142
   x3 = X[3*incX] * y0;
 
143
   goto X_IS_LOADED;
 
144
ALPHA1:
 
145
   x0 = X[0*incX];
 
146
   x1 = X[1*incX];
 
147
   x2 = X[2*incX];
 
148
   x3 = X[3*incX];
 
149
   goto X_IS_LOADED;
 
150
}
 
151
static void ATL_gerk_Meq5
 
152
   (const int M, const int N, const SCALAR alpha, const TYPE *X, const int incX,
 
153
    const TYPE *Y, const int incY, TYPE *A, const int lda)
 
154
{
 
155
   register int j;
 
156
   register TYPE y0, x0, x1, x2, x3, x4;
 
157
 
 
158
   if (alpha == ATL_rone) goto ALPHA1;
 
159
   else if (alpha != ATL_rnone) goto ALPHAX;
 
160
   x0 = -X[0*incX];
 
161
   x1 = -X[1*incX];
 
162
   x2 = -X[2*incX];
 
163
   x3 = -X[3*incX];
 
164
   x4 = -X[4*incX];
 
165
X_IS_LOADED:
 
166
 
 
167
   for (j=0; j < N; j++, A += lda)
 
168
   {
 
169
      y0 = Y[incY*j];
 
170
      A[0] += x0 * y0;
 
171
      A[1] += x1 * y0;
 
172
      A[2] += x2 * y0;
 
173
      A[3] += x3 * y0;
 
174
      A[4] += x4 * y0;
 
175
   }
 
176
   return;
 
177
 
 
178
ALPHAX:
 
179
   y0 = alpha;
 
180
   x0 = X[0*incX] * y0;
 
181
   y0 = alpha;
 
182
   x1 = X[1*incX] * y0;
 
183
   y0 = alpha;
 
184
   x2 = X[2*incX] * y0;
 
185
   y0 = alpha;
 
186
   x3 = X[3*incX] * y0;
 
187
   y0 = alpha;
 
188
   x4 = X[4*incX] * y0;
 
189
   goto X_IS_LOADED;
 
190
ALPHA1:
 
191
   x0 = X[0*incX];
 
192
   x1 = X[1*incX];
 
193
   x2 = X[2*incX];
 
194
   x3 = X[3*incX];
 
195
   x4 = X[4*incX];
 
196
   goto X_IS_LOADED;
 
197
}
 
198
static void ATL_gerk_Meq6
 
199
   (const int M, const int N, const SCALAR alpha, const TYPE *X, const int incX,
 
200
    const TYPE *Y, const int incY, TYPE *A, const int lda)
 
201
{
 
202
   register int j;
 
203
   register TYPE y0, x0, x1, x2, x3, x4, x5;
 
204
 
 
205
   if (alpha == ATL_rone) goto ALPHA1;
 
206
   else if (alpha != ATL_rnone) goto ALPHAX;
 
207
   x0 = -X[0*incX];
 
208
   x1 = -X[1*incX];
 
209
   x2 = -X[2*incX];
 
210
   x3 = -X[3*incX];
 
211
   x4 = -X[4*incX];
 
212
   x5 = -X[5*incX];
 
213
X_IS_LOADED:
 
214
 
 
215
   for (j=0; j < N; j++, A += lda)
 
216
   {
 
217
      y0 = Y[incY*j];
 
218
      A[0] += x0 * y0;
 
219
      A[1] += x1 * y0;
 
220
      A[2] += x2 * y0;
 
221
      A[3] += x3 * y0;
 
222
      A[4] += x4 * y0;
 
223
      A[5] += x5 * y0;
 
224
   }
 
225
   return;
 
226
 
 
227
ALPHAX:
 
228
   y0 = alpha;
 
229
   x0 = X[0*incX] * y0;
 
230
   y0 = alpha;
 
231
   x1 = X[1*incX] * y0;
 
232
   y0 = alpha;
 
233
   x2 = X[2*incX] * y0;
 
234
   y0 = alpha;
 
235
   x3 = X[3*incX] * y0;
 
236
   y0 = alpha;
 
237
   x4 = X[4*incX] * y0;
 
238
   y0 = alpha;
 
239
   x5 = X[5*incX] * y0;
 
240
   goto X_IS_LOADED;
 
241
ALPHA1:
 
242
   x0 = X[0*incX];
 
243
   x1 = X[1*incX];
 
244
   x2 = X[2*incX];
 
245
   x3 = X[3*incX];
 
246
   x4 = X[4*incX];
 
247
   x5 = X[5*incX];
 
248
   goto X_IS_LOADED;
 
249
}
 
250
static void ATL_gerk_Meq7
 
251
   (const int M, const int N, const SCALAR alpha, const TYPE *X, const int incX,
 
252
    const TYPE *Y, const int incY, TYPE *A, const int lda)
 
253
{
 
254
   register int j;
 
255
   register TYPE y0, x0, x1, x2, x3, x4, x5, x6;
 
256
 
 
257
   if (alpha == ATL_rone) goto ALPHA1;
 
258
   else if (alpha != ATL_rnone) goto ALPHAX;
 
259
   x0 = -X[0*incX];
 
260
   x1 = -X[1*incX];
 
261
   x2 = -X[2*incX];
 
262
   x3 = -X[3*incX];
 
263
   x4 = -X[4*incX];
 
264
   x5 = -X[5*incX];
 
265
   x6 = -X[6*incX];
 
266
X_IS_LOADED:
 
267
 
 
268
   for (j=0; j < N; j++, A += lda)
 
269
   {
 
270
      y0 = Y[incY*j];
 
271
      A[0] += x0 * y0;
 
272
      A[1] += x1 * y0;
 
273
      A[2] += x2 * y0;
 
274
      A[3] += x3 * y0;
 
275
      A[4] += x4 * y0;
 
276
      A[5] += x5 * y0;
 
277
      A[6] += x6 * y0;
 
278
   }
 
279
   return;
 
280
 
 
281
ALPHAX:
 
282
   y0 = alpha;
 
283
   x0 = X[0*incX] * y0;
 
284
   y0 = alpha;
 
285
   x1 = X[1*incX] * y0;
 
286
   y0 = alpha;
 
287
   x2 = X[2*incX] * y0;
 
288
   y0 = alpha;
 
289
   x3 = X[3*incX] * y0;
 
290
   y0 = alpha;
 
291
   x4 = X[4*incX] * y0;
 
292
   y0 = alpha;
 
293
   x5 = X[5*incX] * y0;
 
294
   y0 = alpha;
 
295
   x6 = X[6*incX] * y0;
 
296
   goto X_IS_LOADED;
 
297
ALPHA1:
 
298
   x0 = X[0*incX];
 
299
   x1 = X[1*incX];
 
300
   x2 = X[2*incX];
 
301
   x3 = X[3*incX];
 
302
   x4 = X[4*incX];
 
303
   x5 = X[5*incX];
 
304
   x6 = X[6*incX];
 
305
   goto X_IS_LOADED;
 
306
}
 
307
static void ATL_gerk_Meq8
 
308
   (const int M, const int N, const SCALAR alpha, const TYPE *X, const int incX,
 
309
    const TYPE *Y, const int incY, TYPE *A, const int lda)
 
310
{
 
311
   register int j;
 
312
   register TYPE y0, x0, x1, x2, x3, x4, x5, x6, x7;
 
313
 
 
314
   if (alpha == ATL_rone) goto ALPHA1;
 
315
   else if (alpha != ATL_rnone) goto ALPHAX;
 
316
   x0 = -X[0*incX];
 
317
   x1 = -X[1*incX];
 
318
   x2 = -X[2*incX];
 
319
   x3 = -X[3*incX];
 
320
   x4 = -X[4*incX];
 
321
   x5 = -X[5*incX];
 
322
   x6 = -X[6*incX];
 
323
   x7 = -X[7*incX];
 
324
X_IS_LOADED:
 
325
 
 
326
   for (j=0; j < N; j++, A += lda)
 
327
   {
 
328
      y0 = Y[incY*j];
 
329
      A[0] += x0 * y0;
 
330
      A[1] += x1 * y0;
 
331
      A[2] += x2 * y0;
 
332
      A[3] += x3 * y0;
 
333
      A[4] += x4 * y0;
 
334
      A[5] += x5 * y0;
 
335
      A[6] += x6 * y0;
 
336
      A[7] += x7 * y0;
 
337
   }
 
338
   return;
 
339
 
 
340
ALPHAX:
 
341
   y0 = alpha;
 
342
   x0 = X[0*incX] * y0;
 
343
   y0 = alpha;
 
344
   x1 = X[1*incX] * y0;
 
345
   y0 = alpha;
 
346
   x2 = X[2*incX] * y0;
 
347
   y0 = alpha;
 
348
   x3 = X[3*incX] * y0;
 
349
   y0 = alpha;
 
350
   x4 = X[4*incX] * y0;
 
351
   y0 = alpha;
 
352
   x5 = X[5*incX] * y0;
 
353
   y0 = alpha;
 
354
   x6 = X[6*incX] * y0;
 
355
   y0 = alpha;
 
356
   x7 = X[7*incX] * y0;
 
357
   goto X_IS_LOADED;
 
358
ALPHA1:
 
359
   x0 = X[0*incX];
 
360
   x1 = X[1*incX];
 
361
   x2 = X[2*incX];
 
362
   x3 = X[3*incX];
 
363
   x4 = X[4*incX];
 
364
   x5 = X[5*incX];
 
365
   x6 = X[6*incX];
 
366
   x7 = X[7*incX];
 
367
   goto X_IS_LOADED;
 
368
}
 
369
static void ATL_gerk_Meq9
 
370
   (const int M, const int N, const SCALAR alpha, const TYPE *X, const int incX,
 
371
    const TYPE *Y, const int incY, TYPE *A, const int lda)
 
372
{
 
373
   register int j;
 
374
   register TYPE y0, x0, x1, x2, x3, x4, x5, x6, x7, x8;
 
375
 
 
376
   if (alpha == ATL_rone) goto ALPHA1;
 
377
   else if (alpha != ATL_rnone) goto ALPHAX;
 
378
   x0 = -X[0*incX];
 
379
   x1 = -X[1*incX];
 
380
   x2 = -X[2*incX];
 
381
   x3 = -X[3*incX];
 
382
   x4 = -X[4*incX];
 
383
   x5 = -X[5*incX];
 
384
   x6 = -X[6*incX];
 
385
   x7 = -X[7*incX];
 
386
   x8 = -X[8*incX];
 
387
X_IS_LOADED:
 
388
 
 
389
   for (j=0; j < N; j++, A += lda)
 
390
   {
 
391
      y0 = Y[incY*j];
 
392
      A[0] += x0 * y0;
 
393
      A[1] += x1 * y0;
 
394
      A[2] += x2 * y0;
 
395
      A[3] += x3 * y0;
 
396
      A[4] += x4 * y0;
 
397
      A[5] += x5 * y0;
 
398
      A[6] += x6 * y0;
 
399
      A[7] += x7 * y0;
 
400
      A[8] += x8 * y0;
 
401
   }
 
402
   return;
 
403
 
 
404
ALPHAX:
 
405
   y0 = alpha;
 
406
   x0 = X[0*incX] * y0;
 
407
   y0 = alpha;
 
408
   x1 = X[1*incX] * y0;
 
409
   y0 = alpha;
 
410
   x2 = X[2*incX] * y0;
 
411
   y0 = alpha;
 
412
   x3 = X[3*incX] * y0;
 
413
   y0 = alpha;
 
414
   x4 = X[4*incX] * y0;
 
415
   y0 = alpha;
 
416
   x5 = X[5*incX] * y0;
 
417
   y0 = alpha;
 
418
   x6 = X[6*incX] * y0;
 
419
   y0 = alpha;
 
420
   x7 = X[7*incX] * y0;
 
421
   y0 = alpha;
 
422
   x8 = X[8*incX] * y0;
 
423
   goto X_IS_LOADED;
 
424
ALPHA1:
 
425
   x0 = X[0*incX];
 
426
   x1 = X[1*incX];
 
427
   x2 = X[2*incX];
 
428
   x3 = X[3*incX];
 
429
   x4 = X[4*incX];
 
430
   x5 = X[5*incX];
 
431
   x6 = X[6*incX];
 
432
   x7 = X[7*incX];
 
433
   x8 = X[8*incX];
 
434
   goto X_IS_LOADED;
 
435
}
 
436
static void ATL_gerk_Meq10
 
437
   (const int M, const int N, const SCALAR alpha, const TYPE *X, const int incX,
 
438
    const TYPE *Y, const int incY, TYPE *A, const int lda)
 
439
{
 
440
   register int j;
 
441
   register TYPE y0, x0, x1, x2, x3, x4, x5, x6, x7, x8, x9;
 
442
 
 
443
   if (alpha == ATL_rone) goto ALPHA1;
 
444
   else if (alpha != ATL_rnone) goto ALPHAX;
 
445
   x0 = -X[0*incX];
 
446
   x1 = -X[1*incX];
 
447
   x2 = -X[2*incX];
 
448
   x3 = -X[3*incX];
 
449
   x4 = -X[4*incX];
 
450
   x5 = -X[5*incX];
 
451
   x6 = -X[6*incX];
 
452
   x7 = -X[7*incX];
 
453
   x8 = -X[8*incX];
 
454
   x9 = -X[9*incX];
 
455
X_IS_LOADED:
 
456
 
 
457
   for (j=0; j < N; j++, A += lda)
 
458
   {
 
459
      y0 = Y[incY*j];
 
460
      A[0] += x0 * y0;
 
461
      A[1] += x1 * y0;
 
462
      A[2] += x2 * y0;
 
463
      A[3] += x3 * y0;
 
464
      A[4] += x4 * y0;
 
465
      A[5] += x5 * y0;
 
466
      A[6] += x6 * y0;
 
467
      A[7] += x7 * y0;
 
468
      A[8] += x8 * y0;
 
469
      A[9] += x9 * y0;
 
470
   }
 
471
   return;
 
472
 
 
473
ALPHAX:
 
474
   y0 = alpha;
 
475
   x0 = X[0*incX] * y0;
 
476
   y0 = alpha;
 
477
   x1 = X[1*incX] * y0;
 
478
   y0 = alpha;
 
479
   x2 = X[2*incX] * y0;
 
480
   y0 = alpha;
 
481
   x3 = X[3*incX] * y0;
 
482
   y0 = alpha;
 
483
   x4 = X[4*incX] * y0;
 
484
   y0 = alpha;
 
485
   x5 = X[5*incX] * y0;
 
486
   y0 = alpha;
 
487
   x6 = X[6*incX] * y0;
 
488
   y0 = alpha;
 
489
   x7 = X[7*incX] * y0;
 
490
   y0 = alpha;
 
491
   x8 = X[8*incX] * y0;
 
492
   y0 = alpha;
 
493
   x9 = X[9*incX] * y0;
 
494
   goto X_IS_LOADED;
 
495
ALPHA1:
 
496
   x0 = X[0*incX];
 
497
   x1 = X[1*incX];
 
498
   x2 = X[2*incX];
 
499
   x3 = X[3*incX];
 
500
   x4 = X[4*incX];
 
501
   x5 = X[5*incX];
 
502
   x6 = X[6*incX];
 
503
   x7 = X[7*incX];
 
504
   x8 = X[8*incX];
 
505
   x9 = X[9*incX];
 
506
   goto X_IS_LOADED;
 
507
}
 
508
static void ATL_gerk_Meq11
 
509
   (const int M, const int N, const SCALAR alpha, const TYPE *X, const int incX,
 
510
    const TYPE *Y, const int incY, TYPE *A, const int lda)
 
511
{
 
512
   register int j;
 
513
   register TYPE y0, x0, x1, x2, x3, x4, x5, x6, x7, x8, x9, x10;
 
514
 
 
515
   if (alpha == ATL_rone) goto ALPHA1;
 
516
   else if (alpha != ATL_rnone) goto ALPHAX;
 
517
   x0 = -X[0*incX];
 
518
   x1 = -X[1*incX];
 
519
   x2 = -X[2*incX];
 
520
   x3 = -X[3*incX];
 
521
   x4 = -X[4*incX];
 
522
   x5 = -X[5*incX];
 
523
   x6 = -X[6*incX];
 
524
   x7 = -X[7*incX];
 
525
   x8 = -X[8*incX];
 
526
   x9 = -X[9*incX];
 
527
   x10 = -X[10*incX];
 
528
X_IS_LOADED:
 
529
 
 
530
   for (j=0; j < N; j++, A += lda)
 
531
   {
 
532
      y0 = Y[incY*j];
 
533
      A[0] += x0 * y0;
 
534
      A[1] += x1 * y0;
 
535
      A[2] += x2 * y0;
 
536
      A[3] += x3 * y0;
 
537
      A[4] += x4 * y0;
 
538
      A[5] += x5 * y0;
 
539
      A[6] += x6 * y0;
 
540
      A[7] += x7 * y0;
 
541
      A[8] += x8 * y0;
 
542
      A[9] += x9 * y0;
 
543
      A[10] += x10 * y0;
 
544
   }
 
545
   return;
 
546
 
 
547
ALPHAX:
 
548
   y0 = alpha;
 
549
   x0 = X[0*incX] * y0;
 
550
   y0 = alpha;
 
551
   x1 = X[1*incX] * y0;
 
552
   y0 = alpha;
 
553
   x2 = X[2*incX] * y0;
 
554
   y0 = alpha;
 
555
   x3 = X[3*incX] * y0;
 
556
   y0 = alpha;
 
557
   x4 = X[4*incX] * y0;
 
558
   y0 = alpha;
 
559
   x5 = X[5*incX] * y0;
 
560
   y0 = alpha;
 
561
   x6 = X[6*incX] * y0;
 
562
   y0 = alpha;
 
563
   x7 = X[7*incX] * y0;
 
564
   y0 = alpha;
 
565
   x8 = X[8*incX] * y0;
 
566
   y0 = alpha;
 
567
   x9 = X[9*incX] * y0;
 
568
   y0 = alpha;
 
569
   x10 = X[10*incX] * y0;
 
570
   goto X_IS_LOADED;
 
571
ALPHA1:
 
572
   x0 = X[0*incX];
 
573
   x1 = X[1*incX];
 
574
   x2 = X[2*incX];
 
575
   x3 = X[3*incX];
 
576
   x4 = X[4*incX];
 
577
   x5 = X[5*incX];
 
578
   x6 = X[6*incX];
 
579
   x7 = X[7*incX];
 
580
   x8 = X[8*incX];
 
581
   x9 = X[9*incX];
 
582
   x10 = X[10*incX];
 
583
   goto X_IS_LOADED;
 
584
}
 
585
static void ATL_gerk_Meq12
 
586
   (const int M, const int N, const SCALAR alpha, const TYPE *X, const int incX,
 
587
    const TYPE *Y, const int incY, TYPE *A, const int lda)
 
588
{
 
589
   register int j;
 
590
   register TYPE y0, x0, x1, x2, x3, x4, x5, x6, x7, x8, x9, x10, x11;
 
591
 
 
592
   if (alpha == ATL_rone) goto ALPHA1;
 
593
   else if (alpha != ATL_rnone) goto ALPHAX;
 
594
   x0 = -X[0*incX];
 
595
   x1 = -X[1*incX];
 
596
   x2 = -X[2*incX];
 
597
   x3 = -X[3*incX];
 
598
   x4 = -X[4*incX];
 
599
   x5 = -X[5*incX];
 
600
   x6 = -X[6*incX];
 
601
   x7 = -X[7*incX];
 
602
   x8 = -X[8*incX];
 
603
   x9 = -X[9*incX];
 
604
   x10 = -X[10*incX];
 
605
   x11 = -X[11*incX];
 
606
X_IS_LOADED:
 
607
 
 
608
   for (j=0; j < N; j++, A += lda)
 
609
   {
 
610
      y0 = Y[incY*j];
 
611
      A[0] += x0 * y0;
 
612
      A[1] += x1 * y0;
 
613
      A[2] += x2 * y0;
 
614
      A[3] += x3 * y0;
 
615
      A[4] += x4 * y0;
 
616
      A[5] += x5 * y0;
 
617
      A[6] += x6 * y0;
 
618
      A[7] += x7 * y0;
 
619
      A[8] += x8 * y0;
 
620
      A[9] += x9 * y0;
 
621
      A[10] += x10 * y0;
 
622
      A[11] += x11 * y0;
 
623
   }
 
624
   return;
 
625
 
 
626
ALPHAX:
 
627
   y0 = alpha;
 
628
   x0 = X[0*incX] * y0;
 
629
   y0 = alpha;
 
630
   x1 = X[1*incX] * y0;
 
631
   y0 = alpha;
 
632
   x2 = X[2*incX] * y0;
 
633
   y0 = alpha;
 
634
   x3 = X[3*incX] * y0;
 
635
   y0 = alpha;
 
636
   x4 = X[4*incX] * y0;
 
637
   y0 = alpha;
 
638
   x5 = X[5*incX] * y0;
 
639
   y0 = alpha;
 
640
   x6 = X[6*incX] * y0;
 
641
   y0 = alpha;
 
642
   x7 = X[7*incX] * y0;
 
643
   y0 = alpha;
 
644
   x8 = X[8*incX] * y0;
 
645
   y0 = alpha;
 
646
   x9 = X[9*incX] * y0;
 
647
   y0 = alpha;
 
648
   x10 = X[10*incX] * y0;
 
649
   y0 = alpha;
 
650
   x11 = X[11*incX] * y0;
 
651
   goto X_IS_LOADED;
 
652
ALPHA1:
 
653
   x0 = X[0*incX];
 
654
   x1 = X[1*incX];
 
655
   x2 = X[2*incX];
 
656
   x3 = X[3*incX];
 
657
   x4 = X[4*incX];
 
658
   x5 = X[5*incX];
 
659
   x6 = X[6*incX];
 
660
   x7 = X[7*incX];
 
661
   x8 = X[8*incX];
 
662
   x9 = X[9*incX];
 
663
   x10 = X[10*incX];
 
664
   x11 = X[11*incX];
 
665
   goto X_IS_LOADED;
 
666
}
 
667
static void ATL_gerk_Meq13
 
668
   (const int M, const int N, const SCALAR alpha, const TYPE *X, const int incX,
 
669
    const TYPE *Y, const int incY, TYPE *A, const int lda)
 
670
{
 
671
   register int j;
 
672
   register TYPE y0, x0, x1, x2, x3, x4, x5, x6, x7, x8, x9, x10, x11, x12;
 
673
 
 
674
   if (alpha == ATL_rone) goto ALPHA1;
 
675
   else if (alpha != ATL_rnone) goto ALPHAX;
 
676
   x0 = -X[0*incX];
 
677
   x1 = -X[1*incX];
 
678
   x2 = -X[2*incX];
 
679
   x3 = -X[3*incX];
 
680
   x4 = -X[4*incX];
 
681
   x5 = -X[5*incX];
 
682
   x6 = -X[6*incX];
 
683
   x7 = -X[7*incX];
 
684
   x8 = -X[8*incX];
 
685
   x9 = -X[9*incX];
 
686
   x10 = -X[10*incX];
 
687
   x11 = -X[11*incX];
 
688
   x12 = -X[12*incX];
 
689
X_IS_LOADED:
 
690
 
 
691
   for (j=0; j < N; j++, A += lda)
 
692
   {
 
693
      y0 = Y[incY*j];
 
694
      A[0] += x0 * y0;
 
695
      A[1] += x1 * y0;
 
696
      A[2] += x2 * y0;
 
697
      A[3] += x3 * y0;
 
698
      A[4] += x4 * y0;
 
699
      A[5] += x5 * y0;
 
700
      A[6] += x6 * y0;
 
701
      A[7] += x7 * y0;
 
702
      A[8] += x8 * y0;
 
703
      A[9] += x9 * y0;
 
704
      A[10] += x10 * y0;
 
705
      A[11] += x11 * y0;
 
706
      A[12] += x12 * y0;
 
707
   }
 
708
   return;
 
709
 
 
710
ALPHAX:
 
711
   y0 = alpha;
 
712
   x0 = X[0*incX] * y0;
 
713
   y0 = alpha;
 
714
   x1 = X[1*incX] * y0;
 
715
   y0 = alpha;
 
716
   x2 = X[2*incX] * y0;
 
717
   y0 = alpha;
 
718
   x3 = X[3*incX] * y0;
 
719
   y0 = alpha;
 
720
   x4 = X[4*incX] * y0;
 
721
   y0 = alpha;
 
722
   x5 = X[5*incX] * y0;
 
723
   y0 = alpha;
 
724
   x6 = X[6*incX] * y0;
 
725
   y0 = alpha;
 
726
   x7 = X[7*incX] * y0;
 
727
   y0 = alpha;
 
728
   x8 = X[8*incX] * y0;
 
729
   y0 = alpha;
 
730
   x9 = X[9*incX] * y0;
 
731
   y0 = alpha;
 
732
   x10 = X[10*incX] * y0;
 
733
   y0 = alpha;
 
734
   x11 = X[11*incX] * y0;
 
735
   y0 = alpha;
 
736
   x12 = X[12*incX] * y0;
 
737
   goto X_IS_LOADED;
 
738
ALPHA1:
 
739
   x0 = X[0*incX];
 
740
   x1 = X[1*incX];
 
741
   x2 = X[2*incX];
 
742
   x3 = X[3*incX];
 
743
   x4 = X[4*incX];
 
744
   x5 = X[5*incX];
 
745
   x6 = X[6*incX];
 
746
   x7 = X[7*incX];
 
747
   x8 = X[8*incX];
 
748
   x9 = X[9*incX];
 
749
   x10 = X[10*incX];
 
750
   x11 = X[11*incX];
 
751
   x12 = X[12*incX];
 
752
   goto X_IS_LOADED;
 
753
}
 
754
static void ATL_gerk_Meq14
 
755
   (const int M, const int N, const SCALAR alpha, const TYPE *X, const int incX,
 
756
    const TYPE *Y, const int incY, TYPE *A, const int lda)
 
757
{
 
758
   register int j;
 
759
   register TYPE y0, x0, x1, x2, x3, x4, x5, x6, x7, x8, x9, x10, x11, x12;
 
760
   register TYPE x13;
 
761
 
 
762
   if (alpha == ATL_rone) goto ALPHA1;
 
763
   else if (alpha != ATL_rnone) goto ALPHAX;
 
764
   x0 = -X[0*incX];
 
765
   x1 = -X[1*incX];
 
766
   x2 = -X[2*incX];
 
767
   x3 = -X[3*incX];
 
768
   x4 = -X[4*incX];
 
769
   x5 = -X[5*incX];
 
770
   x6 = -X[6*incX];
 
771
   x7 = -X[7*incX];
 
772
   x8 = -X[8*incX];
 
773
   x9 = -X[9*incX];
 
774
   x10 = -X[10*incX];
 
775
   x11 = -X[11*incX];
 
776
   x12 = -X[12*incX];
 
777
   x13 = -X[13*incX];
 
778
X_IS_LOADED:
 
779
 
 
780
   for (j=0; j < N; j++, A += lda)
 
781
   {
 
782
      y0 = Y[incY*j];
 
783
      A[0] += x0 * y0;
 
784
      A[1] += x1 * y0;
 
785
      A[2] += x2 * y0;
 
786
      A[3] += x3 * y0;
 
787
      A[4] += x4 * y0;
 
788
      A[5] += x5 * y0;
 
789
      A[6] += x6 * y0;
 
790
      A[7] += x7 * y0;
 
791
      A[8] += x8 * y0;
 
792
      A[9] += x9 * y0;
 
793
      A[10] += x10 * y0;
 
794
      A[11] += x11 * y0;
 
795
      A[12] += x12 * y0;
 
796
      A[13] += x13 * y0;
 
797
   }
 
798
   return;
 
799
 
 
800
ALPHAX:
 
801
   y0 = alpha;
 
802
   x0 = X[0*incX] * y0;
 
803
   y0 = alpha;
 
804
   x1 = X[1*incX] * y0;
 
805
   y0 = alpha;
 
806
   x2 = X[2*incX] * y0;
 
807
   y0 = alpha;
 
808
   x3 = X[3*incX] * y0;
 
809
   y0 = alpha;
 
810
   x4 = X[4*incX] * y0;
 
811
   y0 = alpha;
 
812
   x5 = X[5*incX] * y0;
 
813
   y0 = alpha;
 
814
   x6 = X[6*incX] * y0;
 
815
   y0 = alpha;
 
816
   x7 = X[7*incX] * y0;
 
817
   y0 = alpha;
 
818
   x8 = X[8*incX] * y0;
 
819
   y0 = alpha;
 
820
   x9 = X[9*incX] * y0;
 
821
   y0 = alpha;
 
822
   x10 = X[10*incX] * y0;
 
823
   y0 = alpha;
 
824
   x11 = X[11*incX] * y0;
 
825
   y0 = alpha;
 
826
   x12 = X[12*incX] * y0;
 
827
   y0 = alpha;
 
828
   x13 = X[13*incX] * y0;
 
829
   goto X_IS_LOADED;
 
830
ALPHA1:
 
831
   x0 = X[0*incX];
 
832
   x1 = X[1*incX];
 
833
   x2 = X[2*incX];
 
834
   x3 = X[3*incX];
 
835
   x4 = X[4*incX];
 
836
   x5 = X[5*incX];
 
837
   x6 = X[6*incX];
 
838
   x7 = X[7*incX];
 
839
   x8 = X[8*incX];
 
840
   x9 = X[9*incX];
 
841
   x10 = X[10*incX];
 
842
   x11 = X[11*incX];
 
843
   x12 = X[12*incX];
 
844
   x13 = X[13*incX];
 
845
   goto X_IS_LOADED;
 
846
}
 
847
static void ATL_gerk_Meq15
 
848
   (const int M, const int N, const SCALAR alpha, const TYPE *X, const int incX,
 
849
    const TYPE *Y, const int incY, TYPE *A, const int lda)
 
850
{
 
851
   register int j;
 
852
   register TYPE y0, x0, x1, x2, x3, x4, x5, x6, x7, x8, x9, x10, x11, x12;
 
853
   register TYPE x13, x14;
 
854
 
 
855
   if (alpha == ATL_rone) goto ALPHA1;
 
856
   else if (alpha != ATL_rnone) goto ALPHAX;
 
857
   x0 = -X[0*incX];
 
858
   x1 = -X[1*incX];
 
859
   x2 = -X[2*incX];
 
860
   x3 = -X[3*incX];
 
861
   x4 = -X[4*incX];
 
862
   x5 = -X[5*incX];
 
863
   x6 = -X[6*incX];
 
864
   x7 = -X[7*incX];
 
865
   x8 = -X[8*incX];
 
866
   x9 = -X[9*incX];
 
867
   x10 = -X[10*incX];
 
868
   x11 = -X[11*incX];
 
869
   x12 = -X[12*incX];
 
870
   x13 = -X[13*incX];
 
871
   x14 = -X[14*incX];
 
872
X_IS_LOADED:
 
873
 
 
874
   for (j=0; j < N; j++, A += lda)
 
875
   {
 
876
      y0 = Y[incY*j];
 
877
      A[0] += x0 * y0;
 
878
      A[1] += x1 * y0;
 
879
      A[2] += x2 * y0;
 
880
      A[3] += x3 * y0;
 
881
      A[4] += x4 * y0;
 
882
      A[5] += x5 * y0;
 
883
      A[6] += x6 * y0;
 
884
      A[7] += x7 * y0;
 
885
      A[8] += x8 * y0;
 
886
      A[9] += x9 * y0;
 
887
      A[10] += x10 * y0;
 
888
      A[11] += x11 * y0;
 
889
      A[12] += x12 * y0;
 
890
      A[13] += x13 * y0;
 
891
      A[14] += x14 * y0;
 
892
   }
 
893
   return;
 
894
 
 
895
ALPHAX:
 
896
   y0 = alpha;
 
897
   x0 = X[0*incX] * y0;
 
898
   y0 = alpha;
 
899
   x1 = X[1*incX] * y0;
 
900
   y0 = alpha;
 
901
   x2 = X[2*incX] * y0;
 
902
   y0 = alpha;
 
903
   x3 = X[3*incX] * y0;
 
904
   y0 = alpha;
 
905
   x4 = X[4*incX] * y0;
 
906
   y0 = alpha;
 
907
   x5 = X[5*incX] * y0;
 
908
   y0 = alpha;
 
909
   x6 = X[6*incX] * y0;
 
910
   y0 = alpha;
 
911
   x7 = X[7*incX] * y0;
 
912
   y0 = alpha;
 
913
   x8 = X[8*incX] * y0;
 
914
   y0 = alpha;
 
915
   x9 = X[9*incX] * y0;
 
916
   y0 = alpha;
 
917
   x10 = X[10*incX] * y0;
 
918
   y0 = alpha;
 
919
   x11 = X[11*incX] * y0;
 
920
   y0 = alpha;
 
921
   x12 = X[12*incX] * y0;
 
922
   y0 = alpha;
 
923
   x13 = X[13*incX] * y0;
 
924
   y0 = alpha;
 
925
   x14 = X[14*incX] * y0;
 
926
   goto X_IS_LOADED;
 
927
ALPHA1:
 
928
   x0 = X[0*incX];
 
929
   x1 = X[1*incX];
 
930
   x2 = X[2*incX];
 
931
   x3 = X[3*incX];
 
932
   x4 = X[4*incX];
 
933
   x5 = X[5*incX];
 
934
   x6 = X[6*incX];
 
935
   x7 = X[7*incX];
 
936
   x8 = X[8*incX];
 
937
   x9 = X[9*incX];
 
938
   x10 = X[10*incX];
 
939
   x11 = X[11*incX];
 
940
   x12 = X[12*incX];
 
941
   x13 = X[13*incX];
 
942
   x14 = X[14*incX];
 
943
   goto X_IS_LOADED;
 
944
}
 
945
 
 
946
void Mjoin(PATL,gerk_Mlt16)
 
947
   (const int M, const int N, const SCALAR alpha, const TYPE *X, const int incX,
 
948
    const TYPE *Y, const int incY, TYPE *A, const int lda)
 
949
/*
 
950
 * ATLAS's normal GER kernels are optimized for long-M, and loop over rows in
 
951
 * the inner loop.  To avoid this killing us on short, wide matrices, have
 
952
 * special case code for M < 16.  This also allows kernels to assume M >= 16.
 
953
 */
 
954
{
 
955
   static gerk_t gerks[15]={ATL_gerk_Meq1, ATL_gerk_Meq2, ATL_gerk_Meq3,
 
956
                            ATL_gerk_Meq4, ATL_gerk_Meq5, ATL_gerk_Meq6,
 
957
                            ATL_gerk_Meq7, ATL_gerk_Meq8, ATL_gerk_Meq9,
 
958
                            ATL_gerk_Meq10, ATL_gerk_Meq11, ATL_gerk_Meq12,
 
959
                            ATL_gerk_Meq13, ATL_gerk_Meq14, ATL_gerk_Meq15};
 
960
   if (M < 1 || N < 1 || SCALAR_IS_ZERO(alpha))
 
961
      return;
 
962
   #ifdef ATL_GAS_x8664
 
963
   if (M > 14)
 
964
   #elif defined(ATL_GAS_x8632)
 
965
   if (M > 6)
 
966
   #else
 
967
   if (M > 15)
 
968
   #endif
 
969
   {
 
970
      Mjoin(PATL,gerk_axpy)(M, N, alpha, X, incX, Y, incY, A, lda);
 
971
      return;
 
972
   }
 
973
   gerks[M-1](M, N, alpha, X, incX, Y, incY, A, lda);
 
974
}
 
975
 
 
976
#else /* complex type */
 
977
#ifdef Conj_
 
978
void Mjoin(PATL,gerck_Mlt16)
 
979
#else
 
980
void Mjoin(PATL,gerk_Mlt16)
 
981
#endif
 
982
   (const int M, const int N, const SCALAR alpha, const TYPE *X, const int incX,
 
983
    const TYPE *Y, const int incY, TYPE *A, const int lda)
 
984
{
 
985
/*
 
986
 * For now, complex simply calls refblas for short M, axpy-based for large.
 
987
 * Probably not worth additional instruction load for complex to unroll.
 
988
 */
 
989
   #ifndef TUNING
 
990
   if (M < 8)
 
991
#ifdef Conj_
 
992
      Mjoin(PATL,refgerc)(M, N, alpha, X, incX, Y, incY, A, lda);
 
993
#else
 
994
      Mjoin(PATL,refgeru)(M, N, alpha, X, incX, Y, incY, A, lda);
 
995
#endif
 
996
   else
 
997
   #endif
 
998
#ifdef Conj_
 
999
      Mjoin(PATL,gerck_axpy)(M, N, alpha, X, incX, Y, incY, A, lda);
 
1000
#else
 
1001
      Mjoin(PATL,gerk_axpy)(M, N, alpha, X, incX, Y, incY, A, lda);
 
1002
#endif
 
1003
}
 
1004
#endif