Resolve "mixed-precision"

Closes #146 (closed)

alpine/ChargedParticles.hpp

@@ -25,6 +25,7 @@ constexpr unsigned Dim = 3;
 // some typedefs
 typedef ippl::ParticleSpatialLayout<double,Dim>   PLayout_t;
 typedef ippl::UniformCartesian<double, Dim>        Mesh_t;
+typedef ippl::UniformCartesian<float, Dim>        Mesh_st;
 typedef ippl::FieldLayout<Dim> FieldLayout_t;
 typedef ippl::OrthogonalRecursiveBisection<double, Dim, Mesh_t> ORB;
 
@@ -33,16 +34,17 @@ using size_type = ippl::detail::size_type;
 template<typename T, unsigned Dim>
 using Vector = ippl::Vector<T, Dim>;
 
-template<typename T, unsigned Dim>
-using Field = ippl::Field<T, Dim>;
+template<typename T, unsigned Dim, class M=Mesh_t>
+using Field = ippl::Field<T, Dim, M>;
 
 template<typename T>
 using ParticleAttrib = ippl::ParticleAttrib<T>;
 
 typedef Vector<double, Dim>  Vector_t;
+typedef Vector<float, Dim> Vector_st;
 typedef Field<double, Dim>   Field_t;
-typedef Field<Vector_t, Dim> VField_t;
-typedef ippl::FFTPeriodicPoissonSolver<Vector_t, double, Dim> Solver_t;
+typedef Field<Vector_st, Dim, Mesh_st> VField_t;
+typedef ippl::FFTPeriodicPoissonSolver<Vector_st, double, Dim, Mesh_st> Solver_t;
 
 const double pi = std::acos(-1.0);
 
@@ -56,7 +58,9 @@ void dumpVTK(VField_t& E, int nx, int ny, int nz, int iteration,
     typename VField_t::view_type::host_mirror_type host_view = E.getHostMirror();
 
     std::stringstream fname;
-    fname << "data/ef_";
+    fname << "data/";
+    fname << nx+3;
+    fname << "/ef_";
     fname << std::setw(4) << std::setfill('0') << iteration;
     fname << ".vtk";
 
@@ -162,7 +166,7 @@ public:
 public:
     ParticleAttrib<double>     q; // charge
     typename ippl::ParticleBase<PLayout>::particle_position_type P;  // particle velocity
-    typename ippl::ParticleBase<PLayout>::particle_position_type E;  // electric field at particle position
+    ParticleAttrib<Vector_st> E;  // electric field at particle position
 
 
     /*
@@ -399,32 +403,32 @@ public:
 
         const int nghostE = E_m.getNghost();
         auto Eview = E_m.getView();
-        Vector_t normE;
+        Vector_st normE;
         using mdrange_type = Kokkos::MDRangePolicy<Kokkos::Rank<3>>;
 
         for (unsigned d=0; d<Dim; ++d) {
 
-        double temp = 0.0;
+        float temp = 0.0;
         Kokkos::parallel_reduce("Vector E reduce",
                                 mdrange_type({nghostE, nghostE, nghostE},
                                              {Eview.extent(0) - nghostE,
                                               Eview.extent(1) - nghostE,
                                               Eview.extent(2) - nghostE}),
                                 KOKKOS_LAMBDA(const size_t i, const size_t j,
-                                              const size_t k, double& valL)
+                                              const size_t k, float& valL)
                                 {
-                                    double myVal = std::pow(Eview(i, j, k)[d], 2);
+                                    float myVal = std::pow(Eview(i, j, k)[d], 2);
                                     valL += myVal;
-                                }, Kokkos::Sum<double>(temp));
-            double globaltemp = 0.0;
-            MPI_Reduce(&temp, &globaltemp, 1, MPI_DOUBLE, MPI_SUM, 0, Ippl::getComm());
+                                }, Kokkos::Sum<float>(temp));
+            float globaltemp = 0.0;
+            MPI_Reduce(&temp, &globaltemp, 1, MPI_FLOAT, MPI_SUM, 0, Ippl::getComm());
             normE[d] = std::sqrt(globaltemp);
         }
 
         if (Ippl::Comm->rank() == 0) {
             std::stringstream fname;
             fname << "data/ParticleField_";
-            fname << Ippl::Comm->size();
+            fname << nr_m[0];
             fname << ".csv";
 
             Inform csvout(NULL, fname.str().c_str(), Inform::APPEND);
@@ -580,7 +584,7 @@ public:
         Kokkos::deep_copy(P_host, P.getView());
         std::stringstream pname;
         pname << "data/ParticleIC_";
-        pname << Ippl::Comm->rank();
+        pname <<  Ippl::Comm->rank();
         pname << ".csv";
         Inform pcsvout(NULL, pname.str().c_str(), Inform::OVERWRITE, Ippl::Comm->rank());
         pcsvout.precision(10);