cpu/FFT.h

#ifndef PRDC_CPU_FFT_H
#define PRDC_CPU_FFT_H
 
/*
* PSCF - Polymer Self-Consistent Field 
*
* Copyright 2015 - 2025, The Regents of the University of Minnesota
* Distributed under the terms of the GNU General Public License.
*/
 
#include <prdc/cpu/RField.h>
#include <prdc/cpu/RFieldDft.h>
#include <prdc/cpu/CField.h>
#include <pscf/math/IntVec.h>
#include <util/global.h>
 
#include <fftw3.h>
 
namespace Pscf {
namespace Prdc {
namespace Cpu {
 
   using namespace Util;
   using namespace Pscf;

   template <int D>
   class FFT 
   {
 
   public:

      FFT();

      virtual ~FFT();

      void setup(IntVec<D> const & meshDimensions);
 
      // Real Data (Real <-> Complex Transforms)

      void forwardTransform(RField<D> const & in, RFieldDft<D>& out) const;

      void inverseTransformUnsafe(RFieldDft<D>& in, RField<D>& out) const;

      void inverseTransformSafe(RFieldDft<D> const & in, RField<D>& out) 
      const;
 
      // Complex Data (Complex <-> Complex Transforms)

      void forwardTransform(CField<D> const & in, CField<D>& out) const;

      void inverseTransform(CField<D> const & in, CField<D>& out) const;
 
      // Accessors

      IntVec<D> const & meshDimensions() const;

      bool isSetup() const;
 
      // Static function

      static 
      void computeKMesh(IntVec<D> const & rMeshDimensions,
                        IntVec<D> & kMeshDimensions,
                        int & kSize );

      static 
      bool hasImplicitInverse(IntVec<D> const & wavevector,
                              IntVec<D> const & meshDimensions);
 
   private:

      mutable RFieldDft<D> kFieldCopy_;

      IntVec<D> meshDimensions_;

      int rSize_;

      int kSize_;

      fftw_plan rcfPlan_;

      fftw_plan criPlan_;

      fftw_plan ccfPlan_;

      fftw_plan cciPlan_;

      bool isSetup_;

      void makePlans(RField<D>& rField, RFieldDft<D>& kField,
                     CField<D>& cFieldIn, CField<D>& cFieldOut);
 
   };
 
   // Declarations of explicit specializations
 
   template <>
   void FFT<1>::makePlans(RField<1>& rField, RFieldDft<1>& kField,
                          CField<1>& cFieldIn, CField<1>& cFieldOut);
 
   template <>
   void FFT<2>::makePlans(RField<2>& rField, RFieldDft<2>& kField,
                          CField<2>& cFieldIn, CField<2>& cFieldOut);
 
   template <>
   void FFT<3>::makePlans(RField<3>& rField, RFieldDft<3>& kField,
                          CField<3>& cFieldIn, CField<3>& cFieldOut);
 
   // Inline member functions
 
   /*
   * Has this object been setup?
   */
   template <int D>
   inline bool FFT<D>::isSetup() const
   {  return isSetup_; }
 
   /*
   * Return the dimensions of the grid for which this was allocated.
   */
   template <int D>
   inline IntVec<D> const & FFT<D>::meshDimensions() const
   {  return meshDimensions_; }
 
   /*
   * Does this wavevector have an implicit inverse?
   */
   template <int D>
   inline
   bool FFT<D>::hasImplicitInverse(IntVec<D> const & wavevector,
                                   IntVec<D> const & meshDimensions)
   {
      int i = wavevector[D-1];
      int d = meshDimensions[D-1];
      if ((i != 0) && (d - i > d/2 + 1)) {
         return true;
      } else {
         return false;
      }
   }
 
   #ifndef PRDC_CPU_FFT_TPP
   // Suppress implicit instantiation
   extern template class FFT<1>;
   extern template class FFT<2>;
   extern template class FFT<3>;
   #endif
 
} // namespace Pscf::Prdc::Cpu
} // namespace Pscf::Prdc
} // namespace Pscf
#endif