View Javadoc
1   // ******************************************************************************
2   //
3   // Title:       Force Field X.
4   // Description: Force Field X - Software for Molecular Biophysics.
5   // Copyright:   Copyright (c) Michael J. Schnieders 2001-2025.
6   //
7   // This file is part of Force Field X.
8   //
9   // Force Field X is free software; you can redistribute it and/or modify it
10  // under the terms of the GNU General Public License version 3 as published by
11  // the Free Software Foundation.
12  //
13  // Force Field X is distributed in the hope that it will be useful, but WITHOUT
14  // ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
15  // FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
16  // details.
17  //
18  // You should have received a copy of the GNU General Public License along with
19  // Force Field X; if not, write to the Free Software Foundation, Inc., 59 Temple
20  // Place, Suite 330, Boston, MA 02111-1307 USA
21  //
22  // Linking this library statically or dynamically with other modules is making a
23  // combined work based on this library. Thus, the terms and conditions of the
24  // GNU General Public License cover the whole combination.
25  //
26  // As a special exception, the copyright holders of this library give you
27  // permission to link this library with independent modules to produce an
28  // executable, regardless of the license terms of these independent modules, and
29  // to copy and distribute the resulting executable under terms of your choice,
30  // provided that you also meet, for each linked independent module, the terms
31  // and conditions of the license of that module. An independent module is a
32  // module which is not derived from or based on this library. If you modify this
33  // library, you may extend this exception to your version of the library, but
34  // you are not obligated to do so. If you do not wish to do so, delete this
35  // exception statement from your version.
36  //
37  // ******************************************************************************
38  package ffx.numerics.fft;
39  
40  import ffx.utilities.FFXTest;
41  import org.junit.Test;
42  import org.junit.runner.RunWith;
43  import org.junit.runners.Parameterized;
44  import org.junit.runners.Parameterized.Parameters;
45  
46  import java.util.Arrays;
47  import java.util.Collection;
48  import java.util.Random;
49  
50  import static org.junit.Assert.assertEquals;
51  
52  /**
53   * @author Michael J. Schnieders
54   */
55  @RunWith(Parameterized.class)
56  public class ComplexTest extends FFXTest {
57  
58    private final int n;
59    private final String info;
60    private final boolean preferred;
61    private final double[] data;
62    private final double[] dataBlocked;
63    private final double[] orig;
64    private final double[] origBlocked;
65    private final double[] dft;
66  
67    public ComplexTest(String info, int n, boolean preferred) {
68      this.info = info;
69      this.n = n;
70      this.preferred = preferred;
71      data = new double[n * 2];
72      dataBlocked = new double[n * 2];
73      orig = new double[n * 2];
74      origBlocked = new double[n * 2];
75      dft = new double[n * 2];
76      Random r = new Random();
77      for (int i = 0; i < n; i++) {
78        double d = r.nextDouble();
79        data[i * 2] = d;
80        orig[i * 2] = d;
81        dataBlocked[i] = d;
82        origBlocked[i] = d;
83      }
84    }
85  
86    @Parameters
87    public static Collection<Object[]> data() {
88      return Arrays.asList(
89          new Object[][]{
90              {"Test n = 162 with factors [6, 3, 3, 3]", 162, true},
91              {"Test n = 160 with factors [5, 4, 4, 2]", 160, true},
92              {"Test n = 120 with factors [6, 5, 4]", 120, true},
93              {"Test n = 64 with factors [4, 4, 4]", 64, true},
94              {"Test n = 48 with factors [6, 4, 2]", 48, true},
95              {"Test n = 21 with factors [7, 3]", 21, true},
96              {"Test n = 38 with factors [2, 19]", 38, false},
97              {"Test n = 22 with factors [2, 11]", 22, false},
98          });
99    }
100 
101   /**
102    * Test of preferredDimension method, of class Complex.
103    */
104   @Test
105   public void testPreferredDimension() {
106     boolean result = Complex.preferredDimension(n);
107     assertEquals(info, preferred, result);
108   }
109 
110   /**
111    * Test of fft method, of class Complex.
112    */
113   @Test
114   public void testInterleavedScalar() {
115     double tolerance = 1.0e-11;
116 
117     int offset = 0;
118     int stride = 2;
119     Complex complex = new Complex(n);
120     complex.setUseSIMD(false);
121 
122     // System.out.println(info + "\n Factors " + Arrays.toString(complex.getFactors()));
123 
124     long dftTime = System.nanoTime();
125     Complex.dft(data, dft);
126     dftTime = System.nanoTime() - dftTime;
127     String dftString = " DFT Time: " + dftTime * 1.0e-9 + " s\n";
128 
129     long fftTime = System.nanoTime();
130     complex.fft(data, offset, stride);
131     fftTime = System.nanoTime() - fftTime;
132     String fftString = " FFT Time: " + fftTime * 1.0e-9 + " s";
133 
134     // Test the FFT is equals the DFT result.
135     for (int i = 0; i < 2 * n; i++) {
136       assertEquals(" Forward " + info + " at position: " + i, dft[i], data[i], tolerance);
137     }
138 
139     // The FFT is faster than the DFT.
140     String message = fftString + dftString;
141     // assertTrue(message, fftTime < dftTime);
142 
143     // Test that X = IFFT(FFT(X)).
144     complex.inverse(data, offset, stride);
145     for (int i = 0; i < n; i++) {
146       double orig = this.orig[i * 2];
147       double actual = data[i * 2];
148       assertEquals(" IFFT(FFT(X)) " + info + " at position: " + i, orig, actual, tolerance);
149     }
150   }
151 
152   /**
153    * Test of fft method, of class Complex.
154    */
155   @Test
156   public void testBlockedScalar() {
157     if (!this.preferred) {
158       return;
159     }
160     double tolerance = 1.0e-11;
161     int offset = 0;
162     int stride = 1;
163     Complex complex = new Complex(n, DataLayout1D.BLOCKED, n);
164     complex.setUseSIMD(false);
165 
166     // System.out.println(info + "\n Factors " + Arrays.toString(complex.getFactors()));
167 
168     long dftTime = System.nanoTime();
169     Complex.dftBlocked(dataBlocked, dft);
170     dftTime = System.nanoTime() - dftTime;
171     String dftString = " DFT Time: " + dftTime * 1.0e-9 + " s\n";
172 
173     long fftTime = System.nanoTime();
174     complex.fft(dataBlocked, offset, stride);
175     fftTime = System.nanoTime() - fftTime;
176     String fftString = " FFT Time: " + fftTime * 1.0e-9 + " s";
177 
178     // Test the FFT is equals the DFT result.
179     for (int i = 0; i < 2 * n; i++) {
180       assertEquals(" Forward " + info + " at position: " + i, dft[i], dataBlocked[i], tolerance);
181     }
182 
183     // The FFT is faster than the DFT.
184     String message = fftString + dftString;
185     // assertTrue(message, fftTime < dftTime);
186 
187     // Test that X = IFFT(FFT(X)).
188     complex.inverse(dataBlocked, offset, stride);
189     for (int i = 0; i < n; i++) {
190       double orig = this.origBlocked[i];
191       double actual = dataBlocked[i];
192       assertEquals(" IFFT(FFT(X)) " + info + " at position: " + i, orig, actual, tolerance);
193     }
194   }
195 
196   /**
197    * Test of fft method, of class Complex.
198    */
199   @Test
200   public void testInterleavedSIMD() {
201     double tolerance = 1.0e-11;
202 
203     int offset = 0;
204     int stride = 2;
205     Complex complex = new Complex(n);
206     complex.setUseSIMD(true);
207 
208     // System.out.println(info + "\n Factors " + Arrays.toString(complex.getFactors()));
209 
210     long dftTime = System.nanoTime();
211     Complex.dft(data, dft);
212     dftTime = System.nanoTime() - dftTime;
213     String dftString = " DFT Time: " + dftTime * 1.0e-9 + " s\n";
214 
215     long fftTime = System.nanoTime();
216     complex.fft(data, offset, stride);
217     fftTime = System.nanoTime() - fftTime;
218     String fftString = " FFT Time: " + fftTime * 1.0e-9 + " s";
219 
220     // Test the FFT is equals the DFT result.
221     for (int i = 0; i < 2 * n; i++) {
222       assertEquals(" Forward " + info + " at position: " + i, dft[i], data[i], tolerance);
223     }
224 
225     // The FFT is faster than the DFT.
226     String message = fftString + dftString;
227     // assertTrue(message, fftTime < dftTime);
228 
229     // Test that X = IFFT(FFT(X)).
230     complex.inverse(data, offset, stride);
231     for (int i = 0; i < n; i++) {
232       double orig = this.orig[i * 2];
233       double actual = data[i * 2];
234       assertEquals(" IFFT(FFT(X)) " + info + " at position: " + i, orig, actual, tolerance);
235     }
236   }
237 
238   /**
239    * Test of fft method, of class Complex.
240    */
241   @Test
242   public void testBlockedSIMD() {
243     if (!this.preferred) {
244       return;
245     }
246     double tolerance = 1.0e-11;
247     int offset = 0;
248     int stride = 1;
249     Complex complex = new Complex(n, DataLayout1D.BLOCKED, n);
250     complex.setUseSIMD(true);
251 
252     long dftTime = System.nanoTime();
253     Complex.dftBlocked(dataBlocked, dft);
254     dftTime = System.nanoTime() - dftTime;
255     String dftString = " DFT Time: " + dftTime * 1.0e-9 + " s\n";
256 
257     long fftTime = System.nanoTime();
258     complex.fft(dataBlocked, offset, stride);
259     fftTime = System.nanoTime() - fftTime;
260     String fftString = " FFT Time: " + fftTime * 1.0e-9 + " s";
261 
262     // Test the FFT is equals the DFT result.
263     for (int i = 0; i < 2 * n; i++) {
264       assertEquals(" Forward " + info + " at position: " + i, dft[i], dataBlocked[i], tolerance);
265     }
266 
267     // The FFT is faster than the DFT.
268     String message = fftString + dftString;
269     // assertTrue(message, fftTime < dftTime);
270 
271     // Test that X = IFFT(FFT(X)).
272     complex.inverse(dataBlocked, offset, stride);
273     for (int i = 0; i < n; i++) {
274       double orig = this.origBlocked[i];
275       double actual = dataBlocked[i];
276       assertEquals(" IFFT(FFT(X)) " + info + " at position: " + i, orig, actual, tolerance);
277     }
278   }
279 
280   /**
281    * Test of fft method, of class Complex, using N FFTs packed next to each other.
282    * <p>
283    * This test is designed to mimic a 2D transform, so we'll pack N FFTs into a 1D array of size (2*N)^2
284    */
285   @Test
286   public void testInterleavedScalarNFFT() {
287     double tolerance = 1.0e-11;
288 
289     int nFFTs = n;
290     int nextFFT = 2 * n;
291     int offset = 0;
292     int stride = 2;
293     int im = 1;
294     Complex complex = new Complex(n, DataLayout1D.INTERLEAVED, im, nFFTs);
295     complex.setUseSIMD(false);
296 
297     long dftTime = System.nanoTime();
298     Complex.dft(data, dft);
299     dftTime = System.nanoTime() - dftTime;
300     String dftString = " DFT Time: " + dftTime * 1.0e-9 + " s\n";
301 
302     long fftTime = System.nanoTime();
303     double[] data2 = new double[2 * n * nFFTs];
304     int index = 0;
305     for (int ii = 0; ii < 2 * n; ii += stride) {
306       for (int i = 0; i < nFFTs; i++) {
307         data2[index] = data[ii];
308         //data2[index + im] = data[ii + im];
309         index += stride;
310       }
311     }
312 
313     complex.fft(data2, offset, stride, nextFFT);
314     fftTime = System.nanoTime() - fftTime;
315     String fftString = " FFT Time: " + fftTime * 1.0e-9 + " s";
316 
317     // Test the FFT is equals the DFT result.
318     index = 0;
319     for (int ii = 0; ii < n; ii++) {
320       for (int i = 0; i < nFFTs; i++) {
321         assertEquals(" FFT " + i + " Forward " + info + " at position: " + ii,
322             dft[2 * ii], data2[index], tolerance);
323         assertEquals(" FFT " + i + " Forward " + info + " at position: " + ii,
324             dft[2 * ii + 1], data2[index + 1], tolerance);
325         index += 2;
326       }
327     }
328 
329     // The FFT is faster than the DFT.
330     String message = fftString + dftString;
331     // assertTrue(message, fftTime < dftTime);
332 
333     // Test that X = IFFT(FFT(X)).
334     complex.inverse(data2, offset, stride, nextFFT);
335 
336     // Test the FFT is equals the DFT result.
337     index = 0;
338     for (int ii = 0; ii < nextFFT; ii += 2) {
339       for (int i = 0; i < nFFTs; i++) {
340         assertEquals(" IFFT(FFT(X)) " + info + " at position: " + ii + " FFT " + i,
341             orig[ii], data2[index], tolerance);
342         index += 2;
343       }
344     }
345   }
346 
347   /**
348    * Test of fft method, of class Complex, using N FFTs packed next to each other.
349    * <p>
350    * This test is designed to mimic a 2D transform, so we'll pack N FFTs into a 1D array of size (2*N)^2
351    */
352   @Test
353   public void testInterleavedSIMDNFFT() {
354     double tolerance = 1.0e-11;
355 
356     int nFFTs = n;
357     int nextFFT = 2 * n;
358     int offset = 0;
359     int stride = 2;
360     int im = 1;
361     Complex complex = new Complex(n, DataLayout1D.INTERLEAVED, im, nFFTs);
362     complex.setUseSIMD(true);
363 
364     long dftTime = System.nanoTime();
365     Complex.dft(data, dft);
366     dftTime = System.nanoTime() - dftTime;
367     String dftString = " DFT Time: " + dftTime * 1.0e-9 + " s\n";
368 
369     long fftTime = System.nanoTime();
370     double[] data2 = new double[2 * n * nFFTs];
371     int index = 0;
372     for (int ii = 0; ii < 2 * n; ii += stride) {
373       for (int i = 0; i < nFFTs; i++) {
374         data2[index] = data[ii];
375         //data2[index + im] = data[ii + im];
376         index += stride;
377       }
378     }
379 
380     complex.fft(data2, offset, stride, nextFFT);
381     fftTime = System.nanoTime() - fftTime;
382     String fftString = " FFT Time: " + fftTime * 1.0e-9 + " s";
383 
384     // Test the FFT is equals the DFT result.
385     index = 0;
386     for (int ii = 0; ii < n; ii++) {
387       for (int i = 0; i < nFFTs; i++) {
388         assertEquals(" FFT " + i + " Forward " + info + " at position: " + ii,
389             dft[2 * ii], data2[index], tolerance);
390         assertEquals(" FFT " + i + " Forward " + info + " at position: " + ii,
391             dft[2 * ii + 1], data2[index + 1], tolerance);
392         index += 2;
393       }
394     }
395 
396     // The FFT is faster than the DFT.
397     String message = fftString + dftString;
398     // assertTrue(message, fftTime < dftTime);
399 
400     // Test that X = IFFT(FFT(X)).
401     complex.inverse(data2, offset, stride, nextFFT);
402 
403     // Test the FFT is equals the DFT result.
404     index = 0;
405     for (int ii = 0; ii < nextFFT; ii += 2) {
406       for (int i = 0; i < nFFTs; i++) {
407         assertEquals(" IFFT(FFT(X)) " + info + " at position: " + ii + " FFT " + i,
408             orig[ii], data2[index], tolerance);
409         index += 2;
410       }
411     }
412   }
413 
414   /**
415    * Test of fft method, of class Complex, using N FFTs packed next to each other.
416    * <p>
417    * This test is designed to mimic a 2D transform, so we'll pack N FFTs into a 1D array of size (2*N)^2
418    */
419   @Test
420   public void testBlockedScalarNFFT() {
421     if (!this.preferred) {
422       return;
423     }
424     double tolerance = 1.0e-11;
425 
426     int nFFTs = n;
427     int nextFFT = 2 * n;
428     int offset = 0;
429     int stride = 1;
430     int im = n * nFFTs;
431     Complex complex = new Complex(n, DataLayout1D.BLOCKED, im, nFFTs);
432     complex.setUseSIMD(false);
433 
434     long dftTime = System.nanoTime();
435     Complex.dftBlocked(dataBlocked, dft);
436     dftTime = System.nanoTime() - dftTime;
437     String dftString = " DFT Time: " + dftTime * 1.0e-9 + " s\n";
438 
439     long fftTime = System.nanoTime();
440     double[] data2 = new double[2 * n * nFFTs];
441     int index = 0;
442     for (int ii = 0; ii < n; ii += stride) {
443       for (int i = 0; i < nFFTs; i++) {
444         data2[index++] = dataBlocked[ii];
445       }
446     }
447 
448     complex.fft(data2, offset, stride, nextFFT);
449     fftTime = System.nanoTime() - fftTime;
450     String fftString = " FFT Time: " + fftTime * 1.0e-9 + " s";
451 
452     // Test the FFT is equals the DFT result.
453     index = 0;
454     for (int ii = 0; ii < n; ii++) {
455       for (int i = 0; i < nFFTs; i++) {
456         assertEquals(" FFT " + i + " Forward " + info + " at position: " + ii,
457             dft[ii], data2[index], tolerance);
458         assertEquals(" FFT " + i + " Forward " + info + " at position: " + ii,
459             dft[ii + n], data2[index + im], tolerance);
460         index++;
461       }
462     }
463 
464     // The FFT is faster than the DFT.
465     String message = fftString + dftString;
466     // assertTrue(message, fftTime < dftTime);
467 
468     // Test that X = IFFT(FFT(X)).
469     complex.inverse(data2, offset, stride, nextFFT);
470 
471     // Test the FFT is equals the DFT result.
472     index = 0;
473     for (int ii = 0; ii < nextFFT; ii += 2) {
474       for (int i = 0; i < nFFTs; i++) {
475         assertEquals(" IFFT(FFT(X)) " + info + " at position: " + ii + " FFT " + i,
476             orig[ii], data2[index], tolerance);
477         index++;
478       }
479     }
480   }
481 
482   /**
483    * Test of fft method, of class Complex, using N FFTs packed next to each other.
484    * <p>
485    * This test is designed to mimic a 2D transform, so we'll pack N FFTs into a 1D array of size (2*N)^2
486    */
487   @Test
488   public void testBlockedSIMDNFFT() {
489     if (!this.preferred) {
490       return;
491     }
492     double tolerance = 1.0e-11;
493 
494     int nFFTs = n;
495     int nextFFT = 2 * n;
496     int offset = 0;
497     int stride = 1;
498     int im = n * nFFTs;
499     Complex complex = new Complex(n, DataLayout1D.BLOCKED, im, nFFTs);
500     complex.setUseSIMD(true);
501 
502     long dftTime = System.nanoTime();
503     Complex.dftBlocked(dataBlocked, dft);
504     dftTime = System.nanoTime() - dftTime;
505     String dftString = " DFT Time: " + dftTime * 1.0e-9 + " s\n";
506 
507     long fftTime = System.nanoTime();
508     double[] data2 = new double[2 * n * nFFTs];
509     int index = 0;
510     for (int ii = 0; ii < n; ii += stride) {
511       for (int i = 0; i < nFFTs; i++) {
512         data2[index++] = dataBlocked[ii];
513       }
514     }
515 
516     complex.fft(data2, offset, stride, nextFFT);
517     fftTime = System.nanoTime() - fftTime;
518     String fftString = " FFT Time: " + fftTime * 1.0e-9 + " s";
519 
520     // Test the FFT is equals the DFT result.
521     index = 0;
522     for (int ii = 0; ii < n; ii++) {
523       for (int i = 0; i < nFFTs; i++) {
524         assertEquals(" FFT " + i + " Forward " + info + " at position: " + ii,
525             dft[ii], data2[index], tolerance);
526         assertEquals(" FFT " + i + " Forward " + info + " at position: " + ii,
527             dft[ii + n], data2[index + im], tolerance);
528         index++;
529       }
530     }
531 
532     // The FFT is faster than the DFT.
533     String message = fftString + dftString;
534     // assertTrue(message, fftTime < dftTime);
535 
536     // Test that X = IFFT(FFT(X)).
537     complex.inverse(data2, offset, stride, nextFFT);
538 
539     // Test the FFT is equals the DFT result.
540     index = 0;
541     for (int ii = 0; ii < nextFFT; ii += 2) {
542       for (int i = 0; i < nFFTs; i++) {
543         assertEquals(" IFFT(FFT(X)) " + info + " at position: " + ii + " FFT " + i,
544             orig[ii], data2[index], tolerance);
545         index++;
546       }
547     }
548   }
549 
550 }