diff --git a/CMakeLists.txt b/CMakeLists.txt
index 17f55ec81591cb2e99ea0c0ed99420e3b8330950..eb67ee611dd7a61f8c9bea64ad247602510f77ef 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -297,6 +297,8 @@ set(cpp_for_lib
${PROJECT_SOURCE_DIR}/cpp/particles/interpolation/spline_n9.cpp
${PROJECT_SOURCE_DIR}/cpp/particles/interpolation/spline_n10.cpp
${PROJECT_SOURCE_DIR}/cpp/particles/interpolation/Lagrange_polys.cpp
+ ${PROJECT_SOURCE_DIR}/cpp/particles/interpolation/field_tinterpolator.cpp
+ ${PROJECT_SOURCE_DIR}/cpp/particles/interpolation/particle_set.cpp
${PROJECT_SOURCE_DIR}/cpp/scope_timer.cpp
${PROJECT_SOURCE_DIR}/cpp/full_code/test_interpolation.cpp
${PROJECT_SOURCE_DIR}/cpp/full_code/NSVEparticles.cpp
@@ -348,6 +350,11 @@ set(hpp_for_lib
${PROJECT_SOURCE_DIR}/cpp/particles/interpolation/spline_n9.hpp
${PROJECT_SOURCE_DIR}/cpp/particles/interpolation/spline_n10.hpp
${PROJECT_SOURCE_DIR}/cpp/particles/interpolation/Lagrange_polys.hpp
+ ${PROJECT_SOURCE_DIR}/cpp/particles/interpolation/abstract_particle_set.hpp
+ ${PROJECT_SOURCE_DIR}/cpp/particles/interpolation/particle_set.hpp
+ ${PROJECT_SOURCE_DIR}/cpp/particles/interpolation/field_tinterpolator.hpp
+ ${PROJECT_SOURCE_DIR}/cpp/particles/abstract_particle_rhs.hpp
+ ${PROJECT_SOURCE_DIR}/cpp/particles/particle_solver.hpp
${PROJECT_SOURCE_DIR}/cpp/scope_timer.hpp
${PROJECT_SOURCE_DIR}/cpp/full_code/test_interpolation.hpp
${PROJECT_SOURCE_DIR}/cpp/full_code/NSVEparticles.hpp
diff --git a/cpp/particles/abstract_particle_rhs.hpp b/cpp/particles/abstract_particle_rhs.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..92f94477eeebea57e569f253677db42a6d673188
--- /dev/null
+++ b/cpp/particles/abstract_particle_rhs.hpp
@@ -0,0 +1,46 @@
+/******************************************************************************
+* *
+* Copyright 2020 Max Planck Institute for Dynamics and Self-Organization *
+* *
+* This file is part of TurTLE. *
+* *
+* TurTLE 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. *
+* *
+* TurTLE 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 TurTLE. If not, see <http://www.gnu.org/licenses/> *
+* *
+* Contact: Cristian.Lalescu@ds.mpg.de *
+* *
+******************************************************************************/
+
+
+
+#ifndef ABSTRACT_PARTICLE_RHS_HPP
+#define ABSTRACT_PARTICLE_RHS_HPP
+
+#include "particles/interpolation/abstract_particle_set.hpp"
+
+// children of this class must provide a method to compute the rhs for a given particle set
+class abstract_particle_rhs
+{
+ public:
+ // destructor
+ virtual ~abstract_particle_rhs(){}
+
+ // important bit
+ virtual int operator()(
+ double t,
+ const abstract_particle_set &pset,
+ particle_rnumber *result) const = 0;
+};
+
+#endif//ABSTRACT_PARTICLE_RHS_HPP
+
diff --git a/cpp/particles/interpolation/abstract_particle_set.hpp b/cpp/particles/interpolation/abstract_particle_set.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..da1ce30b62f4b32b98c4601564109c3fa2e62970
--- /dev/null
+++ b/cpp/particles/interpolation/abstract_particle_set.hpp
@@ -0,0 +1,98 @@
+/******************************************************************************
+* *
+* Copyright 2020 Max Planck Institute for Dynamics and Self-Organization *
+* *
+* This file is part of TurTLE. *
+* *
+* TurTLE 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. *
+* *
+* TurTLE 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 TurTLE. If not, see <http://www.gnu.org/licenses/> *
+* *
+* Contact: Cristian.Lalescu@ds.mpg.de *
+* *
+******************************************************************************/
+
+
+
+#ifndef ABSTRACT_PARTICLE_SET_HPP
+#define ABSTRACT_PARTICLE_SET_HPP
+
+#include <array>
+#include <cassert>
+#include <vector>
+#include "field.hpp"
+
+// the point here is to have a class that stores particle state data,
+// knows how to distribute this data among processes,
+// and knows how to compute interpolations.
+// the processes of MPI distribution and computing interpolations can't really be separated
+// so it makes sense to merge the functionalities.
+class abstract_particle_set
+{
+ protected:
+ using partsize_t = long long int;
+ using particle_rnumber = double;
+
+ public:
+ // destructor
+ virtual ~abstract_particle_set(){}
+
+ // extract particle information
+ virtual const particle_rnumber* getParticleState() const = 0;
+ virtual const partsize_t* getParticleIndex() const = 0;
+
+ virtual std::unique_ptr<particle_rnumber[]> extractFromParticleState(
+ const int firstComponent,
+ const int lastComponent) const = 0;
+
+ virtual partsize_t getLocalNumberOfParticles() const = 0;
+ virtual partsize_t getTotalNumberOfParticles() const = 0;
+ virtual int getStateSize() const = 0;
+
+ // reorganize particles within MPI domain
+ virtual int redistribute() = 0;
+
+ // sample field values at particle location
+ virtual int sample(
+ const field<float,
+ FFTW,
+ ONE> &field_to_sample,
+ particle_rnumber *result) const = 0;
+ virtual int sample(
+ const field<float,
+ FFTW,
+ THREE> &field_to_sample,
+ particle_rnumber *result) const = 0;
+ virtual int sample(
+ const field<float,
+ FFTW,
+ THREExTHREE> &field_to_sample,
+ particle_rnumber *result) const = 0;
+ virtual int sample(
+ const field<double,
+ FFTW,
+ ONE> &field_to_sample,
+ particle_rnumber *result) const = 0;
+ virtual int sample(
+ const field<double,
+ FFTW,
+ THREE> &field_to_sample,
+ particle_rnumber *result) const = 0;
+ virtual int sample(
+ const field<double,
+ FFTW,
+ THREExTHREE> &field_to_sample,
+ particle_rnumber *result) const = 0;
+};
+
+#endif//ABSTRACT_PARTICLE_SET_HPP
+
diff --git a/cpp/particles/interpolation/field_tinterpolator.cpp b/cpp/particles/interpolation/field_tinterpolator.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..968829b98702f82c7f53a9d954e3abada249ac35
--- /dev/null
+++ b/cpp/particles/interpolation/field_tinterpolator.cpp
@@ -0,0 +1,41 @@
+/******************************************************************************
+* *
+* Copyright 2020 Max Planck Institute for Dynamics and Self-Organization *
+* *
+* This file is part of TurTLE. *
+* *
+* TurTLE 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. *
+* *
+* TurTLE 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 TurTLE. If not, see <http://www.gnu.org/licenses/> *
+* *
+* Contact: Cristian.Lalescu@ds.mpg.de *
+* *
+******************************************************************************/
+
+#include "particles/interpolation/field_tinterpolator.hpp"
+
+template class field_tinterpolator<float, FFTW, ONE, NONE>;
+template class field_tinterpolator<float, FFTW, THREE, NONE>;
+template class field_tinterpolator<float, FFTW, THREExTHREE, NONE>;
+
+template class field_tinterpolator<double, FFTW, ONE, NONE>;
+template class field_tinterpolator<double, FFTW, THREE, NONE>;
+template class field_tinterpolator<double, FFTW, THREExTHREE, NONE>;
+
+template class field_tinterpolator<float, FFTW, ONE, LINEAR>;
+template class field_tinterpolator<float, FFTW, THREE, LINEAR>;
+template class field_tinterpolator<float, FFTW, THREExTHREE, LINEAR>;
+
+template class field_tinterpolator<double, FFTW, ONE, LINEAR>;
+template class field_tinterpolator<double, FFTW, THREE, LINEAR>;
+template class field_tinterpolator<double, FFTW, THREExTHREE, LINEAR>;
+
diff --git a/cpp/particles/interpolation/field_tinterpolator.hpp b/cpp/particles/interpolation/field_tinterpolator.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..5d77e3361b48736bc6c93ee18c5f558bbe46ab38
--- /dev/null
+++ b/cpp/particles/interpolation/field_tinterpolator.hpp
@@ -0,0 +1,124 @@
+/******************************************************************************
+* *
+* Copyright 2020 Max Planck Institute for Dynamics and Self-Organization *
+* *
+* This file is part of TurTLE. *
+* *
+* TurTLE 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. *
+* *
+* TurTLE 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 TurTLE. If not, see <http://www.gnu.org/licenses/> *
+* *
+* Contact: Cristian.Lalescu@ds.mpg.de *
+* *
+******************************************************************************/
+
+
+
+#ifndef FIELD_TINTERPOLATOR_HPP
+#define FIELD_TINTERPOLATOR_HPP
+
+#include <array>
+#include <cassert>
+#include <vector>
+#include "field.hpp"
+#include "particles/particles_distr_mpi.hpp"
+#include "particles/particles_output_hdf5.hpp"
+#include "particles/interpolation/abstract_particle_set.hpp"
+
+enum temporal_interpolation_type {NONE, LINEAR, CUBIC};
+
+template <typename rnumber,
+ field_backend be,
+ field_components fc,
+ temporal_interpolation_type tt>
+class field_tinterpolator
+{
+ private:
+ using particle_rnumber = double;
+ using partsize_t = long long int;
+
+ std::array<field<rnumber, be, fc>*, 4> field_list; // list of fields for temporal interpolation
+
+ public:
+ field_tinterpolator(){
+ this->field_list[0] = NULL;
+ this->field_list[1] = NULL;
+ this->field_list[2] = NULL;
+ this->field_list[3] = NULL;
+ }
+
+ ~field_tinterpolator(){}
+
+ int set_field(
+ field<rnumber, be, fc> *field_src = NULL,
+ const int tindex = 0)
+ {
+ switch(tt)
+ {
+ case NONE:
+ assert(tindex == 0);
+ this->field_list[0] = field_src;
+ break;
+ case LINEAR:
+ assert(tindex == 0 || tindex == 1);
+ this->field_list[tindex] = field_src;
+ break;
+ case CUBIC:
+ assert(false);
+ break;
+ }
+ return EXIT_SUCCESS;
+ }
+
+ // t is a fraction between 0 and 1.
+ int operator()(
+ double t,
+ const abstract_particle_set *pset,
+ particle_rnumber *result)
+ {
+ switch(tt)
+ {
+ case NONE:
+ {
+ pset->sample(*(this->field_list[0]), result);
+ break;
+ }
+ case LINEAR:
+ {
+ particle_rnumber *result0, *result1;
+ result0 = new particle_rnumber[pset->getLocalNumberOfParticles()*ncomp(fc)];
+ result1 = new particle_rnumber[pset->getLocalNumberOfParticles()*ncomp(fc)];
+ pset->sample(*(this->field_list[0]), result0);
+ pset->sample(*(this->field_list[1]), result1);
+ for (partsize_t idx_part = 0; idx_part < pset->getLocalNumberOfParticles(); idx_part++)
+ {
+ for (unsigned int cc = 0; cc < ncomp(fc); cc++)
+ {
+ result[idx_part*ncomp(fc) + cc] = (
+ (1-t)*result0[idx_part*ncomp(fc) + cc] +
+ t *result1[idx_part*ncomp(fc) + cc]);
+ }
+ }
+ delete[] result0;
+ delete[] result1;
+ break;
+ }
+ case CUBIC:
+ assert(false);
+ break;
+ }
+ return EXIT_SUCCESS;
+ }
+};
+
+#endif//FIELD_TINTERPOLATOR_HPP
+
diff --git a/cpp/particles/interpolation/particle_set.cpp b/cpp/particles/interpolation/particle_set.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..656cf72479e5874ef69df0225aa63ee1bce24790
--- /dev/null
+++ b/cpp/particles/interpolation/particle_set.cpp
@@ -0,0 +1,31 @@
+/******************************************************************************
+* *
+* Copyright 2020 Max Planck Institute for Dynamics and Self-Organization *
+* *
+* This file is part of TurTLE. *
+* *
+* TurTLE 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. *
+* *
+* TurTLE 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 TurTLE. If not, see <http://www.gnu.org/licenses/> *
+* *
+* Contact: Cristian.Lalescu@ds.mpg.de *
+* *
+******************************************************************************/
+
+#include "particles/interpolation/particle_set.hpp"
+
+template class particle_set<3, 1, 0>;
+template class particle_set<3, 1, 1>;
+template class particle_set<3, 2, 0>;
+template class particle_set<3, 2, 1>;
+template class particle_set<3, 2, 2>;
+
diff --git a/cpp/particles/interpolation/particle_set.hpp b/cpp/particles/interpolation/particle_set.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..55045e091750387a7bff1fe13e8188120adf97bc
--- /dev/null
+++ b/cpp/particles/interpolation/particle_set.hpp
@@ -0,0 +1,196 @@
+/******************************************************************************
+* *
+* Copyright 2020 Max Planck Institute for Dynamics and Self-Organization *
+* *
+* This file is part of TurTLE. *
+* *
+* TurTLE 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. *
+* *
+* TurTLE 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 TurTLE. If not, see <http://www.gnu.org/licenses/> *
+* *
+* Contact: Cristian.Lalescu@ds.mpg.de *
+* *
+******************************************************************************/
+
+
+
+#ifndef PARTICLE_SET_HPP
+#define PARTICLE_SET_HPP
+
+#include <array>
+#include <cassert>
+#include <vector>
+#include "field.hpp"
+#include "particles/particles_distr_mpi.hpp"
+#include "particles/particles_output_hdf5.hpp"
+#include "particles/interpolation/particles_field_computer.hpp"
+#include "particles/interpolation/particles_generic_interp.hpp"
+#include "particles/interpolation/abstract_particle_set.hpp"
+
+template <int state_size,
+ int neighbours,
+ int smoothness>
+class particle_set: public abstract_particle_set
+{
+ private:
+ // custom names for types
+ using interpolator_class = particles_field_computer<partsize_t,
+ particle_rnumber,
+ particles_generic_interp<particle_rnumber, neighbours, smoothness>,
+ neighbours>;
+ using distributor_class = particles_distr_mpi<partsize_t, particle_rnumber>;
+
+ // MPI
+ MPI_Comm mpi_comm;
+
+ // data
+ const std::pair<int,int> current_partition_interval;
+ const int partition_interval_size;
+
+ std::unique_ptr<partsize_t[]> number_particles_per_partition;
+ std::unique_ptr<partsize_t[]> offset_particles_for_partition;
+
+ std::unique_ptr<particle_rnumber[]> local_state;
+ std::unique_ptr<partsize_t[]> local_index;
+ partsize_t local_number_of_particles;
+ partsize_t total_number_of_particles;
+
+ // for I/O it helps to have custom array shape, we store it here
+ std::vector<hsize_t> particle_file_layout;
+
+ // raw objects for distributing data in MPI, and for interpolating fields
+ distributor_class pdistributor;
+ particles_generic_interp<particle_rnumber, neighbours, smoothness> interpolation_formula;
+ interpolator_class pinterpolator;
+
+ public:
+ particle_set(
+ const MPI_Comm in_mpi_comm,
+ const std::array<particle_rnumber, 3> &in_spatial_box_width,
+ const std::array<particle_rnumber, 3> &in_spatial_box_offset,
+ const std::array<particle_rnumber, 3> &in_spatial_partition_width,
+ const std::array<size_t, 3> &field_grid_dim,
+ const std::array<size_t, 3> &in_local_field_offset,
+ const std::array<size_t, 3> &in_local_field_dims):
+ mpi_comm(in_mpi_comm),
+ current_partition_interval(
+ {in_local_field_offset[IDXC_Z],
+ in_local_field_offset[IDXC_Z] + in_local_field_dims[IDXC_Z]}),
+ partition_interval_size(current_partition_interval.second - current_partition_interval.first),
+ pdistributor(
+ in_mpi_comm,
+ current_partition_interval,
+ field_grid_dim),
+ interpolation_formula(),
+ pinterpolator(
+ field_grid_dim,
+ current_partition_interval,
+ interpolation_formula,
+ in_spatial_box_width,
+ in_spatial_box_offset,
+ in_spatial_partition_width)
+ {
+ this->number_particles_per_partition.reset(new partsize_t[partition_interval_size]);
+ this->offset_particles_for_partition.reset(new partsize_t[partition_interval_size+1]);
+ }
+
+ ~particle_set(){}
+
+ const particle_rnumber* getParticleState() const
+ {
+ return this->local_state.get();
+ }
+
+ const partsize_t* getParticleIndex() const
+ {
+ return this->local_index.get();
+ }
+
+ partsize_t getLocalNumberOfParticles() const
+ {
+ return this->local_number_of_particles;
+ }
+
+ partsize_t getTotalNumberOfParticles() const
+ {
+ return this->total_number_of_particles;
+ }
+
+ int setParticleFileLayout(std::vector<hsize_t> input_layout)
+ {
+ this->particle_file_layout.resize(input_layout.size());
+ for (unsigned int i=0; i<this->particle_file_layout.size(); i++)
+ this->particle_file_layout[i] = input_layout[i];
+ return EXIT_SUCCESS;
+ }
+
+ std::vector<hsize_t> getParticleFileLayout(void)
+ {
+ return std::vector<hsize_t>(this->particle_file_layout);
+ }
+
+ std::unique_ptr<particle_rnumber[]> extractParticlesState(
+ const int firstState,
+ const int lastState) const
+ {
+ const int numberOfStates = std::max(0,(std::min(lastState, state_size)-firstState));
+
+ std::unique_ptr<particle_rnumber[]> stateExtract(new particle_rnumber[local_number_of_particles*numberOfStates]);
+
+ for(partsize_t idx_part = 0 ; idx_part < this->local_number_of_particles ; ++idx_part){
+ for(int idxState = 0 ; idxState < numberOfStates ; ++idxState){
+ stateExtract[idx_part*numberOfStates + idxState] = this->local_state[idx_part*state_size + idxState+firstState];
+ }
+ }
+
+ return stateExtract;
+ }
+
+ template <typename field_rnumber,
+ field_backend be,
+ field_components fc>
+ int sample(const field<field_rnumber, be, fc> &field_to_sample,
+ particle_rnumber *result)
+ {
+ this->pdistributor.template compute_distr<interpolator_class,
+ field<field_rnumber, be, fc>,
+ state_size,
+ ncomp(fc)>(
+ this->pinterpolator,
+ field_to_sample,
+ &this->local_number_of_particles,
+ &this->local_state,
+ result,
+ neighbours);
+ return EXIT_SUCCESS;
+ }
+
+ int redistribute()
+ {
+ // todo: ask why the last template parameter is needed
+ this->pdistributor.template redistribute<interpolator_class,
+ state_size,
+ state_size,
+ 1>(
+ this->pinterpolator,
+ this->number_particles_per_partition.get(),
+ &this->local_number_of_particles,
+ &this->local_state,
+ NULL,
+ 0,
+ &this->local_index);
+ return EXIT_SUCCESS;
+ }
+};
+
+#endif//PARTICLE_SET_HPP
+
diff --git a/cpp/particles/interpolation/particles_field_computer.hpp b/cpp/particles/interpolation/particles_field_computer.hpp
index d03915f7861d579e59d83f0c470923959f9aa24d..620b4685fdc140c999a9bce13784f4dc1d87a773 100644
--- a/cpp/particles/interpolation/particles_field_computer.hpp
+++ b/cpp/particles/interpolation/particles_field_computer.hpp
@@ -28,6 +28,7 @@
#include <array>
#include <utility>
+#include <cmath>
#include "scope_timer.hpp"
#include "particles/particles_utils.hpp"
diff --git a/cpp/particles/interpolation/spline.hpp b/cpp/particles/interpolation/spline.hpp
index ef990088566ec10f0bbf10937980705ffeb570dc..746b21ae13c120ec7442f923b4901fca7a04cd28 100644
--- a/cpp/particles/interpolation/spline.hpp
+++ b/cpp/particles/interpolation/spline.hpp
@@ -2,20 +2,20 @@
* *
* Copyright 2017 Max Planck Institute for Dynamics and Self-Organization *
* *
-* This file is part of bfps. *
+* This file is part of TurTLE. *
* *
-* bfps is free software: you can redistribute it and/or modify *
+* TurTLE 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, *
+* TurTLE 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/> *
+* along with TurTLE. If not, see <http://www.gnu.org/licenses/> *
* *
* Contact: Cristian.Lalescu@ds.mpg.de *
* *
@@ -26,15 +26,15 @@
#ifndef SPLINE_HPP
#define SPLINE_HPP
-#include "spline_n1.hpp"
-#include "spline_n2.hpp"
-#include "spline_n3.hpp"
-#include "spline_n4.hpp"
-#include "spline_n5.hpp"
-#include "spline_n6.hpp"
-#include "spline_n7.hpp"
-#include "spline_n8.hpp"
-#include "spline_n9.hpp"
-#include "spline_n10.hpp"
+#include "particles/interpolation/spline_n1.hpp"
+#include "particles/interpolation/spline_n2.hpp"
+#include "particles/interpolation/spline_n3.hpp"
+#include "particles/interpolation/spline_n4.hpp"
+#include "particles/interpolation/spline_n5.hpp"
+#include "particles/interpolation/spline_n6.hpp"
+#include "particles/interpolation/spline_n7.hpp"
+#include "particles/interpolation/spline_n8.hpp"
+#include "particles/interpolation/spline_n9.hpp"
+#include "particles/interpolation/spline_n10.hpp"
#endif
diff --git a/cpp/particles/particle_solver.hpp b/cpp/particles/particle_solver.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..a8ecff83fff82faec97f670b73e426c72880760a
--- /dev/null
+++ b/cpp/particles/particle_solver.hpp
@@ -0,0 +1,61 @@
+/******************************************************************************
+* *
+* Copyright 2020 Max Planck Institute for Dynamics and Self-Organization *
+* *
+* This file is part of TurTLE. *
+* *
+* TurTLE 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. *
+* *
+* TurTLE 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 TurTLE. If not, see <http://www.gnu.org/licenses/> *
+* *
+* Contact: Cristian.Lalescu@ds.mpg.de *
+* *
+******************************************************************************/
+
+
+
+#ifndef PARTICLE_SOLVER_HPP
+#define PARTICLE_SOLVER_HPP
+
+#include <cassert>
+#include "particles/interpolation/abstract_particle_set.hpp"
+#include "particles/abstract_particle_rhs.hpp"
+
+namespace particle_solver
+{
+ template<int state_size>
+ int Euler(
+ double dt,
+ abstract_particle_set &pset,
+ const abstract_particle_rhs &rhs)
+ {
+ auto *rhs_val = new abstract_particle_set::particle_rnumber[pset.getLocalNumberOfParticles()*state_size];
+ auto *xx = new abstract_particle_set::particle_rnumber[pset.getLocalNumberOfParticles()*state_size];
+ rhs(0, pset, rhs_val);
+ for (abstract_particle_set::partsize_t idx_part = 0; idx_part < pset.getLocalNumberOfParticles(); idx_part++)
+ {
+ for (unsigned int cc = 0; cc < state_size; cc++)
+ {
+ xx[idx_part*state_size + cc] = (
+ pset.getParticleState()[idx_part*state_size + cc] +
+ dt * rhs_val[idx_part*state_size + cc]);
+ }
+ }
+ // TODO: replace pset particle state with xx values
+ delete[] rhs_val;
+ delete[] xx;
+ return EXIT_SUCCESS;
+ }
+};
+
+#endif//PARTICLE_SOLVER_HPP
+
diff --git a/cpp/particles/particles_system.hpp b/cpp/particles/particles_system.hpp
index 6ed1a890f2acc4b9bac5bed7e14e85074bd57d48..2be42334ed636b95fcb0c3343fb67dbabf373219 100644
--- a/cpp/particles/particles_system.hpp
+++ b/cpp/particles/particles_system.hpp
@@ -29,6 +29,7 @@
#include <array>
#include <set>
#include <algorithm>
+#include <cmath>
#include "particles/abstract_particles_system.hpp"
#include "particles/abstract_particles_system_with_p2p.hpp"
@@ -85,10 +86,12 @@ class particles_system : public abstract_particles_system_with_p2p<partsize_t, r
particles_inner_computer_class computer_particules_inner;
public:
- particles_system(const std::array<size_t,3>& field_grid_dim, const std::array<real_number,3>& in_spatial_box_width,
+ particles_system(const std::array<size_t,3>& field_grid_dim,
+ const std::array<real_number,3>& in_spatial_box_width,
const std::array<real_number,3>& in_spatial_box_offset,
const std::array<real_number,3>& in_spatial_partition_width,
- const real_number in_my_spatial_low_limit, const real_number in_my_spatial_up_limit,
+ const real_number in_my_spatial_low_limit,
+ const real_number in_my_spatial_up_limit,
const std::array<size_t,3>& in_local_field_dims,
const std::array<size_t,3>& in_local_field_offset,
const field_class& in_field,
@@ -99,19 +102,36 @@ public:
particles_inner_computer_class in_computer_particules_inner,
const int in_current_iteration = 1)
: mpi_com(in_mpi_com),
- current_partition_interval({in_local_field_offset[IDXC_Z], in_local_field_offset[IDXC_Z] + in_local_field_dims[IDXC_Z]}),
+ current_partition_interval({in_local_field_offset[IDXC_Z],
+ in_local_field_offset[IDXC_Z] + in_local_field_dims[IDXC_Z]}),
partition_interval_size(current_partition_interval.second - current_partition_interval.first),
interpolator(),
- particles_distr(in_mpi_com, current_partition_interval,field_grid_dim),
+ particles_distr(in_mpi_com,
+ current_partition_interval,
+ field_grid_dim),
positions_updater(),
- computer(field_grid_dim, current_partition_interval,
- interpolator, in_spatial_box_width, in_spatial_box_offset, in_spatial_partition_width),
+ computer(field_grid_dim,
+ current_partition_interval,
+ interpolator,
+ in_spatial_box_width,
+ in_spatial_box_offset,
+ in_spatial_partition_width),
default_field(in_field),
- spatial_box_width(in_spatial_box_width), spatial_partition_width(in_spatial_partition_width),
- my_spatial_low_limit(in_my_spatial_low_limit), my_spatial_up_limit(in_my_spatial_up_limit),
- my_nb_particles(0), total_nb_particles(in_total_nb_particles), step_idx(in_current_iteration),
- distr_p2p(in_mpi_com, current_partition_interval,field_grid_dim, spatial_box_width, in_spatial_box_offset, in_cutoff),
- computer_p2p(std::move(in_computer_p2p)), computer_particules_inner(std::move(in_computer_particules_inner)){
+ spatial_box_width(in_spatial_box_width),
+ spatial_partition_width(in_spatial_partition_width),
+ my_spatial_low_limit(in_my_spatial_low_limit),
+ my_spatial_up_limit(in_my_spatial_up_limit),
+ my_nb_particles(0),
+ total_nb_particles(in_total_nb_particles),
+ step_idx(in_current_iteration),
+ distr_p2p(in_mpi_com,
+ current_partition_interval,
+ field_grid_dim,
+ spatial_box_width,
+ in_spatial_box_offset,
+ in_cutoff),
+ computer_p2p(std::move(in_computer_p2p)),
+ computer_particules_inner(std::move(in_computer_particules_inner)){
current_my_nb_particles_per_partition.reset(new partsize_t[partition_interval_size]);
current_offset_particles_for_partition.reset(new partsize_t[partition_interval_size+1]);