add (broken) new field<> based solver

the idea is to modify the fluid_solver class into something that uses
the field class, to the point where if we want to change the backend to
something other than plain FFTW, the new solver class will not need to
be modified in any way.

the "vorticity_equation" class included in this commit is not
functional, it is only a checkpoint.

bfps/cpp/field.hpp

@@ -224,6 +224,14 @@ class field
                 const std::string dset_name,
                 const hsize_t toffset,
                 const double max_estimate);
+        inline void impose_zero_mode()
+        {
+            if (this->clayout->myrank == this->clayout->rank[0][0] &&
+                this->real_space_representation == false)
+            {
+                std::fill_n(this->data, 2*ncomp(fc), 0.0);
+            }
+        }
 };
 
 template <typename rnumber,

bfps/cpp/vorticity_equation.hpp

+/**********************************************************************
+*                                                                     *
+*  Copyright 2015 Max Planck Institute                                *
+*                 for Dynamics and Self-Organization                  *
+*                                                                     *
+*  This file is part of bfps.                                         *
+*                                                                     *
+*  bfps is free software: you can redistribute it and/or modify       *
+*  it under the terms of the GNU General Public License as published  *
+*  by the Free Software Foundation, either version 3 of the License,  *
+*  or (at your option) any later version.                             *
+*                                                                     *
+*  bfps is distributed in the hope that it will be useful,            *
+*  but WITHOUT ANY WARRANTY; without even the implied warranty of     *
+*  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the      *
+*  GNU General Public License for more details.                       *
+*                                                                     *
+*  You should have received a copy of the GNU General Public License  *
+*  along with bfps.  If not, see <http://www.gnu.org/licenses/>       *
+*                                                                     *
+* Contact: Cristian.Lalescu@ds.mpg.de                                 *
+*                                                                     *
+**********************************************************************/
+
+#include <stdio.h>
+#include <stdlib.h>
+#include <iostream>
+
+#include "field.hpp"
+
+#ifndef VORTICITY_EQUATION
+
+#define VORTICITY_EQUATION
+
+extern int myrank, nprocs;
+
+
+/* container for field descriptor, fields themselves, parameters, etc
+ * using the same big macro idea that they're using in fftw3.h
+ * I feel like I should quote:  Ugh.
+ * */
+
+template <typename rnumber,
+          field_backend be>
+class vorticity_equation
+{
+    public:
+        /* field descriptor, for direct binary I/O */
+        field_descriptor<rnumber> *cd, *rd;
+
+        /* name */
+        char name[256];
+
+        /* fields */
+        field<rnumber, be, THREE> *cvorticity, *cvelocity;
+        field<rnumber, be, THREE> *rvorticity, *rvelocity;
+        kspace<be, SMOOTH> *kk;
+
+
+        /* short names for velocity, and 4 vorticity fields */
+        field<rnumber, be, THREE> *u, *v[4];
+
+        /* physical parameters */
+        double nu;
+        int fmode;         // for Kolmogorov flow
+        double famplitude; // both for Kflow and band forcing
+        double fk0, fk1;   // for band forcing
+        char forcing_type[128];
+
+        /* methods */
+        vorticity_equation(
+                const char *NAME,
+                int nx,
+                int ny,
+                int nz,
+                double DKX = 1.0,
+                double DKY = 1.0,
+                double DKZ = 1.0,
+                unsigned FFTW_PLAN_RIGOR = FFTW_MEASURE);
+        ~vorticity_equation(void);
+
+        void compute_gradient_statistics(
+                rnumber (*__restrict__ vec)[2],
+                double *__restrict__ gradu_moments,
+                double *__restrict__ trS2_Q_R_moments,
+                ptrdiff_t *__restrict__ gradu_histograms,
+                ptrdiff_t *__restrict__ trS2_Q_R_histograms,
+                ptrdiff_t *__restrict__ QR2D_histogram,
+                double trS2_Q_R_max_estimates[3],
+                double gradu_max_estimates[9],
+                const int nbins_1D = 256,
+                const int nbins_2D = 64);
+
+        void compute_vorticity(void);
+        void compute_velocity(rnumber (*__restrict__ vorticity)[2]);
+        void compute_pressure(rnumber (*__restrict__ pressure)[2]);
+        void compute_Eulerian_acceleration(rnumber (*__restrict__ dst)[2]);
+        void compute_Lagrangian_acceleration(rnumber (*__restrict__ dst)[2]);
+        void compute_Lagrangian_acceleration(rnumber *__restrict__ dst);
+        void omega_nonlin(int src);
+        void step(double dt);
+        void impose_zero_modes(void);
+        void add_forcing(rnumber (*__restrict__ acc_field)[2], rnumber (*__restrict__ vort_field)[2], rnumber factor);
+
+        int read(char field, char representation);
+        int write(char field, char representation);
+        int write_rTrS2();
+        int write_renstrophy();
+        int write_rpressure();
+};
+
+#endif//VORTICITY_EQUATION
+

setup.py

@@ -88,7 +88,8 @@ print('This is bfps version ' + VERSION)
 
 
 ### lists of files and MANIFEST.in
-src_file_list = ['field',
+src_file_list = [#'vorticity_equation',
+                 'field',
                  'field_descriptor',
                  'rFFTW_distributed_particles',
                  'distributed_particles',