← Back to branch summary

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

« back to all changes in this revision

Viewing changes to tune/blas/gemv/mvktime.c

  • Committer: Package Import Robot
  • Author(s): Sébastien Villemot
  • Date: 2013-06-11 15:58:16 UTC
  • mfrom: (1.1.3 upstream)
  • mto: (2.2.21 experimental)
  • mto: This revision was merged to the branch mainline in revision 26.
  • Revision ID: package-import@ubuntu.com-20130611155816-b72z8f621tuhbzn0
Tags: upstream-3.10.1
Import upstream version 3.10.1

Show diffs side-by-side

added added

removed removed

Lines of Context:
tune/blas/gemv/mvktime.c
 
1
/*
 
2
 *             Automatically Tuned Linear Algebra Software v3.10.1
 
3
 * Copyright (C) 2010, 2009 R. Clint Whaley
 
4
 *
 
5
 * Redistribution and use in source and binary forms, with or without
 
6
 * modification, are permitted provided that the following conditions
 
7
 * are met:
 
8
 *   1. Redistributions of source code must retain the above copyright
 
9
 *      notice, this list of conditions and the following disclaimer.
 
10
 *   2. Redistributions in binary form must reproduce the above copyright
 
11
 *      notice, this list of conditions, and the following disclaimer in the
 
12
 *      documentation and/or other materials provided with the distribution.
 
13
 *   3. The name of the ATLAS group or the names of its contributers may
 
14
 *      not be used to endorse or promote products derived from this
 
15
 *      software without specific written permission.
 
16
 *
 
17
 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
 
18
 * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
 
19
 * TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
 
20
 * PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE ATLAS GROUP OR ITS CONTRIBUTORS
 
21
 * BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
 
22
 * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
 
23
 * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
 
24
 * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
 
25
 * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
 
26
 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
 
27
 * POSSIBILITY OF SUCH DAMAGE.
 
28
 *
 
29
 */
 
30
#include <stdio.h>
 
31
#include <stdlib.h>
 
32
#include <assert.h>
 
33
#include <string.h>
 
34
#include "atlas_misc.h"
 
35
#include Mstr(Mjoin(Mjoin(atlas_,PRE),sysinfo.h))
 
36
 
 
37
static char *resfile=NULL;
 
38
static FILE *fpres=NULL;
 
39
#include <stdio.h>
 
40
#include <stdlib.h>
 
41
#include <assert.h>
 
42
 
 
43
struct FA_allocs
 
44
{
 
45
   void *mem, *memA;
 
46
   struct FA_allocs *next;
 
47
} *allocQ=NULL;
 
48
 
 
49
struct FA_allocs *NewAlloc(size_t size, struct FA_allocs *next,
 
50
                           int align, int misalign)
 
51
/*
 
52
 * Allocates size allocation that is aligned to [align], but not aligned
 
53
 * to [misalign].  Therefore, misalign > align.  Align must minimally be sizeof
 
54
 * while misalign may be 0 if we don't need to avoid a particular alignment.
 
55
 */
 
56
{
 
57
   void *vp;
 
58
   char *cp;
 
59
   struct FA_allocs *ap;
 
60
   int n, i;
 
61
   const int malign = align >= misalign ? align : misalign;
 
62
 
 
63
   n = size + align + align + malign;
 
64
   i = (n >> 3)<<3;
 
65
   if (n != i)
 
66
      n += n - i;
 
67
   cp = malloc(n + sizeof(struct FA_allocs));
 
68
   assert(cp);
 
69
   ap = (struct FA_allocs *) (cp + n);
 
70
   ap->mem = cp;
 
71
/*
 
72
 * Align to min alignment
 
73
 */
 
74
   ap->memA = align ? (void*) ((((size_t) cp)/align)*align + align) : cp;
 
75
/*
 
76
 * Misalign to misalign
 
77
 * We often need to make sure to unaligned addresses share the same modulo
 
78
 * so that they have the *same* degree of misalignment (so that their alignment
 
79
 * can be fixed by simple peeling), and so in this case force the address
 
80
 * modulo the misalign to be the exact align value.
 
81
 */
 
82
   if (misalign)
 
83
      ap->memA = (void*)((((size_t)ap->memA)/malign)*malign + malign + align);
 
84
   ap->next = next;
 
85
   return(ap);
 
86
}
 
87
 
 
88
/*
 
89
 * no-align malloc free retaining system default behavior
 
90
 */
 
91
void *NA_malloc(size_t size)
 
92
{
 
93
   return(malloc(size));
 
94
}
 
95
void *NA_calloc(size_t n, size_t size)
 
96
{
 
97
   return(calloc(n, size));
 
98
}
 
99
void NA_free(void *ptr)
 
100
{
 
101
   free(ptr);
 
102
}
 
103
 
 
104
 
 
105
/*
 
106
 * malloc/free pair that aligns data to align, but not to misalign
 
107
 */
 
108
void *FA_malloc(size_t size, int align, int misalign)
 
109
{
 
110
   if ((!misalign && align <= 8) || !size)
 
111
      return(malloc(size));
 
112
   else
 
113
   {
 
114
      allocQ = NewAlloc(size, allocQ, align, misalign);
 
115
      return(allocQ->memA);
 
116
   }
 
117
}
 
118
void *FA_calloc(size_t n, size_t size, int align, int misalign)
 
119
{
 
120
   char *cp;
 
121
   int *ip;
 
122
   double *dp;
 
123
   size_t i;
 
124
   size_t tsize;
 
125
   tsize = n * size;
 
126
   cp = FA_malloc(tsize, align, misalign);
 
127
   if (size == sizeof(int))
 
128
      for (ip=(int*)cp,i=0; i < n; i++)
 
129
        ip[i] = 0;
 
130
   else if (size == sizeof(double))
 
131
      for (dp=(double*)cp,i=0; i < n; i++)
 
132
        dp[i] = 0.0;
 
133
   else
 
134
      for (i=0; i < tsize; i++)
 
135
        cp[i] = 0;
 
136
   return(cp);
 
137
}
 
138
 
 
139
void FA_free(void *ptr, int align, int misalign)
 
140
/*
 
141
 * Part of malloc/free pair that aligns data to FALIGN
 
142
 */
 
143
{
 
144
   struct FA_allocs *ap, *prev;
 
145
   if (ptr)
 
146
   {
 
147
      if ((!misalign && align <= 8))
 
148
         free(ptr);
 
149
      else
 
150
      {
 
151
         for (ap=allocQ; ap && ap->memA != ptr; ap = ap->next) prev = ap;
 
152
         if (!ap)
 
153
         {
 
154
            fprintf(stderr, "Couldn't find mem=%ld\nmemQ=\n", (size_t)ptr);
 
155
            for (ap=allocQ; ap; ap = ap->next)
 
156
               fprintf(stderr, "   %ld, %ld\n", (size_t)ap->memA,
 
157
                       (size_t)ap->mem);
 
158
         }
 
159
         assert(ap);
 
160
         if (ap == allocQ)
 
161
            allocQ = allocQ->next;
 
162
         else
 
163
            prev->next = ap->next;
 
164
         free(ap->mem);
 
165
      }
 
166
   }
 
167
}
 
168
#include "atlas_mvtesttime.h"
 
169
 
 
170
#ifdef TIME_KERNEL
 
171
   void ATL_UGEMV(ATL_CINT M, ATL_CINT N, const TYPE *A, ATL_CINT lda,
 
172
                  const TYPE *X, TYPE *Y);
 
173
#else
 
174
   #ifdef ATL_TRANS_
 
175
      void Mjoin(PATL,gemvT)
 
176
         (ATL_CINT M, ATL_CINT N, const SCALAR alpha, const TYPE *A,
 
177
          ATL_CINT lda, const TYPE *X, ATL_CINT incX,
 
178
          const SCALAR beta, TYPE *Y, ATL_CINT incY);
 
179
      #ifdef TREAL
 
180
         #define test_gemv(M, N, alpha, A, lda, X, incX, beta, Y, incY) \
 
181
            Mjoin(PATL,gemvT)(M, N, *alpha, A, lda, X, incX, *beta, Y, incY)
 
182
      #else
 
183
         #define test_gemv(M, N, alpha, A, lda, X, incX, beta, Y, incY) \
 
184
            Mjoin(PATL,gemvT)(M, N, alpha, A, lda, X, incX, beta, Y, incY)
 
185
      #endif
 
186
   #else
 
187
      void Mjoin(PATL,gemvN)
 
188
         (ATL_CINT M, ATL_CINT N, const SCALAR alpha, const TYPE *A,
 
189
          ATL_CINT lda, const TYPE *X, ATL_CINT incX,
 
190
          const SCALAR beta, TYPE *Y, ATL_CINT incY);
 
191
      #ifdef TREAL
 
192
         #define test_gemv(M, N, alpha, A, lda, X, incX, beta, Y, incY) \
 
193
            Mjoin(PATL,gemvN)(M, N, *alpha, A, lda, X, incX, *beta, Y, incY)
 
194
      #else
 
195
         #define test_gemv(M, N, alpha, A, lda, X, incX, beta, Y, incY) \
 
196
            Mjoin(PATL,gemvN)(M, N, alpha, A, lda, X, incX, beta, Y, incY)
 
197
      #endif
 
198
   #endif
 
199
#endif
 
200
 
 
201
double Time2Flop(ATL_INT M, ATL_INT N, double time)
 
202
{
 
203
   if (time == 0.0 || time != time)
 
204
      return(0.0);
 
205
   #ifdef TREAL
 
206
      return(((1.0e-6 * M)*(2.0*N+1.0))/time);
 
207
   #else
 
208
      return((((6.0*M)*(N+1.0) + (2.0*M)*N)*1.0e-6)/time);
 
209
   #endif
 
210
}
 
211
 
 
212
void Times2Flops(ATL_INT M, ATL_INT N, ATL_INT ntim, double *mf)
 
213
/*
 
214
 * Converts time to MFLOP
 
215
 */
 
216
{
 
217
   int i;
 
218
 
 
219
   for (i=0; i < ntim; i++)
 
220
      mf[i] = Time2Flop(M, N, mf[i]);
 
221
}
 
222
 
 
223
static double mysum(ATL_CINT N, double *d)
 
224
{
 
225
   int i;
 
226
   double sum;
 
227
 
 
228
   sum = d[0];
 
229
   for (i=1; i < N; i++)
 
230
      sum += d[i];
 
231
   return(sum);
 
232
}
 
233
 
 
234
#ifdef ATL_TRANS_
 
235
   #define NX M
 
236
   #define NY N
 
237
#else
 
238
   #define NX N
 
239
   #define NY M
 
240
#endif
 
241
#ifdef TIME_KERNEL
 
242
double mvtime(
 
243
   int verb,            /* verbosity */
 
244
   int nreps,           /* number of reps to do for one timing sample */
 
245
   size_t flushelts,    /* size of area to flush to avoid cache reuse */
 
246
   ATL_INT M,           /* # of rows of array A */
 
247
   ATL_INT N,           /* # of cols of array A */
 
248
   TYPE *alpha,
 
249
   TYPE *beta,
 
250
   ATL_INT lda,         /* leading dim */
 
251
   int incX,            /* ignored, assumed 1 by definition */
 
252
   int incY,            /* increment on Y; can be non-1 */
 
253
   int FAa,             /* if (FA. = 0) enforce no alignment */
 
254
   int MAa,             /* else force op to be aligned to at least FA bytes */
 
255
   int FAx,             /* if MA. != 0, disallow op to be aligned to MA. byts */
 
256
   int MAx,
 
257
   int FAy,
 
258
   int MAy)
 
259
/*
 
260
 * This function directly calls the kernel routine on data that has been
 
261
 * preloaded (through initialization) to any cache large enough to hold it.
 
262
 */
 
263
{
 
264
   double t0, t1;
 
265
   TYPE *A, *X, *Y;
 
266
   ATL_INT i;
 
267
   int k;
 
268
 
 
269
   A = FA_malloc(ATL_MulBySize(lda)*N, FAa, MAa);
 
270
   X = FA_malloc(ATL_MulBySize(NX), FAx, MAx);
 
271
   Y = FA_malloc(ATL_MulBySize(NY*incY), FAy, MAy);
 
272
   ATL_assert(A && X && Y);
 
273
 
 
274
   Mjoin(PATL,gegen)(1, NY, Y, incY, M);
 
275
   Mjoin(PATL,gegen)(NX, 1, X, NX, N+127*50+77);
 
276
   Mjoin(PATL,gegen)(M, N, A, lda, N*M+513*7+90);
 
277
/*
 
278
 * NOTE: if nreps too high this could lead to under/overflow
 
279
 */
 
280
   for (k=0; k < 8; k++)  /* loop until results are believable, or give up */
 
281
   {
 
282
      t0 = time00();
 
283
      for (i=nreps; i; i--)
 
284
      {
 
285
         ATL_UGEMV(M, N, A, lda, X, Y);
 
286
      }
 
287
      t1 = time00();
 
288
      if (t1 > t0) break;
 
289
      nreps = (nreps) ? nreps+nreps : 1;
 
290
   }
 
291
   t1 = (t1 - t0)/(1.0*nreps);
 
292
   if (verb)
 
293
      fprintf(stdout, "   M=%d, N=%d, lda=%d, nreps=%d, time=%e, mflop=%.2f\n",
 
294
              M, N, lda, nreps, t1, Time2Flop(M, N, t1));
 
295
   FA_free(A, FAa, MAa);
 
296
   FA_free(X, FAx, MAx);
 
297
   FA_free(Y, FAy, MAy);
 
298
   return(t1);
 
299
}
 
300
#else
 
301
double mvtime(
 
302
   int verb,            /* verbosity */
 
303
   int nreps,           /* number of reps to do for one timing sample */
 
304
   ATL_INT flushelts,   /* size of area to flush to avoid cache reuse */
 
305
   ATL_INT M,           /* # of rows of array A */
 
306
   ATL_INT N,           /* # of cols of array A */
 
307
   TYPE *alpha,
 
308
   TYPE *beta,
 
309
   ATL_INT lda,         /* leading dim */
 
310
   int mu,              /* unrolling on M */
 
311
   int nu,              /* unrolling on N */
 
312
   int FAa,             /* if (FA. = 0) enforce no alignment */
 
313
   int MAa,             /* else force op to be aligned to at least FA bytes */
 
314
   int FAx,             /* if MA. != 0, disallow op to be aligned to MA. bytes*/
 
315
   int MAx,
 
316
   int FAy,
 
317
   int MAy)
 
318
/*
 
319
 * Times the kernel for out-of-cache (where flushelts sets the cache that it
 
320
 * is not allowed to be in) use.
 
321
 * RETURNS: elapsed time in seconds to average repitition of indicated problem.
 
322
 */
 
323
{
 
324
   #ifdef TREAL
 
325
      TYPE NONE = -1.0;
 
326
   TYPE nbeta[1] = {-beta[0]}, *be=beta;
 
327
   #else
 
328
      TYPE NONE[2] = {-1.0, 0.0};
 
329
   TYPE nbeta[2] = {-beta[0], -beta[1]}, *be=beta;
 
330
   #endif
 
331
   double t0, t1;
 
332
   TYPE *A, *X, *Y, *a, *x, *y;
 
333
   void *vmem;
 
334
   size_t Aelts, Xelts, Yelts, setspan, ygap, xgap, agap, pregap, setsz, nsets;
 
335
   size_t ptr_st;
 
336
   ATL_INT i, j;
 
337
   int k, maxalign;
 
338
 
 
339
   if (MAx)
 
340
      assert(MAx != FAx);
 
341
   if (MAy)
 
342
      assert(MAy != FAy);
 
343
   if (MAa)
 
344
      assert(MAa != FAa);
 
345
/*
 
346
 * Find basic length of each operand in elements
 
347
 */
 
348
   Aelts = lda * N;
 
349
   #ifdef ATL_TRANS_
 
350
      Xelts = M;
 
351
      Yelts = N;
 
352
   #else
 
353
      Xelts = N;
 
354
      Yelts = M;
 
355
   #endif
 
356
/*
 
357
 * Map memory so that we can enforce all required alignments while moving
 
358
 * through memory; mem starts with maxalign-aligned memory, so that we can
 
359
 * guarantee all further alignments
 
360
 */
 
361
   maxalign = (MAx) ? MAx : 1;
 
362
   if (MAy)
 
363
      maxalign = ATL_lcm(MAy,maxalign);
 
364
   if (MAa)
 
365
      maxalign = ATL_lcm(MAa,maxalign);
 
366
   if (FAx)
 
367
      maxalign = ATL_lcm(FAx,maxalign);
 
368
   if (FAy)
 
369
      maxalign = ATL_lcm(FAy,maxalign);
 
370
   if (FAa)
 
371
      maxalign = ATL_lcm(FAa,maxalign);
 
372
   if (maxalign == 1)
 
373
      maxalign = 0;
 
374
   j = (FAx) ? FAx : sizeof(TYPE);
 
375
   if (MAx)
 
376
      for (i=0; (i % j != 0 || i%MAx == 0); i += sizeof(TYPE));
 
377
   else if (FAx)
 
378
      for (i=0; i % j != 0 ; i += sizeof(TYPE));
 
379
   else
 
380
      i = 0;
 
381
   pregap = i;
 
382
   xgap = ATL_MulBySize(Xelts);
 
383
   if (FAy || MAy)
 
384
   {
 
385
      j = (FAy) ? FAy : sizeof(TYPE);
 
386
      if (MAy)
 
387
         for (i=pregap+xgap; (i%j != 0 || i%MAy == 0); i += sizeof(TYPE));
 
388
      else
 
389
         for (i=pregap+xgap; (i%j != 0); i += sizeof(TYPE));
 
390
      xgap = i - pregap;
 
391
   }
 
392
   ygap = ATL_MulBySize(Yelts);
 
393
   if (FAa || MAa)
 
394
   {
 
395
      j = (FAa) ? FAa : sizeof(TYPE);
 
396
      if (MAa)
 
397
         for (i=pregap+xgap+ygap; (i%j != 0 || i%MAa == 0); i += sizeof(TYPE));
 
398
      else
 
399
         for (i=pregap+xgap+ygap; (i%j != 0); i += sizeof(TYPE));
 
400
      ygap = i - pregap - xgap;
 
401
   }
 
402
   agap = ATL_MulBySize(Aelts);
 
403
 
 
404
   if (maxalign)
 
405
   {
 
406
      j = pregap;
 
407
      for (i=pregap+xgap+ygap+agap; i%maxalign != 0; i++);
 
408
      agap = i-pregap-xgap-ygap;
 
409
   }
 
410
   setspan = pregap + xgap + ygap + agap;
 
411
   assert(setspan%sizeof(TYPE) == 0);
 
412
   setsz = ATL_MulBySize(M+N+M*N);
 
413
   nsets = (ATL_MulBySize(flushelts)+setsz-1)/setsz;
 
414
   if (!nsets)
 
415
      nsets = 1;
 
416
   vmem = malloc(maxalign + nsets*setspan);
 
417
   assert(vmem);
 
418
   if (maxalign)   /* start maxaligned to guarantee all alignments */
 
419
      for (ptr_st = (size_t)vmem; ptr_st%maxalign; ptr_st++);
 
420
   else ptr_st = (size_t) vmem;
 
421
   X = (TYPE*) (ptr_st + pregap);
 
422
   Y = (TYPE*) (ptr_st + pregap + xgap);
 
423
   A = (TYPE*) (ptr_st + pregap + xgap + ygap);
 
424
/*
 
425
 * Set ptrs to last set in memory
 
426
 */
 
427
   setspan /= sizeof(TYPE);
 
428
   a = A += (nsets-1) * setspan;
 
429
   x = X += (nsets-1) * setspan;
 
430
   y = Y += (nsets-1) * setspan;
 
431
   for (i=nsets; i; i--)
 
432
   {
 
433
      #define DEBUG_FA
 
434
      #ifdef DEBUG_FA
 
435
         if (FAa)
 
436
            assert(((size_t)a)%FAa == 0);
 
437
         if (FAx)
 
438
            assert(((size_t)x)%FAx == 0);
 
439
         if (FAy)
 
440
            assert(((size_t)y)%FAy == 0);
 
441
         if (MAa)
 
442
            assert(((size_t)a)%MAa != 0);
 
443
         if (MAx)
 
444
            assert(((size_t)x)%MAx != 0);
 
445
         if (MAy)
 
446
            assert(((size_t)y)%MAy != 0);
 
447
      #endif
 
448
      Mjoin(PATL,gegen)(Yelts, 1, y, Yelts, M);
 
449
      Mjoin(PATL,gegen)(Xelts, 1, x, Xelts, N+127*50+77);
 
450
      if (i&1)
 
451
         Mjoin(PATL,scal)(Xelts, NONE, x, 1);
 
452
      Mjoin(PATL,gegen)(M, N, A, lda, N*M+513*7+90);
 
453
      a -= setspan; x -= setspan; y -= setspan;
 
454
   }
 
455
   a = A; x = X; y = Y;
 
456
 
 
457
   j=0;
 
458
   for (k=0; k < 8; k++) /* loop until good timing or too many trips */
 
459
   {
 
460
      t0 = time00();
 
461
      for (i=nreps; i; i--)
 
462
      {
 
463
         test_gemv(M, N, alpha, A, lda, x, 1, be, y, 1);
 
464
         if (++j < nsets) { a -= setspan; x -= setspan; y -= setspan; }
 
465
         else  { a = A; x = X; y = Y; j=0; be = (be == beta) ? nbeta : beta;}
 
466
      }
 
467
      t1 = time00();
 
468
      if (t1 > t0)
 
469
         break;
 
470
      nreps = (nreps) ? nreps+nreps : 1;
 
471
   }
 
472
   free(vmem);
 
473
   t1 = (t1-t0) / (1.0*nreps);
 
474
   if (verb)
 
475
      fprintf(stdout, "   M=%d, N=%d, lda=%d, nreps=%d, time=%e, mflop=%.2f\n",
 
476
              M, N, lda, nreps, t1, Time2Flop(M, N, t1));
 
477
   return(t1);
 
478
}
 
479
#endif
 
480
 
 
481
void DoTimes(int verb, size_t flshelts, ATL_INT ntim, ATL_INT nrep,
 
482
             ATL_INT mu, ATL_INT nu, ATL_INT M, ATL_INT N, TYPE *alpha,
 
483
             TYPE *beta, ATL_INT lda,
 
484
             int FAa, int MAa, int FAx, int MAx, int FAy, int MAy)
 
485
{
 
486
   double *times;
 
487
   int i;
 
488
 
 
489
   times = malloc(ntim * sizeof(double));
 
490
   assert(times);
 
491
 
 
492
#ifdef TREAL
 
493
   fprintf(stdout, "GEMV: M=%d, N=%d, lda=%d, AF=[%d,%d,%d], AM=[%d,%d,%d], beta=%e, alpha=%e:\n",
 
494
           M, N, lda, FAa, FAx, FAy, MAa, MAx, MAy, *beta, *alpha);
 
495
#else
 
496
   fprintf(stdout, "GEMV: M=%d, N=%d, lda=%d, AF=[%d,%d,%d], AM=[%d,%d,%d], beta=[%e,%e], alpha=[%e,%e]:\n",
 
497
           M, N, lda, FAa, FAx, FAy, MAa, MAx, MAy,
 
498
           *beta, beta[1], *alpha, alpha[1]);
 
499
#endif
 
500
   for (i=0; i < ntim; i++)
 
501
      times[i] = mvtime(verb, nrep, flshelts, M, N, alpha, beta, lda,
 
502
                        mu, nu, FAa, MAa, FAx, MAx, FAy, MAy);
 
503
   SortDoubles(ntim, times);
 
504
   Times2Flops(M, N, ntim, times);
 
505
   if (fpres)
 
506
   {
 
507
      #if defined(PentiumCPS) || defined(WALL)
 
508
         fprintf(fpres, "%d 1\n", ntim);
 
509
      #else
 
510
         fprintf(fpres, "%d 0\n", ntim);
 
511
      #endif
 
512
      for (i=0; i < ntim; i++)
 
513
         fprintf(fpres, "%le\n", times[i]);
 
514
      fclose(fpres);
 
515
   }
 
516
   fprintf(stdout, "NREPS=%d, MAX=%.2f, MIN=%.2f, AVG=%.2f, MED=%.2f\n",
 
517
           ntim, times[0], times[ntim-1], mysum(ntim, times)/ntim,
 
518
           times[ntim>>1]);
 
519
   free(times);
 
520
}
 
521
 
 
522
void PrintUsage(char *name, char *arg, int i)
 
523
{
 
524
   if (i > 0)
 
525
      fprintf(stderr, "BAD ARG '%s' on %dth FLAG\n", arg, i);
 
526
   fprintf(stderr, "USAGE: %s [flags], where flags are:\n", name);
 
527
   fprintf(stderr, "   -v <#> : set verbosity level\n");
 
528
   fprintf(stderr, "   -C <#> : set flushsz = # (kbytes)\n");
 
529
   fprintf(stderr, "   -x <#> : unrolling for X in kernel is #\n");
 
530
   fprintf(stderr, "   -y <#> : unrolling for Y in kernel is #\n");
 
531
   fprintf(stderr, "   -m <#> : set # of rows of matrix to #\n");
 
532
   fprintf(stderr, "   -n <#> : set # of cols of matrix to #\n");
 
533
   fprintf(stderr, "   -l <#> : set leading dimension of array A to #\n");
 
534
   fprintf(stderr, "   -F <#> : do at least # MFLOPS for each timing interval\n");
 
535
   fprintf(stderr, "   -f <file> : output timing summary in <file>; if file exists read & report\n");
 
536
   fprintf(stderr,
 
537
           "   -r <#> : do # repetitions of the call for each timing interval\n");
 
538
   fprintf(stderr,
 
539
      "   -# <#> : report # timings (each interval may have multiple calls)\n");
 
540
   fprintf(stderr,
 
541
"   -F[x,y,a] <#> : if(# > 0) -> force op to be aligned to at least # bytes\n");
 
542
   fprintf(stderr,
 
543
"                   if(# < 0) -> force op to be aligned to < # bytes.\n");
 
544
   fprintf(stderr, "   -b <beta> : 2 floats for complex, one for real.\n");
 
545
   exit(i ? i : -1);
 
546
}
 
547
 
 
548
void GetFlags(int nargs, char **args, int *verb,
 
549
              size_t *flushelts, ATL_INT *celts, ATL_INT *pgelts,
 
550
              ATL_INT *mu, ATL_INT *nu, ATL_INT *ntim, ATL_INT *nrep,
 
551
              enum ATLAS_TRANS *TA, ATL_INT *m, ATL_INT *n, ATL_INT *lda,
 
552
              TYPE *beta,
 
553
              int *FAa, int *MAa, int *FAx, int *MAx, int *FAy, int *MAy)
 
554
{
 
555
   double mfF=ATL_nkflop/1000.0, flops;
 
556
   ATL_INT j, h;
 
557
   size_t il;
 
558
   int i;
 
559
   char ch;
 
560
 
 
561
   *verb = 1;
 
562
   #ifdef ATL_PAGESZ
 
563
      *pgelts = ATL_DivBySize(ATL_PAGESZ);
 
564
   #else
 
565
      *pgelts = 4*ATL_DivBySize(1024);
 
566
   #endif
 
567
   *celts = 0.75*ATL_L1elts;
 
568
   #ifdef L2SIZE
 
569
      *flushelts = L2SIZE;
 
570
   #else
 
571
      *flushelts = 8*1024*ATL_DivBySize(1024);
 
572
   #endif
 
573
   *mu = *nu = 1;
 
574
   *m = 800;
 
575
   *n = 200;
 
576
   *nrep = *lda = 0;
 
577
   *ntim = 3;
 
578
   *FAa = *MAa = *FAx = *MAx = *FAy = *MAy = 0;
 
579
   *beta = 1.0;
 
580
   #ifdef TCPLX
 
581
      beta[1] = 0.0;
 
582
   #endif
 
583
 
 
584
   for (i=1; i < nargs; i++)
 
585
   {
 
586
      if (args[i][0] != '-')
 
587
         PrintUsage(args[0], "No '-' preceeding flag!", i);
 
588
      switch(args[i][1])
 
589
      {
 
590
      case 'f' :  /* set resfile output */
 
591
         if (++i >= nargs)
 
592
            PrintUsage(args[0], "out of flags in -f ", i-1);
 
593
         resfile = args[i];
 
594
         break;
 
595
      case 'v' :  /* set verbosity */
 
596
         if (++i >= nargs)
 
597
            PrintUsage(args[0], "out of flags in -v ", i-1);
 
598
         *verb = atoi(args[i]);
 
599
         break;
 
600
      case 'G' :  /* set GEMV blocking cache size in KB */
 
601
         if (++i >= nargs)
 
602
            PrintUsage(args[0], "out of flags in -G ", i-1);
 
603
         j = atoi(args[i]);
 
604
         *celts = j*ATL_DivBySize(1024);
 
605
         break;
 
606
      case 'A' :  /* set transpose */
 
607
         if (++i >= nargs)
 
608
            PrintUsage(args[0], "out of flags in -A ", i-1);
 
609
         ch = args[i][0];
 
610
         if (ch == 't' || ch == 'T')
 
611
            *TA = AtlasTrans;
 
612
         else if (ch == 'c' || ch == 'C')
 
613
            *TA = AtlasConjTrans;
 
614
         else if (ch == 'z' || ch == 'Z')
 
615
            *TA = AtlasConj;
 
616
         else
 
617
            *TA = AtlasNoTrans;
 
618
         break;
 
619
      case 'C' :  /* set flushsz in KB */
 
620
         if (++i >= nargs)
 
621
            PrintUsage(args[0], "out of flags in -C ", i-1);
 
622
         il = atol(args[i]);
 
623
         if (il >= 0)
 
624
            *flushelts = il*ATL_DivBySize(1024);
 
625
         break;
 
626
      case 'p' :  /* set pagesz in KB */
 
627
         if (++i >= nargs)
 
628
            PrintUsage(args[0], "out of flags in -p ", i-1);
 
629
         j = atoi(args[i]);
 
630
         *pgelts = j*ATL_DivBySize(1024);
 
631
         break;
 
632
      case 'x' :  /* set mu */
 
633
         if (++i >= nargs)
 
634
            PrintUsage(args[0], "out of flags in -x ", i-1);
 
635
         *mu = atoi(args[i]);
 
636
         break;
 
637
      case 'y' :  /* set nu */
 
638
         if (++i >= nargs)
 
639
            PrintUsage(args[0], "out of flags in -y ", i-1);
 
640
         *nu = atoi(args[i]);
 
641
         break;
 
642
      case 'm' :  /* set M */
 
643
         if (++i >= nargs)
 
644
            PrintUsage(args[0], "out of flags in -m ", i-1);
 
645
         *m = atoi(args[i]);
 
646
         break;
 
647
      case 'n' :  /* set N */
 
648
         if (++i >= nargs)
 
649
            PrintUsage(args[0], "out of flags in -n ", i-1);
 
650
         *n = atoi(args[i]);
 
651
         break;
 
652
      case 'l' :  /* set lda */
 
653
         if (++i >= nargs)
 
654
            PrintUsage(args[0], "out of flags in -l ", i-1);
 
655
         *lda = atoi(args[i]);
 
656
         break;
 
657
      case 'a' : /* alias for setting alpha in r1ktime */
 
658
      case 'b' : /* set beta */
 
659
         if (++i >= nargs)
 
660
            PrintUsage(args[0], "out of flags in -b ", i-1);
 
661
         *beta = atof(args[i]);
 
662
         #ifdef TCPLX
 
663
            if (++i >= nargs)
 
664
               PrintUsage(args[0], "out of flags in -b ", i-1);
 
665
            beta[1] = atof(args[i]);
 
666
         #endif
 
667
         break;
 
668
      case 'F' :  /* set nrep by specifying MFLOPS, or force alignment */
 
669
         ch = args[i][2];
 
670
         if (ch == '\0')   /* specifying MFLOPS */
 
671
         {
 
672
            if (++i >= nargs)
 
673
               PrintUsage(args[0], "out of flags in -F ", i-1);
 
674
            j = atoi(args[i]);
 
675
            mfF = j;
 
676
         }
 
677
         else
 
678
         {
 
679
            if (ch != 'a' && ch != 'y' && ch != 'x')
 
680
               PrintUsage(args[0], args[i], i);
 
681
            if (++i >= nargs)
 
682
               PrintUsage(args[0], args[i-1], i-1);
 
683
            j = atoi(args[i]);
 
684
            if (j < 0)
 
685
            {
 
686
               if (ch == 'a')
 
687
                  *MAa = -j;
 
688
               else if (ch == 'y')
 
689
                  *MAy = -j;
 
690
               else if (ch == 'x')
 
691
                  *MAx = -j;
 
692
            }
 
693
            else
 
694
            {
 
695
               if (ch == 'a')
 
696
                  *FAa = j;
 
697
               else if (ch == 'y')
 
698
                  *FAy = j;
 
699
               else if (ch == 'x')
 
700
                  *FAx = j;
 
701
            }
 
702
         }
 
703
         break;
 
704
      case 'r' :  /* set nrep directly as integer */
 
705
         if (++i >= nargs)
 
706
            PrintUsage(args[0], "out of flags in -r ", i-1);
 
707
         *nrep = atoi(args[i]);
 
708
         break;
 
709
      case '#' :  /* set number of timings to report */
 
710
         if (++i >= nargs)
 
711
            PrintUsage(args[0], "out of flags in -# ", i-1);
 
712
         *ntim = atoi(args[i]);
 
713
         break;
 
714
      default:
 
715
         PrintUsage(args[0], args[i], i);
 
716
      }
 
717
   }
 
718
   if (!(*nrep))
 
719
   {
 
720
      flops = Time2Flop(*m, *n, 1.0) * 1000.0;  /* Get kiloFLOPS in GEMV */
 
721
      *nrep = (mfF+flops-1)/flops;
 
722
      if (*nrep < 1) *nrep = 1;
 
723
   }
 
724
   if (!(*lda))
 
725
      *lda = *m + 8;
 
726
}
 
727
int main(int nargs, char **args)
 
728
{
 
729
   size_t flushelts;
 
730
   ATL_INT celts, pgelts, mu, nu, ntim, nrep, m, n, lda;
 
731
   int FAa, MAa, FAx, MAx, FAy, MAy;    /* Force & Max align for ops */
 
732
   int verb;
 
733
   enum ATLAS_TRANS TA;
 
734
   double *dres;
 
735
   #ifdef TREAL
 
736
      TYPE beta;
 
737
   #else
 
738
      TYPE beta[2];
 
739
   #endif
 
740
 
 
741
   GetFlags(nargs, args, &verb, &flushelts, &celts, &pgelts, &mu, &nu, &ntim,
 
742
            &nrep, &TA, &m, &n, &lda, SADD beta,
 
743
            &FAa, &MAa, &FAx, &MAx, &FAy, &MAy);
 
744
   if (resfile)
 
745
   {
 
746
      dres = ReadResultsFile(1, ntim, resfile);
 
747
      if (dres)
 
748
      {
 
749
         fprintf(stdout, "TIMINGS READ IN FROM '%s':\n", resfile);
 
750
         PrintResultsFromFile(stdout, dres);
 
751
         free(dres);
 
752
         exit(0);
 
753
      }
 
754
      fpres = fopen(resfile, "w");
 
755
      assert(fpres);
 
756
   }
 
757
   DoTimes(verb, flushelts, ntim, nrep, mu, nu, m, n, SADD beta, SADD beta, lda,
 
758
           FAa, MAa, FAx, MAx, FAy, MAy);
 
759
   exit(0);
 
760
}