23 const std::string FFT_3d_local::class_name =
"FFT_3d_local";
25 #ifdef USE_FACTORY_AUTOREGISTER
27 bool init = FFT_3d_local::register_factory();
32 FFT_3d_local::FFT_3d_local()
50 FFT_3d_local::FFT_3d_local(
const Parameters& params)
64 set_parameters(params);
69 FFT_3d_local::~FFT_3d_local()
76 void FFT_3d_local::set_parameters(
const Parameters& params)
83 std::string direction;
85 set_parameters(direction);
91 void FFT_3d_local::get_parameters(
Parameters& params)
const
93 if (m_direction == FORWARD) {
95 }
else if (m_direction == BACKWARD) {
96 params.
set_string(
"FFT_direction",
"Backward");
106 void FFT_3d_local::set_parameters(
const std::string& direction)
108 if (direction ==
"Forward") {
109 m_direction = FORWARD;
110 }
else if (direction ==
"Backward") {
111 m_direction = BACKWARD;
115 vout.
crucial(m_vl,
"Error at %s: unsupported FFT direction \"%s\"\n", class_name.c_str(), direction.c_str());
122 bool FFT_3d_local::check_ok()
127 vout.
crucial(m_vl,
"%s: incompatible with xyz parallelization.\n", class_name.c_str());
136 void FFT_3d_local::initialize()
139 int thread_ok = fftw_init_threads();
149 void FFT_3d_local::initialize_plan(
const Field& src)
152 vout.
detailed(m_vl,
"%s: plan recycled.\n", class_name.c_str());
155 vout.
detailed(m_vl,
"%s: create plan.\n", class_name.c_str());
170 for (
int i = 0; i < m_ndim; ++i) {
175 m_nv = src.
nin() / 2;
177 m_buf_in = fftw_alloc_complex(m_nv * m_vol);
178 m_buf_out = fftw_alloc_complex(m_nv * m_vol);
180 if ((!m_buf_in) || (!m_buf_out)) {
181 vout.
crucial(m_vl,
"%s: memory allocation failed.\n", class_name.c_str());
185 m_plan_fw = fftw_plan_many_dft(m_ndim, m_nsize, m_nv,
186 m_buf_in, m_nsize, m_nv, 1,
187 m_buf_out, m_nsize, m_nv, 1,
188 FFTW_FORWARD, FFTW_ESTIMATE);
190 m_plan_bw = fftw_plan_many_dft(m_ndim, m_nsize, m_nv,
191 m_buf_in, m_nsize, m_nv, 1,
192 m_buf_out, m_nsize, m_nv, 1,
193 FFTW_BACKWARD, FFTW_ESTIMATE);
195 if ((!m_plan_fw) || (!m_plan_bw)) {
196 vout.
crucial(m_vl,
"%s: create plan failed.\n", class_name.c_str());
203 void FFT_3d_local::clear_plan()
205 if (m_buf_in) fftw_free(m_buf_in);
206 if (m_buf_out) fftw_free(m_buf_out);
208 if (m_plan_fw) fftw_destroy_plan(m_plan_fw);
209 if (m_plan_bw) fftw_destroy_plan(m_plan_bw);
214 void FFT_3d_local::finalize()
223 if (not ((dir == FORWARD) || (dir == BACKWARD))) {
224 vout.
crucial(m_vl,
"%s: unsupported direction. %d\n", class_name.c_str(), dir);
228 initialize_plan(src);
235 size_t count = m_nv * m_vol;
237 for (
int iex = 0; iex < nex; ++iex) {
238 for (
int it = 0; it < nt; ++it) {
239 memcpy(m_buf_in, src.
ptr(0, index.
site(0, 0, 0, it), iex),
sizeof(
double) * count * 2);
241 if (dir == FORWARD) {
242 fftw_execute(m_plan_fw);
243 }
else if (dir == BACKWARD) {
244 fftw_execute(m_plan_bw);
246 vout.
crucial(m_vl,
"%s: unsupported direction. %d\n", class_name.c_str(), dir);
250 memcpy(dst.
ptr(0, index.
site(0, 0, 0, it), iex), m_buf_out,
sizeof(
double) * count * 2);
254 if (dir == BACKWARD) {
255 scal(dst, 1.0 / m_vol);
261 void FFT_3d_local::fft(
Field& dst,
const Field& src)
263 return fft(dst, src, m_direction);
268 void FFT_3d_local::fft(
Field& field)
271 vout.
crucial(m_vl,
"Error at %s: fft on-the-fly unsupported.\n", class_name.c_str());