TBTK
FourierTransform.h
Go to the documentation of this file.
1 /* Copyright 2016 Kristofer Bj√∂rnson
2  *
3  * Licensed under the Apache License, Version 2.0 (the "License");
4  * you may not use this file except in compliance with the License.
5  * You may obtain a copy of the License at
6  *
7  * http://www.apache.org/licenses/LICENSE-2.0
8  *
9  * Unless required by applicable law or agreed to in writing, software
10  * distributed under the License is distributed on an "AS IS" BASIS,
11  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
12  * See the License for the specific language governing permissions and
13  * limitations under the License.
14  */
15 
23 #ifndef COM_DAFER45_TBTK_FOURIER_TRANSFORM
24 #define COM_DAFER45_TBTK_FOURIER_TRANSFORM
25 
26 #include "TBTK/Index.h"
27 
28 #include <fftw3.h>
29 
30 #include <complex>
31 
32 namespace TBTK{
33 
35 public:
37  template<typename DataType>
38  class Plan{
39  public:
41  Plan(
42  DataType *in,
43  DataType *out,
44  int sizeX,
45  int sign
46  );
47 
49  Plan(
50  DataType *in,
51  DataType *out,
52  int sizeX,
53  int sizeY,
54  int sign
55  );
56 
58  Plan(
59  DataType *in,
60  DataType *out,
61  int sizeX,
62  int sizeY,
63  int sizeZ,
64  int sign
65  );
66 
68  Plan(const Plan &plan) = delete;
69 
71  Plan(Plan &&plan);
72 
74  ~Plan();
75 
77  Plan& operator=(const Plan &plan) = delete;
78 
80  Plan& operator=(Plan &&plan);
81 
83  void setNormalizationFactor(double normalizationFactor);
84 
86  double getNormalizationFactor() const;
87  private:
89  fftw_plan *plan;
90 
92  double normalizationFactor;
93 
95  unsigned int size;
96 
98  DataType *input;
99 
101  DataType *output;
102 
104  fftw_plan& getFFTWPlan();
105 
107  unsigned int getSize() const;
108 
110  DataType* getInput();
111 
113  DataType* getOutput();
114 
116  friend class FourierTransform;
117  };
118 
120  template<typename DataType>
121  class ForwardPlan : public Plan<DataType>{
122  public:
125  DataType *in,
126  DataType *out,
127  int sizeX
128  ) : Plan<DataType>(
129  in,
130  out,
131  sizeX,
132  -1
133  ){}
134 
137  DataType *in,
138  DataType *out,
139  int sizeX,
140  int sizeY
141  ) : Plan<DataType>(
142  in,
143  out,
144  sizeX,
145  sizeY,
146  -1
147  ){}
148 
151  DataType *in,
152  DataType *out,
153  int sizeX,
154  int sizeY,
155  int sizeZ
156  ) : Plan<DataType>(
157  in,
158  out,
159  sizeX,
160  sizeY,
161  sizeZ,
162  -1
163  ){}
164  };
165 
167  template<typename DataType>
168  class InversePlan : public Plan<DataType>{
169  public:
172  DataType *in,
173  DataType *out,
174  int sizeX
175  ) : Plan<DataType>(
176  in,
177  out,
178  sizeX,
179  1
180  ){}
181 
184  DataType *in,
185  DataType *out,
186  int sizeX,
187  int sizeY
188  ) : Plan<DataType>(
189  in,
190  out,
191  sizeX,
192  sizeY,
193  1
194  ){}
195 
198  DataType *in,
199  DataType *out,
200  int sizeX,
201  int sizeY,
202  int sizeZ
203  ) : Plan<DataType>(
204  in,
205  out,
206  sizeX,
207  sizeY,
208  sizeZ,
209  1
210  ){}
211  };
212 
214  static void transform(
215  std::complex<double> *in,
216  std::complex<double> *out,
217  int sizeX,
218  int sign
219  );
220 
222  static void transform(
223  std::complex<double> *in,
224  std::complex<double> *out,
225  int sizeX,
226  int sizeY,
227  int sign
228  );
229 
231  static void transform(
232  std::complex<double> *in,
233  std::complex<double> *out,
234  int sizeX,
235  int sizeY,
236  int sizeZ,
237  int sign
238  );
239 
241  template<typename DataType>
242  static void transform(Plan<DataType> &plan);
243 
245  static void forward(
246  std::complex<double> *in,
247  std::complex<double> *out,
248  int sizeX
249  );
250 
252  static void forward(
253  std::complex<double> *in,
254  std::complex<double> *out,
255  int sizeX,
256  int sizeY
257  );
258 
260  static void forward(
261  std::complex<double> *in,
262  std::complex<double> *out,
263  int sizeX,
264  int sizeY,
265  int sizeZ
266  );
267 
269  static void inverse(
270  std::complex<double> *in,
271  std::complex<double> *out,
272  int sizeX
273  );
274 
276  static void inverse(
277  std::complex<double> *in,
278  std::complex<double> *out,
279  int sizeX,
280  int sizeY
281  );
282 
284  static void inverse(
285  std::complex<double> *in,
286  std::complex<double> *out,
287  int sizeX,
288  int sizeY,
289  int sizeZ
290  );
291 private:
292 };
293 
294 template<typename DataType>
296  fftw_execute(plan.getFFTWPlan());
297 
298  double normalizationFactor = plan.getNormalizationFactor();
299  if(normalizationFactor != 1.){
300  Streams::out << "Normalizing\n";
301  DataType *output = plan.getOutput();
302  for(unsigned int n = 0; n < plan.getSize(); n++)
303  output[n] /= normalizationFactor;
304  }
305 }
306 
308  std::complex<double> *in,
309  std::complex<double> *out,
310  int sizeX
311 ){
312  transform(in, out, sizeX, -1);
313 }
314 
316  std::complex<double> *in,
317  std::complex<double> *out,
318  int sizeX,
319  int sizeY
320 ){
321  transform(in, out, sizeX, sizeY, -1);
322 }
323 
325  std::complex<double> *in,
326  std::complex<double> *out,
327  int sizeX,
328  int sizeY,
329  int sizeZ
330 ){
331  transform(in, out, sizeX, sizeY, sizeZ, -1);
332 }
333 
335  std::complex<double> *in,
336  std::complex<double> *out,
337  int sizeX
338 ){
339  transform(in, out, sizeX, 1);
340 }
341 
343  std::complex<double> *in,
344  std::complex<double> *out,
345  int sizeX,
346  int sizeY
347 ){
348  transform(in, out, sizeX, sizeY, 1);
349 }
350 
352  std::complex<double> *in,
353  std::complex<double> *out,
354  int sizeX,
355  int sizeY,
356  int sizeZ
357 ){
358  transform(in, out, sizeX, sizeY, sizeZ, 1);
359 }
360 
361 template<typename DataType>
363  this->plan = plan.plan;
364  plan.plan = nullptr;
365 
366  normalizationFactor = plan.normalizationFactor;
367  size = plan.size;
368  input = plan.input;
369  output = plan.output;
370 }
371 
372 template<typename DataType>
374  if(plan != nullptr){
375  #pragma omp critical (TBTK_FOURIER_TRANSFORM)
376  fftw_destroy_plan(*plan);
377 
378  delete plan;
379  }
380 
381 }
382 
383 template<typename DataType>
385  DataType
386 >::operator=(Plan &&rhs){
387  if(this != &rhs){
388  if(this->plan != nullptr){
389  #pragma omp critical (TBTK_FOURIER_TRANSFORM)
390  fftw_destroy_plan(*this->plan);
391 
392  delete this->plan;
393 
394  this->plan = rhs.plan;
395 
396  normalizationFactor = rhs.normalizationFactor;
397  size = rhs.size;
398  input = rhs.input;
399  output = rhs.output;
400  }
401  }
402 
403  return *this;
404 }
405 
406 template<typename DataType>
408  double normalizationFactor
409 ){
410  this->normalizationFactor = normalizationFactor;
411 }
412 
413 template<typename DataType>
415  return normalizationFactor;
416 }
417 
418 template<typename DataType>
420  return *plan;
421 }
422 
423 template<typename DataType>
424 inline unsigned int FourierTransform::Plan<DataType>::getSize() const{
425  return size;
426 }
427 
428 template<typename DataType>
430  return input;
431 }
432 
433 template<typename DataType>
435  return output;
436 }
437 
438 }; //End of namespace TBTK
439 
440 #endif
Plan & operator=(const Plan &plan)=delete
static void forward(std::complex< double > *in, std::complex< double > *out, int sizeX)
Definition: FourierTransform.h:307
Flexible physical index.
~Plan()
Definition: FourierTransform.h:373
Definition: FourierTransform.h:38
ForwardPlan(DataType *in, DataType *out, int sizeX, int sizeY)
Definition: FourierTransform.h:136
static void transform(std::complex< double > *in, std::complex< double > *out, int sizeX, int sign)
InversePlan(DataType *in, DataType *out, int sizeX)
Definition: FourierTransform.h:171
InversePlan(DataType *in, DataType *out, int sizeX, int sizeY)
Definition: FourierTransform.h:183
Plan(DataType *in, DataType *out, int sizeX, int sign)
Definition: FourierTransform.h:168
ForwardPlan(DataType *in, DataType *out, int sizeX)
Definition: FourierTransform.h:124
static std::ostream out
Definition: Streams.h:36
double getNormalizationFactor() const
Definition: FourierTransform.h:414
void setNormalizationFactor(double normalizationFactor)
Definition: FourierTransform.h:407
Definition: FourierTransform.h:121
ForwardPlan(DataType *in, DataType *out, int sizeX, int sizeY, int sizeZ)
Definition: FourierTransform.h:150
Definition: ModelFactory.h:35
Definition: FourierTransform.h:34
InversePlan(DataType *in, DataType *out, int sizeX, int sizeY, int sizeZ)
Definition: FourierTransform.h:197
static void inverse(std::complex< double > *in, std::complex< double > *out, int sizeX)
Definition: FourierTransform.h:334