Performance/bin/kernel.c

0001 #include <stdio.h>
0002 #include <stdlib.h>
0003 #include "LU.h"
0004 #include "FFT.h"
0005 #include "SOR.h"
0006 #include "MonteCarlo.h"
0007 #include "LU.h"
0008 #include "Random.h"
0009 #include "Stopwatch.h"
0010 #include "SparseCompRow.h"
0011 #include "array.h"
0012
0013 double kernel_measureFFT(int N, double mintime, Random R) {
0014   /* initialize FFT data as complex (N real/img pairs) */
0015
0016   int twoN = 2 * N;
0017   double *x = RandomVector(twoN, R);
0018   long cycles = 1;
0019   Stopwatch Q = new_Stopwatch();
0020   int i = 0;
0021   double result = 0.0;
0022
0023   while (1) {
0024     Stopwatch_start(Q);
0025     for (i = 0; i < cycles; i++) {
0026       FFT_transform(twoN, x); /* forward transform */
0027       FFT_inverse(twoN, x);   /* backward transform */
0028     }
0029     Stopwatch_stop(Q);
0030     if (Stopwatch_read(Q) >= mintime)
0031       break;
0032
0033     cycles *= 2;
0034   }
0035   /* approx Mflops */
0036
0037   result = FFT_num_flops(N) * cycles / Stopwatch_read(Q) * 1.0e-6;
0038   Stopwatch_delete(Q);
0039   free(x);
0040   return result;
0041 }
0042
0043 double kernel_measureSOR(int N, double min_time, Random R) {
0044   double **G = RandomMatrix(N, N, R);
0045   double result = 0.0;
0046
0047   Stopwatch Q = new_Stopwatch();
0048   int cycles = 1;
0049   while (1) {
0050     Stopwatch_start(Q);
0051     SOR_execute(N, N, 1.25, G, cycles);
0052     Stopwatch_stop(Q);
0053
0054     if (Stopwatch_read(Q) >= min_time)
0055       break;
0056
0057     cycles *= 2;
0058   }
0059   /* approx Mflops */
0060
0061   result = SOR_num_flops(N, N, cycles) / Stopwatch_read(Q) * 1.0e-6;
0062   Stopwatch_delete(Q);
0063   Array2D_double_delete(N, N, G);
0064   return result;
0065 }
0066
0067 double kernel_measureMonteCarlo(double min_time, Random R) {
0068   double result = 0.0;
0069   Stopwatch Q = new_Stopwatch();
0070
0071   int cycles = 1;
0072   while (1) {
0073     Stopwatch_start(Q);
0074     MonteCarlo_integrate(cycles);
0075     Stopwatch_stop(Q);
0076     if (Stopwatch_read(Q) >= min_time)
0077       break;
0078
0079     cycles *= 2;
0080   }
0081   /* approx Mflops */
0082   result = MonteCarlo_num_flops(cycles) / Stopwatch_read(Q) * 1.0e-6;
0083   Stopwatch_delete(Q);
0084   return result;
0085 }
0086
0087 double kernel_measureSparseMatMult(int N, int nz, double min_time, Random R) {
0088   /* initialize vector multipliers and storage for result */
0089   /* y = A*y;  */
0090
0091   double *x = RandomVector(N, R);
0092   double *y = (double *)malloc(sizeof(double) * N);
0093
0094   double result = 0.0;
0095
0096 #if 0
0097         // initialize square sparse matrix
0098         //
0099         // for this test, we create a sparse matrix with M/nz nonzeros
0100         // per row, with spaced-out evenly between the begining of the
0101         // row to the main diagonal.  Thus, the resulting pattern looks
0102         // like
0103         //             +-----------------+
0104         //             +*                +
0105         //             +***              +
0106         //             +* * *            +
0107         //             +** *  *          +
0108         //             +**  *   *        +
0109         //             +* *   *   *      +
0110         //             +*  *   *    *    +
0111         //             +*   *    *    *  +
0112         //             +-----------------+
0113         //
0114         // (as best reproducible with integer artihmetic)
0115         // Note that the first nr rows will have elements past
0116         // the diagonal.
0117 #endif
0118
0119   int nr = nz / N;  /* average number of nonzeros per row  */
0120   int anz = nr * N; /* _actual_ number of nonzeros         */
0121
0122   double *val = RandomVector(anz, R);
0123   int *col = (int *)malloc(sizeof(int) * nz);
0124   int *row = (int *)malloc(sizeof(int) * (N + 1));
0125   int r = 0;
0126   int cycles = 1;
0127
0128   Stopwatch Q = new_Stopwatch();
0129
0130   row[0] = 0;
0131   for (r = 0; r < N; r++) {
0132     /* initialize elements for row r */
0133
0134     int rowr = row[r];
0135     int step = r / nr;
0136     int i = 0;
0137
0138     row[r + 1] = rowr + nr;
0139     if (step < 1)
0140       step = 1; /* take at least unit steps */
0141
0142     for (i = 0; i < nr; i++)
0143       col[rowr + i] = i * step;
0144   }
0145
0146   while (1) {
0147     Stopwatch_start(Q);
0148     SparseCompRow_matmult(N, y, val, row, col, x, cycles);
0149     Stopwatch_stop(Q);
0150     if (Stopwatch_read(Q) >= min_time)
0151       break;
0152
0153     cycles *= 2;
0154   }
0155   /* approx Mflops */
0156   result = SparseCompRow_num_flops(N, nz, cycles) / Stopwatch_read(Q) * 1.0e-6;
0157
0158   Stopwatch_delete(Q);
0159   free(row);
0160   free(col);
0161   free(val);
0162   free(y);
0163   free(x);
0164
0165   return result;
0166 }
0167
0168 double kernel_measureLU(int N, double min_time, Random R) {
0169   double **A = NULL;
0170   double **lu = NULL;
0171   int *pivot = NULL;
0172
0173   Stopwatch Q = new_Stopwatch();
0174   double result = 0.0;
0175   int i = 0;
0176   int cycles = 1;
0177
0178   if ((A = RandomMatrix(N, N, R)) == NULL)
0179     exit(1);
0180   if ((lu = new_Array2D_double(N, N)) == NULL)
0181     exit(1);
0182   if ((pivot = (int *)malloc(N * sizeof(int))) == NULL)
0183     exit(1);
0184
0185   while (1) {
0186     Stopwatch_start(Q);
0187     for (i = 0; i < cycles; i++) {
0188       Array2D_double_copy(N, N, lu, A);
0189       LU_factor(N, N, lu, pivot);
0190     }
0191     Stopwatch_stop(Q);
0192     if (Stopwatch_read(Q) >= min_time)
0193       break;
0194
0195     cycles *= 2;
0196   }
0197   /* approx Mflops */
0198   result = LU_num_flops(N) * cycles / Stopwatch_read(Q) * 1.0e-6;
0199
0200   Stopwatch_delete(Q);
0201   free(pivot);
0202   Array2D_double_delete(N, N, lu);
0203   Array2D_double_delete(N, N, A);
0204
0205   return result;
0206 }