diff --git a/cpp/particles/particles_generic_interp.hpp b/cpp/particles/particles_generic_interp.hpp
index 01f90569d1a9a61d6a8bb0572fef64e9510004f5..fbcc714896131ec7519815e8d6650ae00b69c330 100644
--- a/cpp/particles/particles_generic_interp.hpp
+++ b/cpp/particles/particles_generic_interp.hpp
@@ -26,7 +26,7 @@
 #ifndef PARTICLES_GENERIC_INTERP_HPP
 #define PARTICLES_GENERIC_INTERP_HPP
 
-template <class real_number, int interp_neighbours, int mode>
+template <class real_number, int interp_neighbours, int smoothness>
 class particles_generic_interp;
 
 #include "Lagrange_polys.hpp"
@@ -350,6 +350,28 @@ public:
 /*****************************************************************************/
 /* spline C3 */
 
+template <>
+class particles_generic_interp<double, 1,3>{
+    // no such thing as C3 interpolant with just a 4-point kernel, so we just use C2
+    // this specialisation must exist, otherwise the template_double_for_if will fail
+public:
+    using real_number = double;
+
+    void compute_beta(const int in_derivative, const double in_part_val, double poly_val[]) const {
+        beta_n1_m2(in_derivative, in_part_val, poly_val);
+    }
+};
+
+template <>
+class particles_generic_interp<double, 2,3>{
+public:
+    using real_number = double;
+
+    void compute_beta(const int in_derivative, const double in_part_val, double poly_val[]) const {
+        beta_n2_m3(in_derivative, in_part_val, poly_val);
+    }
+};
+
 template <>
 class particles_generic_interp<double, 3,3>{
 public:
@@ -435,6 +457,48 @@ public:
 /*****************************************************************************/
 /* spline C4 */
 
+template <>
+class particles_generic_interp<double, 1,4>{
+    // no such thing as C4 interpolant with just a 4-point kernel, so we just use C2
+    // this specialisation must exist, otherwise the template_double_for_if will fail
+public:
+    using real_number = double;
+
+    void compute_beta(const int in_derivative, const double in_part_val, double poly_val[]) const {
+        beta_n1_m2(in_derivative, in_part_val, poly_val);
+    }
+};
+
+template <>
+class particles_generic_interp<double, 2,4>{
+public:
+    using real_number = double;
+
+    void compute_beta(const int in_derivative, const double in_part_val, double poly_val[]) const {
+        beta_n2_m4(in_derivative, in_part_val, poly_val);
+    }
+};
+
+template <>
+class particles_generic_interp<double, 3,4>{
+public:
+    using real_number = double;
+
+    void compute_beta(const int in_derivative, const double in_part_val, double poly_val[]) const {
+        beta_n3_m4(in_derivative, in_part_val, poly_val);
+    }
+};
+
+template <>
+class particles_generic_interp<double, 4,4>{
+public:
+    using real_number = double;
+
+    void compute_beta(const int in_derivative, const double in_part_val, double poly_val[]) const {
+        beta_n4_m4(in_derivative, in_part_val, poly_val);
+    }
+};
+
 template <>
 class particles_generic_interp<double, 5,4>{
 public:
diff --git a/cpp/particles/particles_system_builder.hpp b/cpp/particles/particles_system_builder.hpp
index dfc20dc6cac10ffcf98faf682026ea2e4cd745ee..9cc8178e8d743e4374bbb5422c5c27761af24c89 100644
--- a/cpp/particles/particles_system_builder.hpp
+++ b/cpp/particles/particles_system_builder.hpp
@@ -291,7 +291,7 @@ inline std::unique_ptr<abstract_particles_system<partsize_t, particles_rnumber>>
         const int in_current_iteration){
     return Template_double_for_if::evaluate<std::unique_ptr<abstract_particles_system<partsize_t, particles_rnumber>>,
                        int, 1, 11, 1, // interpolation_size
-                       int, 0, 3, 1, // spline_mode
+                       int, 0, 5, 1, // spline_mode
                        particles_system_build_container<partsize_t, field_rnumber,be,fc,particles_rnumber,
                                                         p2p_computer_empty<particles_rnumber,partsize_t>,
                                                         particles_inner_computer_empty<particles_rnumber,partsize_t>,
@@ -321,7 +321,7 @@ inline std::unique_ptr<abstract_particles_system_with_p2p<partsize_t, particles_
         const particles_rnumber cutoff){
     return Template_double_for_if::evaluate<std::unique_ptr<abstract_particles_system_with_p2p<partsize_t, particles_rnumber, p2p_computer_class>>,
                        int, 1, 11, 1, // interpolation_size
-                       int, 0, 3, 1, // spline_mode
+                       int, 0, 5, 1, // spline_mode
                        particles_system_build_container<partsize_t, field_rnumber,be,fc,particles_rnumber,
                                                         p2p_computer_class,
                                                         particles_inner_computer_class,