1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38 package ffx.numerics.math;
39
40 import static org.apache.commons.math3.util.FastMath.floor;
41
42 import java.util.Random;
43
44
45
46
47
48
49
50 public class SSETest {
51
52 public final double[][] A;
53 public final double[] x;
54 private final double[] flatA;
55
56
57
58
59
60
61
62 private SSETest(int m, int n) {
63 A = new double[m][n];
64 flatA = new double[m * n];
65 x = new double[n];
66 }
67
68
69
70
71
72
73 public static void main(String[] args) {
74
75 int m = 500;
76 int n = 500;
77 if (args != null && args.length > 0) {
78 m = Integer.parseInt(args[0]);
79 n = Integer.parseInt(args[1]);
80 }
81
82 SSETest test = new SSETest(m, n);
83
84 int nLoops = 100000;
85
86 test.init(m, n);
87 Random r = new Random(0);
88 double temp = 0;
89 long time = 0;
90 for (int i = 1; i <= nLoops; i++) {
91 time -= System.nanoTime();
92 double[] y = test.matVec(test.A, test.x, m, n);
93 time += System.nanoTime();
94 int j = (int) floor(m * r.nextDouble());
95 temp += y[j];
96 if (i % 10000 == 0) {
97 System.out.println(" Nested: " + temp + " " + time * 1.0e-9);
98 time = 0;
99 }
100 }
101
102 test.init(m, n);
103 r = new Random(0);
104 temp = 0;
105 time = 0;
106 for (int i = 1; i <= nLoops; i++) {
107 time -= System.nanoTime();
108 double[] y = test.matVec(test.flatA, test.x, m, n);
109 time += System.nanoTime();
110 int j = (int) floor(m * r.nextDouble());
111 temp += y[j];
112 if (i % 10000 == 0) {
113 System.out.println(" Flat: " + temp + " " + time * 1.0e-9);
114 time = 0;
115 }
116 }
117 }
118
119 private void init(int n, int m) {
120 Random r = new Random(0);
121
122
123 for (int i = 0; i < m; i++) {
124 int idx = i * n;
125 for (int j = 0; j < n; j++) {
126 A[i][j] = r.nextDouble();
127 flatA[idx + j] = A[i][j];
128 }
129 }
130
131 for (int j = 0; j < n; j++) {
132 x[j] = r.nextDouble();
133 }
134 }
135
136
137
138
139
140
141
142
143
144
145 private double[] matVec(final double[][] A, final double[] x, final int m, final int n) {
146 double[] y = new double[m];
147 for (int i = 0; i < m; i++) {
148 for (int j = 0; j < n; j++) {
149
150 y[i] += A[i][j] * x[j];
151 }
152 }
153 return y;
154 }
155
156
157
158
159
160
161
162
163
164
165 private double[] matVec(final double[] A, final double[] x, final int m, final int n) {
166 double[] y = new double[m];
167 final int extra = n - n % 8;
168 final int ub = ((n / 8) * 8) - 1;
169 int idx = 0;
170 for (int i = 0; i < m; i++) {
171 double acc = 0;
172 for (int j = 0; j < ub; j += 8) {
173 int ptr = idx + j;
174 y[i] +=
175 A[ptr] * x[j]
176 + A[ptr + 1] * x[j + 1]
177 + A[ptr + 2] * x[j + 2]
178 + A[ptr + 3] * x[j + 3];
179 acc += A[ptr + 4] * x[j + 4]
180 + A[ptr + 5] * x[j + 5]
181 + A[ptr + 6] * x[j + 6]
182 + A[ptr + 7] * x[j + 7];
183 }
184 y[i] += acc;
185 for (int j = extra; j < n; j++) {
186 y[i] += A[idx + j] * x[j];
187 }
188 idx += n;
189 }
190 return y;
191 }
192 }