File fft.hxx#
FFT routines
Copyright 2010 B.D.Dudson, S.Farley, M.V.Umansky, X.Q.Xu
Contact: Ben Dudson, bd512@york.ac.uk
This file is part of BOUT++.
BOUT++ is free software: you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.
BOUT++ is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.
You should have received a copy of the GNU Lesser General Public License along with BOUT++. If not, see http://www.gnu.org/licenses/.
Functions
-
BOUT_ENUM_CLASS(FFT_MEASUREMENT_FLAG, estimate, measure, exhaustive)#
-
inline void rfft(const BoutReal *in, int length, dcomplex *out)#
Returns the fft of a real signal using fftw_forward
The fftw_forward returns out_k = sum_{j=0}^(length-1) in_j*exp(-2*pi*j*k*sqrt(-1)/length)
Thus, out_k must be divided by ‘length’ in order for DFT[IDFT[in]] = in where IDFT is the inverse fourier transform. See the the fftw user manual for details.
- Parameters:
in – [in] Pointer to the 1D array to take the fourier transform of
length – [in] Number of points in the input array
out – [out] Pointer to the complex 1D array which is the FFT of in
-
inline void irfft(const dcomplex *in, int length, BoutReal *out)#
Take the inverse fft of signal where the outputs are only reals.
This is done through a call to fftw_plan_dft_c2r_1d which is calling fftw_backwards.
That is out_k = sum_{j=0}^(length-1) in_j*exp(2*pi*j*k*sqrt(-1)/length)
See the the fftw user manual for details.
- Parameters:
in – [in] Pointer to the 1D array to take the inverse fourier transform of
length – [in] Number of points in the input array
out – [out] Pointer to the complex 1D array which is IFFTed
-
namespace bout
Provides access to the Hypre library, handling initialisation and finalisation.
Usage
#include <bout/hyprelib.hxx>
class MyClass { public:
private: HypreLib lib; };
This will then automatically initialise Hypre the first time an object is created, and finalise it when the last object is destroyed.
Copyright 2012 B.D.Dudson, S.Farley, M.V.Umansky, X.Q.Xu
Contact: Ben Dudson, bd512@york.ac.uk
This file is part of BOUT++.
BOUT++ is free software: you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.
BOUT++ is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.
You should have received a copy of the GNU Lesser General Public License along with BOUT++. If not, see http://www.gnu.org/licenses/.
Information about the version of BOUT++
The build system will update this file on every commit, which may result in files that include it getting rebuilt. Therefore it should be included in as few places as possible
Information about the version of BOUT++
The build system will update this file at configure-time
SNB model
-
namespace fft#
Functions
-
inline void rfft(const BoutReal *in, int length, dcomplex *out)#
Returns the fft of a real signal using fftw_forward
The fftw_forward returns out_k = sum_{j=0}^(length-1) in_j*exp(-2*pi*j*k*sqrt(-1)/length)
Thus, out_k must be divided by ‘length’ in order for DFT[IDFT[in]] = in where IDFT is the inverse fourier transform. See the the fftw user manual for details.
- Parameters:
in – [in] Pointer to the 1D array to take the fourier transform of
length – [in] Number of points in the input array
out – [out] Pointer to the complex 1D array which is the FFT of in
-
inline void irfft(const dcomplex *in, int length, BoutReal *out)#
Take the inverse fft of signal where the outputs are only reals.
This is done through a call to fftw_plan_dft_c2r_1d which is calling fftw_backwards.
That is out_k = sum_{j=0}^(length-1) in_j*exp(2*pi*j*k*sqrt(-1)/length)
See the the fftw user manual for details.
- Parameters:
in – [in] Pointer to the 1D array to take the inverse fourier transform of
length – [in] Number of points in the input array
out – [out] Pointer to the complex 1D array which is IFFTed
-
inline void DST(const BoutReal *in, int length, dcomplex *out)#
Discrete Sine Transform
in
andout
arrays must both be of the samelength
-
inline void DST_rev(dcomplex *in, int length, BoutReal *out)#
Inverse Discrete Sine Transform
in
andout
arrays must both be of the samelength
-
void fft_init(bool fft_measure)#
Should the FFT functions find and use an optimised plan?
-
void fft_init(FFT_MEASUREMENT_FLAG fft_flag)#
Should the FFT functions find and use an optimised plan?
-
void fft_init(Options *options = nullptr)#
Should the FFT functions find and use an optimised plan?
If
options
is not nullptr, it should contain a bool called “fftw_measure”. If it is nullptr, use the globalOptions
root
-
Array<dcomplex> rfft(const Array<BoutReal> &in)#
Returns the fft of a real signal
in
using fftw_forward.
-
Array<BoutReal> irfft(const Array<dcomplex> &in, int length)#
Take the inverse fft of signal
in
where the outputs are only reals. Requires thelength
of the original real signallength
is required because input signals to the forward transform of lengthn
andn + 1
both produce ffts of length(n / 2) + 1
— i.e. it’s not possible to recover the length of the original signal from the fft alone.Expects that
in.size() == (length / 2) + 1
-
inline void rfft(const BoutReal *in, int length, dcomplex *out)#
-
namespace fft#