1/* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2
3/*
4 Copyright (C) 2006 Joseph Wang
5 Copyright (C) 2009 Liquidnet Holdings, Inc.
6
7 This file is part of QuantLib, a free-software/open-source library
8 for financial quantitative analysts and developers - http://quantlib.org/
9
10 QuantLib is free software: you can redistribute it and/or modify it
11 under the terms of the QuantLib license. You should have received a
12 copy of the license along with this program; if not, please email
13 <quantlib-dev@lists.sf.net>. The license is also available online at
14 <http://quantlib.org/license.shtml>.
15
16 This program is distributed in the hope that it will be useful, but WITHOUT
17 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
18 FOR A PARTICULAR PURPOSE. See the license for more details.
19*/
20
21#include "fastfouriertransform.hpp"
22#include "utilities.hpp"
23#include <ql/math/fastfouriertransform.hpp>
24#include <ql/math/array.hpp>
25#include <complex>
26#include <vector>
27#include <functional>
28
29using namespace QuantLib;
30using namespace boost::unit_test_framework;
31
32void FastFourierTransformTest::testSimple() {
33 BOOST_TEST_MESSAGE("Testing complex direct FFT...");
34 typedef std::complex<Real> cx;
35 cx a[] = { cx(0,0), cx(1,1), cx(3,3), cx(4,4),
36 cx(4,4), cx(3,3), cx(1,1), cx(0,0) };
37 cx b[8];
38 FastFourierTransform fft(3);
39 fft.transform(inBegin: a, inEnd: a+8, out: b);
40 cx expected[] = { cx(16,16), cx(-4.8284,-11.6569),
41 cx(0,0), cx(-0.3431,0.8284),
42 cx(0,0), cx(0.8284, -0.3431),
43 cx(0,0), cx(-11.6569,-4.8284) };
44 for (size_t i = 0; i<8; i++) {
45 if ((std::fabs(x: b[i].real() - expected[i].real()) > 1.0e-2) ||
46 (std::fabs(x: b[i].imag() - expected[i].imag()) > 1.0e-2))
47 BOOST_ERROR("Convolution(" << i << ")\n"
48 << std::setprecision(4) << std::scientific
49 << " calculated: " << b[i] << "\n"
50 << " expected: " << expected[i]);
51 }
52}
53
54void FastFourierTransformTest::testInverse() {
55 BOOST_TEST_MESSAGE("Testing convolution via inverse FFT...");
56 Array x(3);
57 x[0] = 1;
58 x[1] = 2;
59 x[2] = 3;
60
61 size_t order = FastFourierTransform::min_order(inputSize: x.size())+1;
62 FastFourierTransform fft(order);
63 size_t nFrq = fft.output_size();
64 std::vector< std::complex<Real> > ft (nFrq);
65 std::vector< Real > tmp (nFrq);
66 std::complex<Real> z = std::complex<Real>();
67
68 fft.inverse_transform(inBegin: x.begin(), inEnd: x.end(), out: ft.begin());
69 for (Size i=0; i<nFrq; ++i) {
70 tmp[i] = std::norm(z: ft[i]);
71 ft[i] = z;
72 }
73 fft.inverse_transform(inBegin: tmp.begin(), inEnd: tmp.end(), out: ft.begin());
74
75 // 0
76 Real calculated = ft[0].real() / nFrq;
77 Real expected = x[0]*x[0] + x[1]*x[1] + x[2]*x[2];
78 if (fabs (x: calculated - expected) > 1.0e-10)
79 BOOST_ERROR("Convolution(0)\n"
80 << std::setprecision(16) << std::scientific
81 << " calculated: " << calculated << "\n"
82 << " expected: " << expected);
83
84 // 1
85 calculated = ft[1].real() / nFrq;
86 expected = x[0]*x[1] + x[1]*x[2];
87 if (fabs (x: calculated - expected) > 1.0e-10)
88 BOOST_ERROR("Convolution(1)\n"
89 << std::setprecision(16) << std::scientific
90 << " calculated: " << calculated << "\n"
91 << " expected: " << expected);
92
93 // 2
94 calculated = ft[2].real() / nFrq;
95 expected = x[0]*x[2];
96 if (fabs (x: calculated - expected) > 1.0e-10)
97 BOOST_ERROR("Convolution(1)\n"
98 << std::setprecision(16) << std::scientific
99 << " calculated: " << calculated << "\n"
100 << " expected: " << expected);
101
102}
103
104
105
106test_suite* FastFourierTransformTest::suite() {
107 auto* suite = BOOST_TEST_SUITE("fast fourier transform tests");
108 suite->add(QUANTLIB_TEST_CASE(&FastFourierTransformTest::testSimple));
109 suite->add(QUANTLIB_TEST_CASE(&FastFourierTransformTest::testInverse));
110 return suite;
111}
112
113

source code of quantlib/test-suite/fastfouriertransform.cpp