diff --git a/bfps/NavierStokes.py b/bfps/NavierStokes.py index a5a0abd2350241b2604f22725c30aceb3b44b827..f3fa39d04c76ac5b400b481ecf423095b4135bcd 100644 --- a/bfps/NavierStokes.py +++ b/bfps/NavierStokes.py @@ -408,8 +408,7 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base): update_fields += 'fs->low_pass_Fourier(fs->cvelocity, 3, {0});\n'.format(kcut) update_fields += ('fs->ift_velocity();\n' + 'vel_{0}->read_rFFTW(fs->rvelocity);\n' + - 'fs->compute_Lagrangian_acceleration(acc_{0}->temp);\n' + - 'acc_{0}->read_rFFTW(acc_{0}->temp);\n').format(name) + 'fs->compute_Lagrangian_acceleration(acc_{0}->field);\n').format(name) self.fluid_start += update_fields self.fluid_loop += update_fields return None diff --git a/bfps/cpp/rFFTW_interpolator.cpp b/bfps/cpp/rFFTW_interpolator.cpp index 4bd47c8c00211aecfc520674b63bece0d159b068..be8546f67e0f3968e814b3f5941e325573fc4460 100644 --- a/bfps/cpp/rFFTW_interpolator.cpp +++ b/bfps/cpp/rFFTW_interpolator.cpp @@ -38,11 +38,9 @@ rFFTW_interpolator<rnumber, interp_neighbours>::rFFTW_interpolator( this->field_size = 2*fs->cd->local_size; this->compute_beta = BETA_POLYS; if (sizeof(rnumber) == 4) - this->f0 = (rnumber*)((void*)fftwf_alloc_real(this->field_size)); + this->field = (rnumber*)((void*)fftwf_alloc_real(this->field_size)); else if (sizeof(rnumber) == 8) - this->f0 = (rnumber*)((void*)fftw_alloc_real(this->field_size)); - this->f1 = this->f0; - this->temp = this->f1; + this->field = (rnumber*)((void*)fftw_alloc_real(this->field_size)); // compute dx, dy, dz; this->dx = 4*acos(0) / (fs->dkx*fs->rd->sizes[2]); @@ -84,9 +82,9 @@ template <class rnumber, int interp_neighbours> rFFTW_interpolator<rnumber, interp_neighbours>::~rFFTW_interpolator() { if (sizeof(rnumber) == 4) - fftwf_free((float*)((void*)this->f0)); + fftwf_free((float*)((void*)this->field)); else if (sizeof(rnumber) == 8) - fftw_free((double*)((void*)this->f0)); + fftw_free((double*)((void*)this->field)); delete[] this->lbound; delete[] this->ubound; } @@ -95,11 +93,10 @@ template <class rnumber, int interp_neighbours> int rFFTW_interpolator<rnumber, interp_neighbours>::read_rFFTW(void *void_src) { rnumber *src = (rnumber*)void_src; - rnumber *dst = this->f1; /* do big copy of middle stuff */ std::copy(src, src + this->field_size, - dst); + this->field); return EXIT_SUCCESS; } @@ -130,7 +127,6 @@ template <class rnumber, int interp_neighbours> void rFFTW_interpolator<rnumber, interp_neighbours>::sample( const int nparticles, const int pdimension, - const double t, const double *__restrict__ x, double *__restrict__ y, const int *deriv) @@ -143,7 +139,7 @@ void rFFTW_interpolator<rnumber, interp_neighbours>::sample( this->get_grid_coordinates(nparticles, pdimension, x, xg, xx); /* perform interpolation */ for (int p=0; p<nparticles; p++) - this->operator()(t, xg + p*3, xx + p*3, yy + p*3, deriv); + this->operator()(xg + p*3, xx + p*3, yy + p*3, deriv); MPI_Allreduce( yy, y, @@ -158,7 +154,6 @@ void rFFTW_interpolator<rnumber, interp_neighbours>::sample( template <class rnumber, int interp_neighbours> void rFFTW_interpolator<rnumber, interp_neighbours>::operator()( - const double t, const int *xg, const double *xx, double *dest, @@ -194,9 +189,9 @@ void rFFTW_interpolator<rnumber, interp_neighbours>::operator()( bigiy)*(this->descriptor->sizes[2]+2) + bigix)*3; for (int c=0; c<3; c++) - dest[c] += this->f1[tindex+c]*(bz[iz+interp_neighbours]* - by[iy+interp_neighbours]* - bx[ix+interp_neighbours]); + dest[c] += this->field[tindex+c]*(bz[iz+interp_neighbours]* + by[iy+interp_neighbours]* + bx[ix+interp_neighbours]); } } } diff --git a/bfps/cpp/rFFTW_interpolator.hpp b/bfps/cpp/rFFTW_interpolator.hpp index 1f0d3065ea3d3f25d045ab6355c2534910e2b28d..1b07bcb25d51fc7d884bc342e0842e555e8018f0 100644 --- a/bfps/cpp/rFFTW_interpolator.hpp +++ b/bfps/cpp/rFFTW_interpolator.hpp @@ -55,7 +55,7 @@ class rFFTW_interpolator /* pointers to fields that are to be interpolated * */ - rnumber *f0, *f1, *temp; + rnumber *field; /* physical parameters of field */ double dx, dy, dz; @@ -80,13 +80,11 @@ class rFFTW_interpolator void sample( const int nparticles, const int pdimension, - const double t, const double *__restrict__ x, double *__restrict__ y, const int *deriv = NULL); /* interpolate 1 point */ void operator()( - const double t, const int *__restrict__ xg, const double *__restrict__ xx, double *__restrict__ dest, diff --git a/bfps/cpp/rFFTW_particles.cpp b/bfps/cpp/rFFTW_particles.cpp index 151d057a7958260c50302f08c28b708eb58ab822..4cd79c1259940861c1b0a8f07d147bc5a1356250 100644 --- a/bfps/cpp/rFFTW_particles.cpp +++ b/bfps/cpp/rFFTW_particles.cpp @@ -87,12 +87,12 @@ rFFTW_particles<particle_type, rnumber, interp_neighbours>::~rFFTW_particles() } template <int particle_type, class rnumber, int interp_neighbours> -void rFFTW_particles<particle_type, rnumber, interp_neighbours>::get_rhs(double t, double *x, double *y) +void rFFTW_particles<particle_type, rnumber, interp_neighbours>::get_rhs(double *x, double *y) { switch(particle_type) { case VELOCITY_TRACER: - this->vel->sample(this->nparticles, this->ncomponents, t, x, y); + this->vel->sample(this->nparticles, this->ncomponents, x, y); break; } } @@ -112,7 +112,7 @@ template <int particle_type, class rnumber, int interp_neighbours> void rFFTW_particles<particle_type, rnumber, interp_neighbours>::AdamsBashforth(int nsteps) { ptrdiff_t ii; - this->get_rhs(1, this->state, this->rhs[0]); + this->get_rhs(this->state, this->rhs[0]); switch(nsteps) { case 1: diff --git a/bfps/cpp/rFFTW_particles.hpp b/bfps/cpp/rFFTW_particles.hpp index e2ba8a87b63d7463f021e0ab991f97f161fac4c3..68f8badea0359436d3aaf3f06b02543ee9ff0d89 100644 --- a/bfps/cpp/rFFTW_particles.hpp +++ b/bfps/cpp/rFFTW_particles.hpp @@ -79,11 +79,10 @@ class rFFTW_particles ~rFFTW_particles(); void get_rhs(double *__restrict__ x, double *__restrict__ rhs); - void get_rhs(double t, double *__restrict__ x, double *__restrict__ rhs); inline void sample_vec_field(rFFTW_interpolator<rnumber, interp_neighbours> *field, double *vec_values) { - field->sample(this->nparticles, this->ncomponents, 1.0, this->state, vec_values, NULL); + field->sample(this->nparticles, this->ncomponents, this->state, vec_values, NULL); } /* input/output */