double buffered real to complex transform (idaholab#401)

Author: dschwen

include/userobjects/FFTBufferBase.h

@@ -8,18 +8,24 @@
 
 #pragma once
 
+#include "ComplexTypes.h"
 #include "ElementUserObject.h"
+#include "FFTData.h"
 
 template <typename T>
 class FFTBufferBase;
 
 #define usingFFTBufferBaseMembers                                                                  \
   using ElementUserObject::_perf_graph;                                                            \
   using FFTBufferBase<T>::_dim;                                                                    \
-  using FFTBufferBase<T>::_grid;                                                                   \
-  using FFTBufferBase<T>::_buffer;                                                                 \
-  using FFTBufferBase<T>::_start;                                                                  \
-  using FFTBufferBase<T>::_stride;                                                                 \
+  using FFTBufferBase<T>::_real_space_grid;                                                        \
+  using FFTBufferBase<T>::_reciprocal_space_grid;                                                  \
+  using FFTBufferBase<T>::_real_space_data;                                                        \
+  using FFTBufferBase<T>::_reciprocal_space_data;                                                  \
+  using FFTBufferBase<T>::_real_space_data_start;                                                  \
+  using FFTBufferBase<T>::_reciprocal_space_data_start;                                            \
+  using FFTBufferBase<T>::_real_space_data_stride;                                                 \
+  using FFTBufferBase<T>::_reciprocal_space_data_stride;                                           \
   using FFTBufferBase<T>::_how_many
 
 /**
@@ -28,6 +34,8 @@ class FFTBufferBase;
 template <typename T>
 class FFTBufferBase : public ElementUserObject
 {
+  using ComplexT = typename ComplexType<T>::type;
+
 public:
   static InputParameters validParams();
 
@@ -43,31 +51,20 @@ class FFTBufferBase : public ElementUserObject
   virtual void backward() = 0;
   ///@}
 
-  ///@{ data access by index
-  const T & operator[](std::size_t i) const { return _buffer[i]; }
-  T & operator[](std::size_t i) { return _buffer[i]; }
-  ///@}
+  ///@{ buffer access
+  FFTData<T> & realSpace() { return _real_space_data; }
+  FFTData<ComplexT> & reciprocalSpace() { return _reciprocal_space_data; }
+  const FFTData<T> & realSpace() const { return _real_space_data; }
+  const FFTData<ComplexT> & reciprocalSpace() const { return _reciprocal_space_data; }
 
-  ///@{ data access by location
+  ///@{ real space data access by location
   const T & operator()(const Point & p) const;
   T & operator()(const Point & p);
   ///@}
 
-  ///@{ convenience math operators
-  FFTBufferBase<T> & operator+=(FFTBufferBase<T> const & rhs);
-  FFTBufferBase<T> & operator-=(FFTBufferBase<T> const & rhs);
-  FFTBufferBase<T> & operator*=(FFTBufferBase<Real> const & rhs);
-  FFTBufferBase<T> & operator/=(FFTBufferBase<Real> const & rhs);
-  FFTBufferBase<T> & operator*=(Real rhs);
-  FFTBufferBase<T> & operator/=(Real rhs);
-  ///@}
-
   /// return the number of grid cells along each dimension without padding
   const std::vector<int> & grid() const { return _grid; }
 
-  /// return the number of proper grid cells without padding
-  const std::size_t & size() const { return _grid_size; }
-
   /// return the buffer dimension
   const unsigned int & dim() const { return _dim; }
 
@@ -79,12 +76,6 @@ class FFTBufferBase : public ElementUserObject
   }
 
 protected:
-  /// get the addres of the first data element of the ith object in the bufefr
-  Real * start(std::size_t i);
-
-  /// get the number of transforms required for type T
-  std::size_t howMany() const;
-
   ///@{ mesh data
   MooseMesh & _mesh;
   unsigned int _dim;
@@ -103,15 +94,17 @@ class FFTBufferBase : public ElementUserObject
   Real _cell_volume;
 
   ///@{ FFT data buffer and unpadded number of grid cells
-  std::vector<T> _buffer;
-  std::size_t _grid_size;
+  FFTData<T> _real_space_data;
+  FFTData<ComplexT> _reciprocal_space_data;
   ///@}
 
   /// pointer to the start of the data
-  Real * _start;
+  Real * _real_space_data_start;
+  Complex * _reciprocal_space_data_start;
 
   /// stride in units of double size
-  std::ptrdiff_t _stride;
+  std::ptrdiff_t _real_space_data_stride;
+  std::ptrdiff_t _reciprocal_space_data_stride;
 
   /// optional moose sister variabe (to obtain IC from)
   std::vector<const VariableValue *> _moose_variable;

src/executioners/SpectralExecutionerBase.C

@@ -57,16 +57,22 @@ SpectralExecutionerBase::execute()
   _fe_problem.advanceState();
 
   // back and forth test
-  auto & c = getFFTBuffer<Real>("c");
-  c.forward();
-  auto & R = getFFTBuffer<RealVectorValue>("R");
-  R.forward();
+  auto & c_buffer = getFFTBuffer<Real>("c");
+  auto c = c_buffer.realSpace();
+  c_buffer.forward();
+
+  auto & R_buffer = getFFTBuffer<RealVectorValue>("R");
+  auto R = R_buffer.realSpace();
+  R_buffer.forward();
 
   // gradient test
-  auto & u = getFFTBuffer<Real>("u");
-  u.forward();
-  auto & grad_u = getFFTBuffer<RealVectorValue>("grad_u");
-  kVectorMultiply(u, grad_u);
+  auto & u_buffer = getFFTBuffer<Real>("u");
+  auto u = u_buffer.realSpace();
+  u_buffer.forward();
+
+  auto & grad_u_buffer = getFFTBuffer<RealVectorValue>("grad_u");
+  auto grad_u = grad_u_buffer.realSpace();
+  kVectorMultiply(u_buffer, grad_u_buffer);
 
   _time_step = 1;
   _fe_problem.execute(EXEC_FINAL);
@@ -75,14 +81,15 @@ SpectralExecutionerBase::execute()
   _fe_problem.advanceState();
 
   // back and forth test
-  c.backward();
-  R.backward();
+  c_buffer.backward();
+  R_buffer.backward();
   R /= 10000.0;
   c /= 10000.0;
 
   // gradient test
-  u.backward();
-  grad_u.backward();
+  u_buffer.backward();
+  grad_u_buffer.backward();
+
   u /= 10000.0;
   grad_u /= 100.0;
 
@@ -93,23 +100,24 @@ SpectralExecutionerBase::execute()
 }
 
 void
-SpectralExecutionerBase::kVectorMultiply(const FFTBufferBase<Real> & in,
-                                         FFTBufferBase<RealVectorValue> & out) const
+SpectralExecutionerBase::kVectorMultiply(const FFTBufferBase<Real> & in_buffer,
+                                         FFTBufferBase<RealVectorValue> & out_buffer) const
 {
+  mooseAssert(in_buffer.dim() == out_buffer.dim(), "Buffer dimensions must be equal");
+
+  const FFTData<Complex> & in = in_buffer.reciprocalSpace();
+  FFTData<ComplexVectorValue> & out = out_buffer.reciprocalSpace();
   mooseAssert(in.size() == out.size(), "Buffer sizes must be equal");
-  mooseAssert(in.dim() == out.dim(), "Buffer dimensions must be equal");
 
-  const auto & grid = in.grid();
-  switch (in.dim())
+  const auto & grid = in_buffer.grid();
+  switch (in_buffer.dim())
   {
     case 1:
     {
       const int ni = grid[0];
       for (int i = 0; i < ni; ++i)
       {
         out[i](0) = in[i] * i;
-        out[i](1) = 0.0;
-        out[i](2) = 0.0;
       }
       return;
     }
@@ -124,7 +132,6 @@ SpectralExecutionerBase::kVectorMultiply(const FFTBufferBase<Real> & in,
         {
           out[index](0) = in[index] * i;
           out[index](1) = in[index] * j;
-          out[index](2) = 0.0;
           index++;
         }
       return;