▶ 第三章，逐步优化了一个二维卷积计算的过程
● 基准代码
1 #include <stdio.h>
2 #include <stdlib.h>
3 #include <string.h>
4 #include <math.h>
5 #include <time.h>
6 #include <sys/time.h>
7 #include <omp.h>
8 #include <assert.h>
9 #include <sys/mman.h>
10
11 #define REAL double
12 #define WIDTH 1030
13 #define HEIGHT 2048
14 #define COUNT 1000
15 #define PAD64 0
16 #if PAD64 // 是否对齐到 64 Byte
17 #define WIDTHP ((WIDTH * sizeof(REAL) + 63) / 64 * 64 / sizeof(REAL))
18 #else
19 #define WIDTHP WIDTH
20 #endif
21
22 void initbuf(REAL *fbuf, const int width, const int height) // 初始化矩阵
23 {
24 for (int y = 0; y < height; y++)
25 {
26 REAL val = (y % 10) ? 0 : 1.0;
27 for (int x = 0; x < width; x++)
28 fbuf[y * WIDTHP + x] = val;
29 }
30 return;
31 }
32
33 // 重复计算模糊 count 次
34 void stencil9pt(REAL *finp, REAL *foutp, const int width, const int height,
35 const REAL ctr, const REAL next, const REAL diag, const int count)
36 {
37 REAL *fin = finp, *fout = foutp;
38 for (int i = 0; i < count; i++)
39 {
40 for (int y = 1; y < height - 1; y++) // 不处理光环元素
41 {
42 int c = 1 + y * WIDTHP; // ？为什么要加两次 1
43 int n = c - WIDTHP, s = c + WIDTHP, w = c - 1, e = c + 1;
44 int nw = n - 1, ne = n + 1, sw = s - 1, se = s + 1;
45 for (int x = 1; x < width - 1; x++)
46 {
47 fout[c] = diag * fin[nw] + diag * fin[ne] + diag * fin[sw] + diag * fin[se] +
48 next * fin[w] + next * fin[e] + next * fin[n] + next * fin[s] + ctr * fin[c];
49 c++; n++; s++; e++; w++; nw++; ne++; sw++; se++;
50 }
51 }
52 REAL *ftmp = fin;
53 fin = fout;
54 fout = ftmp;
55 }
56 return;
57 }
58
59 static double dtime()
60 {
61 double tseconds = 0.0;
62 struct timeval mytime;
63 gettimeofday(&mytime, (struct timezone *) 0);
64 tseconds = (double)(mytime.tv_sec + (double)mytime.tv_usec * 1.0e-6);
65 return (tseconds);
66 }
67
68 int main(int argc, char *argv[])
69 {
70 REAL *fa = (REAL *)malloc(sizeof(REAL)*WIDTHP*HEIGHT), *fb = (REAL *)malloc(sizeof(REAL)*WIDTHP*HEIGHT);
71 assert(fa != MAP_FAILED);
72 assert(fb != MAP_FAILED);
73
74 printf("Initializing..%d Threads, %d x %d, PAD=%d..\n\n", omp_get_num_threads(), WIDTH, HEIGHT, WIDTHP);
75 initbuf(fa, WIDTHP, HEIGHT);
76 initbuf(fb, WIDTHP, HEIGHT);
77
78 printf("Running stencil kernel %d times\n", COUNT);
79 const REAL stendiag = 0.00125, stennext = 0.00125, stenctr = 0.99;
80 double time_b, time_e;
81 time_b = dtime();
82 stencil9pt(fa, fb, WIDTHP, HEIGHT, stenctr, stennext, stendiag, COUNT);
83 time_e = dtime();
84 printf("Elapsed time : %.3f (s)\n", time_e - time_b);
85 printf("FLOPS: %.3f (MFlops)\n", (WIDTHP * HEIGHT) * 17.0 * COUNT / (time_e - time_b) * 1.0e-06);// 计算一个元素需要 17 次乘法或加法
86
87 free(fa);
88 free(fb);
89 return 0;
90 }
1 #include <stdio.h>
