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

« back to all changes in this revision

Viewing changes to tune/blas/level1/dottest.c

  • Committer: Bazaar Package Importer
  • Author(s): Sylvestre Ledru
  • Date: 2009-09-17 23:31:54 UTC
  • mto: (2.2.1 experimental)
  • mto: This revision was merged to the branch mainline in revision 10.
  • Revision ID: james.westby@ubuntu.com-20090917233154-9esw88ub02twbuab
Tags: upstream-3.8.3
ImportĀ upstreamĀ versionĀ 3.8.3

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
1
#include "atlas_misc.h"
2
2
#include <assert.h>
3
3
 
 
4
int FAx=0, MAx=0, FAy=0, MAy=0;
 
5
#include <stdio.h>
 
6
#include <stdlib.h>
 
7
#include <assert.h>
 
8
 
 
9
struct FA_allocs
 
10
{
 
11
   void *mem, *memA;
 
12
   struct FA_allocs *next;
 
13
} *allocQ=NULL;
 
14
 
 
15
struct FA_allocs *NewAlloc(size_t size, struct FA_allocs *next,
 
16
                           int align, int misalign)
 
17
/*
 
18
 * Allocates size allocation that is aligned to [align], but not aligned
 
19
 * to [misalign].  Therefore, misalign > align.  Align must minimally be sizeof
 
20
 * while misalign may be 0 if we don't need to avoid a particular alignment.
 
21
 */
 
22
{
 
23
   void *vp;
 
24
   char *cp;
 
25
   struct FA_allocs *ap;
 
26
   int n, i;
 
27
   const int malign = align >= misalign ? align : misalign;
 
28
 
 
29
   n = size + align + malign;
 
30
   i = (n >> 3)<<3;
 
31
   if (n != i)
 
32
      n += n - i;
 
33
   cp = malloc(n + sizeof(struct FA_allocs));
 
34
   assert(cp);
 
35
   ap = (struct FA_allocs *) (cp + n);
 
36
   ap->mem = cp;
 
37
/*
 
38
 * Align to min alignment
 
39
 */
 
40
   ap->memA = align ? (void*) ((((size_t) cp)/align)*align + align) : cp;
 
41
/*
 
42
 * Misalign to misalign
 
43
 */
 
44
   if (misalign)
 
45
   {
 
46
      if (((size_t)ap->memA)%misalign == 0)
 
47
         ap->memA = ((char*)ap->memA) + align;
 
48
   }
 
49
   ap->next = next;
 
50
   return(ap);
 
51
}
 
52
 
 
53
/*
 
54
 * no-align malloc free retaining system default behavior
 
55
 */
 
56
void *NA_malloc(size_t size)
 
57
{
 
58
   return(malloc(size));
 
59
}
 
60
void *NA_calloc(size_t n, size_t size)
 
61
{
 
62
   return(calloc(n, size));
 
63
}
 
64
void NA_free(void *ptr)
 
65
{
 
66
   free(ptr);
 
67
}
 
68
 
 
69
 
 
70
/*
 
71
 * malloc/free pair that aligns data to align, but not to misalign
 
72
 */
 
73
void *FA_malloc(size_t size, int align, int misalign)
 
74
{
 
75
   if ((!misalign && align <= 8) || !size)
 
76
      return(malloc(size));
 
77
   else
 
78
   {
 
79
      allocQ = NewAlloc(size, allocQ, align, misalign);
 
80
      return(allocQ->memA);
 
81
   }
 
82
}
 
83
void *FA_calloc(size_t n, size_t size, int align, int misalign)
 
84
{
 
85
   char *cp;
 
86
   int *ip;
 
87
   double *dp;
 
88
   size_t i;
 
89
   size_t tsize;
 
90
   tsize = n * size;
 
91
   cp = FA_malloc(tsize, align, misalign);
 
92
   if (size == sizeof(int))
 
93
      for (ip=(int*)cp,i=0; i < n; i++)
 
94
        ip[i] = 0;
 
95
   else if (size == sizeof(double))
 
96
      for (dp=(double*)cp,i=0; i < n; i++)
 
97
        dp[i] = 0.0;
 
98
   else
 
99
      for (i=0; i < tsize; i++)
 
100
        cp[i] = 0;
 
101
   return(cp);
 
102
}
 
103
 
 
104
void FA_free(void *ptr, int align, int misalign)
 
105
/*
 
106
 * Part of malloc/free pair that aligns data to FALIGN
 
107
 */
 
108
{
 
109
   struct FA_allocs *ap, *prev;
 
110
   if (ptr)
 
111
   {
 
112
      if ((!misalign && align <= 8))
 
113
         free(ptr);
 
114
      else
 
115
      {
 
116
         for (ap=allocQ; ap && ap->memA != ptr; ap = ap->next) prev = ap;
 
117
         if (!ap)
 
118
         {
 
119
            fprintf(stderr, "Couldn't find mem=%ld\nmemQ=\n", ptr);
 
120
            for (ap=allocQ; ap; ap = ap->next)
 
121
               fprintf(stderr, "   %ld, %ld\n", ap->memA, ap->mem);
 
122
         }
 
123
         assert(ap);
 
124
         if (ap == allocQ)
 
125
            allocQ = allocQ->next;
 
126
         else
 
127
            prev->next = ap->next;
 
128
         free(ap->mem);
 
129
      }
 
130
   }
 
131
}
 
132
 
4
133
#define dumb_seed(iseed_) srand(iseed_)
5
134
#ifndef RAND_MAX  /* rather dangerous non-ansi workaround */
6
135
   #define RAND_MAX ((unsigned long)(1<<30))
56
185
   for (i=0; i < N; i++) X[i] = alpha;
57
186
}
58
187
 
59
 
static TYPE *getvec(int npad, TYPE padval, int N, int incX)
 
188
static TYPE *getvec(int npad, TYPE padval, int N, int incX, int FA, int MA)
60
189
{
61
190
   TYPE *X, *x;
62
191
   int i, n;
64
193
   if (N <= 0) return(NULL);
65
194
   incX = Mabs(incX);
66
195
   n = 2*npad + 1+(N-1)*incX;
67
 
   X = malloc( ATL_sizeof*n );
 
196
   X = FA_malloc(ATL_sizeof*n, FA, MA);
68
197
   assert(X);
69
198
   vecset(n, padval, X);
70
199
   #ifdef TCPLX
95
224
   #endif
96
225
}
97
226
 
98
 
static TYPE *dupvec(int npad, int N, TYPE *X, int incX)
 
227
static TYPE *dupvec(int npad, int N, TYPE *X, int incX, int FA, int MA)
99
228
{
100
229
   int i, n;
101
230
   TYPE *y;
102
231
 
103
232
   incX = Mabs(incX);
104
233
   n = 1+(N-1)*incX + 2*npad;
105
 
   y = malloc(ATL_sizeof*n);
 
234
   y = FA_malloc(ATL_sizeof*n, FA, MA);
106
235
   assert(y);
107
236
   #ifdef TCPLX
108
237
      n *= 2;
112
241
}
113
242
 
114
243
static TYPE *gen_dupvec(int N, TYPE padval, int npadX, TYPE *X, int incX,
115
 
                        int npadY, int incY)
 
244
                        int npadY, int incY, int FA, int MA)
116
245
{
117
246
   int i, n;
118
247
   TYPE *y, *yy, *xx=X+(npadX SHIFT);
119
248
 
120
 
   y = getvec(npadY, padval, N, incY);
 
249
   y = getvec(npadY, padval, N, incY, FA, MA);
121
250
   yy = y + (npadY SHIFT);
122
251
   if (incY < 1) yy -= ((N-1)SHIFT) * incY;
123
252
   if (incX < 1) xx -= ((N-1)SHIFT) * incX;
139
268
   TYPE Mjoin(PATL, epsilon)(void);
140
269
   eps = Mjoin(PATL,epsilon)();
141
270
 
142
 
   x = X = getvec(0, padval, N, incX);
143
 
   y = Y = getvec(0, padval, N, incY);
 
271
   x = X = getvec(0, padval, N, incX, FAx, MAx);
 
272
   y = Y = getvec(0, padval, N, incY, FAy, MAy);
144
273
   if (incX < 1) x -= ((N-1)SHIFT) * incX;
145
274
   if (incY < 1) y -= ((N-1)SHIFT) * incY;
146
275
 
147
276
   #ifdef TREAL
 
277
      dotT = TEST_DOT(N, x, incX, y, incY);
148
278
      dotG = good_dot(N, x, incX, y, incY);
149
 
      dotT = TEST_DOT(N, x, incX, y, incY);
150
279
      diff = dotG - dotT;
151
280
      diff = Mabs(diff);
152
 
      if (diff > 2*N*eps) /* could be double this, but shouldn't in practice */
153
 
      {
 
281
      if (diff > 2*N*eps || dotT != dotT)
 
282
      {                     /* diff could 4*N*eps, but isn't in practice */
154
283
         fprintf(stderr,
155
284
                 "   DOT ERROR: N=%d, correct=%e, computed=%e, diff=%e!!\n",
156
285
                 N, dotG, dotT, diff);
179
308
      }
180
309
   #endif
181
310
 
182
 
   free(X);
183
 
   free(Y);
 
311
   FA_free(X, FAx, MAx);
 
312
   FA_free(Y, FAy, MAy);
184
313
   return(iret);
185
314
}
186
315
int DoAllTests(int nN, int *Ns, int nX, int *Xs, int nY, int *Ys)
224
353
void GetFlags(int nargs, char **args, int *nN, int **Ns, int *nX, int **incXs,
225
354
              int *nY, int **incYs)
226
355
{
227
 
   int i, j;
 
356
   int i, j, k, ig;
228
357
 
229
358
   *nY = *nX = *nN = -1;
230
359
 
234
363
      if (i == nargs-1) PrintUsage(args[0]);
235
364
      switch(args[i][1])
236
365
      {
 
366
      case 'F':
 
367
         j = args[i][2] != 'y';
 
368
         k = atoi(args[++i]);
 
369
         if (j)
 
370
         {
 
371
            if (k < 0)
 
372
               MAx = -k;
 
373
            else
 
374
               FAx = k;
 
375
         }
 
376
         else
 
377
         {
 
378
            if (k < 0)
 
379
               MAy = -k;
 
380
            else
 
381
               FAy = k;
 
382
         }
 
383
         break;
 
384
      case 'b':
 
385
      case 'a':
 
386
         ig = atoi(args[++i]);
 
387
         if (ig > nargs-i) PrintUsage(args[0]);
 
388
         i += ig SHIFT;
 
389
         break;
237
390
      case 'Y':
238
391
         *nY = atoi(args[++i]);
239
392
         if (*nY > nargs-i) PrintUsage(args[0]);
301
454
      assert(*incYs);
302
455
      **incYs = 1;
303
456
   }
 
457
   if (FAx < sizeof(TYPE))
 
458
      FAx = sizeof(TYPE);
 
459
   if (FAy < sizeof(TYPE))
 
460
      FAy = sizeof(TYPE);
304
461
}
305
462
 
306
463
main(int nargs, char **args)