1 /*
2 * Copyright (c) 2026, Oracle and/or its affiliates. All rights reserved.
3 * DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER.
4 *
5 * This code is free software; you can redistribute it and/or modify it
6 * under the terms of the GNU General Public License version 2 only, as
7 * published by the Free Software Foundation. Oracle designates this
8 * particular file as subject to the "Classpath" exception as provided
9 * by Oracle in the LICENSE file that accompanied this code.
10 *
11 * This code is distributed in the hope that it will be useful, but WITHOUT
12 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
14 * version 2 for more details (a copy is included in the LICENSE file that
15 * accompanied this code).
16 *
17 * You should have received a copy of the GNU General Public License version
18 * 2 along with this work; if not, write to the Free Software Foundation,
19 * Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
20 *
21 * Please contact Oracle, 500 Oracle Parkway, Redwood Shores, CA 94065 USA
22 * or visit www.oracle.com if you need additional information or have any
23 * questions.
24 */
25 package dft;
26
27 import dft.Main.ComplexArray.Complex;
28 import hat.Accelerator;
29 import hat.Accelerator.Compute;
30 import hat.ComputeContext;
31 import hat.HATMath;
32
33 import static hat.KernelContext.*;
34 import hat.NDRange;
35 import hat.backend.Backend;
36 import hat.buffer.F32Array;
37 import hat.examples.common.ParseArgs;
38 import jdk.incubator.code.Reflect;
39 import optkl.ifacemapper.BoundSchema;
40 import optkl.ifacemapper.Buffer;
41 import optkl.ifacemapper.Schema;
42
43 import java.lang.invoke.MethodHandles;
44 import java.util.ArrayList;
45 import java.util.List;
46 import java.util.Random;
47 import java.util.stream.IntStream;
48
49 import static hat.examples.common.StatUtils.computeAverage;
50 import static hat.examples.common.StatUtils.computeSpeedup;
51 import static hat.examples.common.StatUtils.dumpStatsToCSVFile;
52 import static hat.examples.common.StatUtils.printCheckResult;
53
54 /**
55 * How to run?
56 *
57 * <p>
58 * With the OpenCL Backend:
59 * <code>
60 * java @.ffi-opencl-example dft.Main --size=<size> --iterations=<iterations> --verbose
61 * </code>
62 * </p>
63 *
64 * <p>
65 * With the CUDA Backend:
66 * <code>
67 * java @.ffi-cuda-example dft.Main --size=<size> --iterations=<iterations> --verbose
68 * </code>
69 *
70 * <p>
71 * Link to DFT: <a href="https://en.wikipedia.org/wiki/Discrete_Fourier_transform">link</a>
72 * </p>
73 *
74 * </p>
75 */
76 public class Main {
77
78 public static final float DELTA = 0.001f;
79
80 // Use a custom data structure for dealing with Array of Complex Numbers
81 public interface ComplexArray extends Buffer {
82 int length();
83
84 interface Complex extends Struct {
85 float real();
86 float imag();
87 void real(float real);
88 void imag(float imag);
89 }
90
91 Complex complex(long index);
92
93 Schema<ComplexArray> schema = Schema.of(ComplexArray.class, complex ->
94 complex.arrayLen("length")
95 .array("complex",
96 array -> array.fields("real", "imag")));
97
98 static ComplexArray create(Accelerator accelerator, int length) {
99 return BoundSchema.of(accelerator, schema, length).allocate();
100 }
101 }
102
103 @Reflect
104 private static void dftKernel(ComplexArray input, ComplexArray output) {
105 int size = input.length();
106 int idx = GIX();
107 if (idx < GSX()) {
108 float sumReal = 0.0f;
109 float sumImag = 0.0f;
110 for (int k = 0; k < size; k++) {
111 float angle = -2 * HATMath.PI * ((k * idx) % size) / size;
112 Complex complexInput = input.complex(k);
113 float cReal = HATMath.native_cosf(angle);
114 float cImag = HATMath.native_sinf(angle);
115 sumReal += (complexInput.real() * cReal) - (complexInput.imag() * cImag);
116 sumImag += (complexInput.real() * cImag) + (complexInput.imag() * cReal);
117 }
118 Complex complexOutput = output.complex(idx);
119 complexOutput.real(sumReal);
120 complexOutput.imag(sumImag);
121 }
122 }
123
124 @Reflect
125 private static void dftCompute(ComputeContext cc, ComplexArray input, ComplexArray output) {
126 var range = NDRange.of1D(input.length(), 256);
127 cc.dispatchKernel(range, () -> dftKernel( input, output));
128 }
129
130 @Reflect
131 private static void dftPlainKernel( F32Array inReal, F32Array inImag, F32Array outReal, F32Array outImag) {
132 int size = inReal.length();
133 int idx = GIX();
134 if (idx < GSX()) {
135 float sumReal = 0.0f;
136 float sumImag = 0.0f;
137 for (int k = 0; k < size; k++) {
138 float angle = -2 * HATMath.PI * ((idx * k) % size) / size;
139 float cReal = HATMath.native_cosf(angle);
140 float cImag = HATMath.native_sinf(angle);
141 sumReal += (inReal.array(k) * cReal) - (inImag.array(k) * cImag);
142 sumImag += (inReal.array(k) * cImag) + (inImag.array(k) * cReal);
143 }
144 outReal.array(idx, sumReal);
145 outImag.array(idx, sumImag);
146 }
147 }
148
149 @Reflect
150 private static void dftPlainCompute(ComputeContext cc, F32Array inReal, F32Array inImag, F32Array outReal, F32Array outImag) {
151 var range = NDRange.of1D(inReal.length(), 256);
152 cc.dispatchKernel(range, ()-> dftPlainKernel( inReal, inImag, outReal, outImag));
153 }
154
155 private static void dftJava(ComplexArray input, ComplexArray output) {
156 int size = input.length();
157 for (int k = 0; k < size; k++) {
158 Complex complexOutput = output.complex(k);
159 complexOutput.real(0.0f);
160 complexOutput.imag(0.0f);
161 float sumReal = 0.0f;
162 float sumImag = 0.0f;
163 for (int j = 0; j < size; j++) {
164 float angle = -2 * HATMath.PI * ((j * k) % size) / size;
165 Complex complexInput = input.complex(j);
166 float cReal = HATMath.cosf(angle);
167 float cImag = HATMath.sinf(angle);
168 sumReal += (complexInput.real() * cReal) - (complexInput.imag() * cImag);
169 sumImag += (complexInput.real() * cImag) + (complexInput.imag() * cReal);
170 }
171 complexOutput.real(sumReal);
172 complexOutput.imag(sumImag);
173 }
174 }
175
176 private static void dftJavaStreams(ComplexArray input, ComplexArray output) {
177 int size = input.length();
178 IntStream.range(0, size).parallel().forEach(idx -> {
179 float sumReal = 0.0f;
180 float sumImag = 0.0f;
181 for (int k = 0; k < size; k++) {
182 float angle = -2 * HATMath.PI * ((idx * k) % size) / size;
183 Complex complexInput = input.complex(k);
184 float cReal = HATMath.cosf(angle);
185 float cImag = HATMath.sinf(angle);
186 sumReal += (complexInput.real() * cReal) - (complexInput.imag() * cImag);
187 sumImag += (complexInput.real() * cImag) + (complexInput.imag() * cReal);
188 }
189 Complex complexOutput = output.complex(idx);
190 complexOutput.real(sumReal);
191 complexOutput.imag(sumImag);
192 });
193 }
194
195 private static boolean checkResult(ComplexArray expected, ComplexArray obtained) {
196 for (int i = 0; i < expected.length(); i++) {
197 if (Math.abs(expected.complex(i).real() - obtained.complex(i).real()) > DELTA) {
198 IO.println(expected.complex(i).real() + " vs " + obtained.complex(i).real());
199 return false;
200 }
201 if (Math.abs(expected.complex(i).imag() - obtained.complex(i).imag()) > DELTA) {
202 IO.println(expected.complex(i).imag() + " vs " + obtained.complex(i).imag());
203 return false;
204 }
205 }
206 return true;
207 }
208
209 // Just for debugging
210 private static void printSignal(ComplexArray signal) {
211 for (int i = 0; i < signal.length(); i++) {
212 IO.println(signal.complex(i).real() + "," + signal.complex(i).imag());
213 }
214 }
215
216 private static boolean checkResult(ComplexArray outputSeq, F32Array outReal, F32Array outImag) {
217 for (int i = 0; i < outputSeq.length(); i++) {
218 if (Math.abs(outputSeq.complex(i).real() - outReal.array(i)) > DELTA) {
219 IO.println(outputSeq.complex(i).real() + " vs " + outReal.array(i));
220 return false;
221 }
222 if (Math.abs(outputSeq.complex(i).imag() - outImag.array(i)) > DELTA) {
223 IO.println(outputSeq.complex(i).imag() + " vs " + outImag.array(i));
224 return false;
225 }
226 }
227 return true;
228 }
229
230 static void main(String[] args) {
231 IO.println("=========================================");
232 IO.println("Example: Discrete Fourier Transform (DFT)");
233 IO.println("=========================================");
234
235 int size = 32768;
236 int iterations = 100;
237 ParseArgs parseArgs = new ParseArgs(args);
238 ParseArgs.Options options = parseArgs.parseWithDefaults(size, iterations);
239
240 boolean verbose = options.verbose();
241 size = options.size();
242 iterations = options.iterations();
243 boolean skipSequential = options.skipSequential();
244 IO.println("Input Size = " + size);
245 IO.println("Num Iterations = " + iterations);
246
247 var lookup = MethodHandles.lookup();
248 var accelerator = new Accelerator(lookup, Backend.FIRST);
249 ComplexArray input = ComplexArray.create(accelerator, size);
250 ComplexArray outputSeq = ComplexArray.create(accelerator, size);
251 ComplexArray.create(accelerator, size);
252 ComplexArray outputStreams = ComplexArray.create(accelerator, size);
253 ComplexArray outputHAT = ComplexArray.create(accelerator, size);
254
255 F32Array inReal = F32Array.create(accelerator, size);
256 F32Array inImag = F32Array.create(accelerator, size);
257 F32Array outReal = F32Array.create(accelerator, size);
258 F32Array outImag = F32Array.create(accelerator, size);
259
260 // Initialize
261 Random r = new Random(71);
262 for (int i = 0; i < size; i++) {
263 Complex c = input.complex(i);
264 c.real(r.nextFloat());
265 c.imag(r.nextFloat());
266 inReal.array(i, c.real());
267 inImag.array(i, c.imag());
268 }
269
270 List<Long> timersJavaDFT = new ArrayList<>();
271 List<Long> timersStreams = new ArrayList<>();
272 List<Long> timersDFTHat = new ArrayList<>();
273 List<Long> timersDFTHatFlatten = new ArrayList<>();
274
275 // Run Java sequential version, DFT
276 if (!skipSequential) {
277 for (int i = 0; i < iterations; i++) {
278 long start = System.nanoTime();
279 dftJava(input, outputSeq);
280 long end = System.nanoTime();
281 timersJavaDFT.add((end - start));
282 if (verbose) {
283 IO.println("[Timer] Java DFT: " + (end - start));
284 }
285 }
286 }
287
288 // Run Java Parallel Stream Version of the DFT
289 for (int i = 0; i < iterations; i++) {
290 long start = System.nanoTime();
291 dftJavaStreams(input, outputStreams);
292 long end = System.nanoTime();
293 timersStreams.add((end-start));
294 if (verbose) {
295 IO.println("[Timer] Parallel Stream: " + (end-start));
296 }
297 }
298
299 // HAT: Initial version (DFT)
300 for (int i = 0; i < iterations; i++) {
301 long start = System.nanoTime();
302 accelerator.compute((@Reflect Compute) computeContext -> dftCompute(computeContext, input, outputHAT));
303 long end = System.nanoTime();
304 timersDFTHat.add((end-start));
305 if (verbose) {
306 IO.println("[Timer] HAT with CDS: " + (end-start));
307 }
308 }
309
310 // HAT: DFT using plain arrays instead of a custom data structure
311 for (int i = 0; i < iterations; i++) {
312 long start = System.nanoTime();
313 accelerator.compute((@Reflect Compute) computeContext -> dftPlainCompute(computeContext, inReal, inImag, outReal, outImag));
314 long end = System.nanoTime();
315 timersDFTHatFlatten.add((end - start));
316 if (verbose) {
317 IO.println("[Timer] HAT-Plain: " + (end - start));
318 }
319 }
320
321 // Check results
322 ComplexArray baseline = skipSequential ? outputStreams : outputSeq;
323 boolean isStreamCorrect = checkResult(baseline, outputStreams);
324 boolean isHATCorrect = checkResult(baseline, outputHAT);
325 boolean isHATPlainCorrect = checkResult(baseline, outReal, outImag);
326 printCheckResult(isStreamCorrect, "Java-Stream");
327 printCheckResult(isHATCorrect, "HAT-Naive");
328 printCheckResult(isHATPlainCorrect, "HAT-Plain");
329
330 // Print Performance Metrics
331 final int skip = iterations / 2;
332 double averageJavaTimer = computeAverage(timersJavaDFT, skip);
333 double averageJavaStreamTimer = computeAverage(timersStreams, skip);
334 double averageHATTimers = computeAverage(timersDFTHat, skip);
335 double averageHATTimersFlatten = computeAverage(timersDFTHatFlatten, skip);
336
337 IO.println("\nAverage elapsed time:");
338 IO.println("Average Java-Seq DFT : " + averageJavaTimer);
339 IO.println("Average Java-Streams DFT : " + averageJavaStreamTimer);
340 IO.println("Average HAT DFT : " + averageHATTimers);
341 IO.println("Average HAT DFT (Flatten) : " + averageHATTimersFlatten);
342
343 if (!skipSequential) {
344 IO.println("\nSpeedups vs Java:");
345 IO.println("Java / Java Parallel Stream : " + computeSpeedup(averageJavaTimer, averageJavaStreamTimer) + "x");
346 IO.println("Java / HAT : " + computeSpeedup(averageJavaTimer, averageHATTimers) + "x");
347 IO.println("Java / HAT-Flatten : " + computeSpeedup(averageJavaTimer, averageHATTimersFlatten) + "x");
348 }
349
350 IO.println("\nSpeedups vs Java Parallel Streams:");
351 IO.println("Java / HAT : " + computeSpeedup(averageJavaStreamTimer, averageHATTimers) + "x");
352 IO.println("Java / HAT-Flatten : " + computeSpeedup(averageJavaStreamTimer, averageHATTimersFlatten) + "x");
353
354 IO.println("\nSpeedups vs HAT-DFT:");
355 IO.println("HAT / HAT-Flatten : " + computeSpeedup(averageHATTimers, averageHATTimersFlatten) + "x");
356
357 // Write CSV table with all results
358 List<List<Long>> timers = List.of(timersJavaDFT, timersStreams, timersDFTHat, timersDFTHatFlatten);
359 List<String> header = List.of("Java-fp32-" + size, "Streams-fp32-" + size, "HAT-UDT-fp32-" + size, "HAT-Plain-fp32-" + size);
360 String fileName = "table-results-dft-" + size + ".csv";
361 if (skipSequential) {
362 timers = List.of(timersStreams, timersDFTHat, timersDFTHatFlatten);
363 header = List.of("Streams-fp32-" + size, "HAT-UDT-fp32-" + size, "HAT-Plain-fp32-" + size);
364 }
365 dumpStatsToCSVFile(timers, header, fileName);
366 }
367 }