Merge branch 'update/split_fftw_mpi_many' into develop

bfps/cpp/field.cpp

@@ -77,7 +77,7 @@ field<rnumber, be, fc>::field(
             ptrdiff_t local_n0, local_0_start;
             ptrdiff_t local_n1, local_1_start;
             //tmp_local_size = fftw_mpi_local_size_many_transposed(
-            fftw_mpi_local_size_many_transposed(
+            fftw_interface<rnumber>::mpi_local_size_many_transposed(
                     3, nfftw, ncomp(fc),
                     FFTW_MPI_DEFAULT_BLOCK, FFTW_MPI_DEFAULT_BLOCK, this->comm,
                     &local_n0, &local_0_start,

bfps/cpp/field.hpp

@@ -72,8 +72,8 @@ class field
         field_layout<fc> *clayout, *rlayout, *rmemlayout;
 
         /* FFT plans */
-        typename fftw_interface<rnumber>::plan c2r_plan;
-        typename fftw_interface<rnumber>::plan r2c_plan;
+        typename fftw_interface<rnumber>::many_plan c2r_plan;
+        typename fftw_interface<rnumber>::many_plan r2c_plan;
         unsigned fftw_plan_rigor;
 
         /* HDF5 data types for arrays */

bfps/cpp/field_descriptor.cpp

@@ -62,7 +62,7 @@ field_descriptor<rnumber>::field_descriptor(
     ptrdiff_t local_n0, local_0_start;
     for (int i = 0; i < this->ndims; i++)
         nfftw[i] = n[i];
-    this->local_size = fftw_mpi_local_size_many(
+    this->local_size = fftw_interface<rnumber>::mpi_local_size_many(
                 this->ndims,
                 &nfftw.front(),
                 1,

bfps/cpp/fluid_solver.cpp

@@ -86,10 +86,10 @@ fluid_solver<rnumber>::fluid_solver(
     this->rv[1] = fftw_interface<rnumber>::alloc_real(this->cd->local_size*2);
     this->rv[2] = this->rv[1];
 
-    this->c2r_vorticity = new typename fftw_interface<rnumber>::plan;
-    this->r2c_vorticity = new typename fftw_interface<rnumber>::plan;
-    this->c2r_velocity  = new typename fftw_interface<rnumber>::plan;
-    this->r2c_velocity  = new typename fftw_interface<rnumber>::plan;
+    this->c2r_vorticity = new typename fftw_interface<rnumber>::many_plan;
+    this->r2c_vorticity = new typename fftw_interface<rnumber>::many_plan;
+    this->c2r_velocity  = new typename fftw_interface<rnumber>::many_plan;
+    this->r2c_velocity  = new typename fftw_interface<rnumber>::many_plan;
 
     ptrdiff_t sizes[] = {nz,
                          ny,
@@ -120,10 +120,10 @@ fluid_solver<rnumber>::fluid_solver(
     this->vc2r[0] = this->c2r_vorticity;
     this->vr2c[0] = this->r2c_vorticity;
 
-    this->vc2r[1] = new typename fftw_interface<rnumber>::plan;
-    this->vr2c[1] = new typename fftw_interface<rnumber>::plan;
-    this->vc2r[2] = new typename fftw_interface<rnumber>::plan;
-    this->vr2c[2] = new typename fftw_interface<rnumber>::plan;
+    this->vc2r[1] = new typename fftw_interface<rnumber>::many_plan;
+    this->vr2c[1] = new typename fftw_interface<rnumber>::many_plan;
+    this->vc2r[2] = new typename fftw_interface<rnumber>::many_plan;
+    this->vr2c[2] = new typename fftw_interface<rnumber>::many_plan;
 
     *(this->vc2r[1]) = fftw_interface<rnumber>::mpi_plan_many_dft_c2r(
                 3, sizes, 3, FFTW_MPI_DEFAULT_BLOCK, FFTW_MPI_DEFAULT_BLOCK,

bfps/cpp/fluid_solver.hpp

@@ -55,12 +55,12 @@ class fluid_solver:public fluid_solver_base<rnumber>
         typename fluid_solver_base<rnumber>::cnumber *cu, *cv[4];
 
         /* plans */
-        typename fftw_interface<rnumber>::plan *c2r_vorticity;
-        typename fftw_interface<rnumber>::plan *r2c_vorticity;
-        typename fftw_interface<rnumber>::plan *c2r_velocity;
-        typename fftw_interface<rnumber>::plan *r2c_velocity;
-        typename fftw_interface<rnumber>::plan *uc2r, *ur2c;
-        typename fftw_interface<rnumber>::plan *vr2c[3], *vc2r[3];
+        typename fftw_interface<rnumber>::many_plan *c2r_vorticity;
+        typename fftw_interface<rnumber>::many_plan *r2c_vorticity;
+        typename fftw_interface<rnumber>::many_plan *c2r_velocity;
+        typename fftw_interface<rnumber>::many_plan *r2c_velocity;
+        typename fftw_interface<rnumber>::many_plan *uc2r, *ur2c;
+        typename fftw_interface<rnumber>::many_plan *vr2c[3], *vc2r[3];
 
         /* physical parameters */
         double nu;

bfps/cpp/fluid_solver_base.cpp

@@ -52,7 +52,7 @@ void fluid_solver_base<rnumber>::clean_up_real_space(rnumber *a, int howmany)
 template <class rnumber>
 double fluid_solver_base<rnumber>::autocorrel(cnumber *a)
 {
-    double *spec = fftw_alloc_real(this->nshells*9);
+    double *spec = fftw_interface<double>::alloc_real(this->nshells*9);
     double sum_local;
     this->cospectrum(a, a, spec);
     sum_local = 0.0;
@@ -60,7 +60,7 @@ double fluid_solver_base<rnumber>::autocorrel(cnumber *a)
     {
         sum_local += spec[n*9] + spec[n*9 + 4] + spec[n*9 + 8];
     }
-    fftw_free(spec);
+    fftw_interface<double>::free(spec);
     return sum_local;
 }
 
@@ -427,7 +427,7 @@ template <class rnumber>
 void fluid_solver_base<rnumber>::write_spectrum(const char *fname, cnumber *a, const double k2exponent)
 {
     TIMEZONE("fluid_solver_base::write_spectrum");
-    double *spec = fftw_alloc_real(this->nshells);
+    double *spec = fftw_interface<double>::alloc_real(this->nshells);
     this->cospectrum(a, a, spec, k2exponent);
     if (this->cd->myrank == 0)
     {
@@ -439,7 +439,7 @@ void fluid_solver_base<rnumber>::write_spectrum(const char *fname, cnumber *a, c
         fwrite((void*)spec, sizeof(double), this->nshells, spec_file);
         fclose(spec_file);
     }
-    fftw_free(spec);
+    fftw_interface<double>::free(spec);
 }
 
 /*****************************************************************************/

bfps/cpp/slab_field_particles.cpp

@@ -69,11 +69,11 @@ slab_field_particles<rnumber>::slab_field_particles(
     this->buffer_width = this->interp_neighbours+1;
     this->buffer_size = this->buffer_width*this->fs->rd->slice_size;
     this->array_size = this->nparticles * this->ncomponents;
-    this->state = fftw_alloc_real(this->array_size);
+    this->state = fftw_interface<rnumber>::alloc_real(this->array_size);
     std::fill_n(this->state, this->array_size, 0.0);
     for (int i=0; i < this->integration_steps; i++)
     {
-        this->rhs[i] = fftw_alloc_real(this->array_size);
+        this->rhs[i] = fftw_interface<rnumber>::alloc_real(this->array_size);
         std::fill_n(this->rhs[i], this->array_size, 0.0);
     }
     this->watching = new bool[this->fs->rd->nprocs*nparticles];
@@ -131,10 +131,10 @@ slab_field_particles<rnumber>::~slab_field_particles()
 {
     delete[] this->computing;
     delete[] this->watching;
-    fftw_free(this->state);
+    fftw_interface<rnumber>::free(this->state);
     for (int i=0; i < this->integration_steps; i++)
     {
-        fftw_free(this->rhs[i]);
+        fftw_interface<rnumber>::free(this->rhs[i]);
     }
     delete[] this->lbound;
     delete[] this->ubound;
@@ -193,7 +193,7 @@ void slab_field_particles<rnumber>::synchronize_single_particle_state(int p, dou
 template <class rnumber>
 void slab_field_particles<rnumber>::synchronize()
 {
-    double *tstate = fftw_alloc_real(this->array_size);
+    double *tstate = fftw_interface<double>::alloc_real(this->array_size);
     // first, synchronize state and jump across CPUs
     std::fill_n(tstate, this->array_size, 0.0);
     for (int p=0; p<this->nparticles; p++)
@@ -236,14 +236,14 @@ void slab_field_particles<rnumber>::synchronize()
                     this->fs->rd->comm);
         }
     }
-    fftw_free(tstate);
+    fftw_interface<double>::free(tstate);
     // assignment of particles
     for (int p=0; p<this->nparticles; p++)
     {
         this->computing[p] = this->get_rank(this->state[p*this->ncomponents + 2]);
         //DEBUG_MSG("synchronizing particles, particle %d computing is %d\n", p, this->computing[p]);
     }
-    double *jump = fftw_alloc_real(this->nparticles);
+    double *jump = fftw_interface<double>::alloc_real(this->nparticles);
     this->jump_estimate(jump);
     // now, see who needs to watch
     bool *local_watching = new bool[this->fs->rd->nprocs*this->nparticles];
@@ -255,7 +255,7 @@ void slab_field_particles<rnumber>::synchronize()
             local_watching[this->get_rank(this->state[this->ncomponents*p+2]-jump[p])*this->nparticles+p] = true;
             local_watching[this->get_rank(this->state[this->ncomponents*p+2]+jump[p])*this->nparticles+p] = true;
         }
-    fftw_free(jump);
+    fftw_interface<double>::free(jump);
     MPI_Allreduce(
             local_watching,
             this->watching,
@@ -389,7 +389,7 @@ void slab_field_particles<rnumber>::step()
 template <class rnumber>
 void slab_field_particles<rnumber>::Euler()
 {
-    double *y = fftw_alloc_real(this->array_size);
+    double *y = fftw_interface<double>::alloc_real(this->array_size);
     this->get_rhs(this->state, y);
     for (int p=0; p<this->nparticles; p++) if (this->fs->rd->myrank == this->computing[p])
     {
@@ -399,7 +399,7 @@ void slab_field_particles<rnumber>::Euler()
         //        "particle %d state is %lg %lg %lg\n",
         //        p, this->state[p*this->ncomponents], this->state[p*this->ncomponents+1], this->state[p*this->ncomponents+2]);
     }
-    fftw_free(y);
+    fftw_interface<double>::free(y);
 }