diff --git a/bench/fmm-computation.hpp b/bench/fmm-computation.hpp
index 19f922a4f4ac3597a96f27815c0a844c6e996dee..c7ca40f5cb1db1df896d9cac69217815de88ce17 100644
--- a/bench/fmm-computation.hpp
+++ b/bench/fmm-computation.hpp
@@ -37,6 +37,25 @@
 #include <random>
 #include <vector>
 
+namespace local_utils
+{
+    // Primary template handles matrix kernels that are not gaussian kernels
+    template<typename MatrixKernelType>
+    struct is_gaussian_kernel : std::false_type
+    {
+    };
+
+    // Partial specialization recognizes gaussian kernels
+    template<typename MatrixKernelType>
+    struct is_gaussian_kernel<scalfmm::matrix_kernels::gaussian<MatrixKernelType>> : std::true_type
+    {
+    };
+
+    // Helper alias
+    template<typename MatrixKernelType>
+    inline static constexpr bool is_gaussian_kernel_v = is_gaussian_kernel<MatrixKernelType>::value;
+}   // namespace local_utils
+
 namespace local_args
 {
     struct interp_settings : cpp_tools::cl_parser::required_tag
@@ -192,7 +211,22 @@ namespace local_args
             flagged
         };
     };
-
+    struct gauss_coeff
+    {
+        cpp_tools::cl_parser::str_vec flags = {"--gauss-coeff", "-gc"};
+        std::string description = "Coefficient of the gaussian.";
+        using type = double;
+        std::string input_hint = "double"; /*!< The input hint */
+        type def = type(1.0);
+    };
+    struct gauss_regul
+    {
+        cpp_tools::cl_parser::str_vec flags = {"--gauss-regul", "-gr"};
+        std::string description = "Regularisation Coefficient of the gaussian on the diagonal term.";
+        using type = double;
+        std::string input_hint = "double"; /*!< The input hint */
+        type def = type(0.0);
+    };
     template<typename T>
     std::ostream& operator<<(std::ostream& os, const std::vector<T>& vec)
     {
@@ -430,11 +464,24 @@ auto run(ParserType const& parser) -> void
     // construct the fmm operator
     // construct the near field
     near_field_type near_field;
-    // a reference on the matrix_kernel of the near_field
     auto& near_mk = near_field.matrix_kernel();
 
+    far_matrix_kernel_type far_mk{};
+    if constexpr(local_utils::is_gaussian_kernel_v<near_matrix_kernel_type>)
+    {
+        value_type coeff{parser.template get<local_args::gauss_coeff>()};
+        value_type reg{parser.template get<local_args::gauss_regul>()};
+        // if we deal with a gaussian kernel, we set the different parameters
+        // a reference on the matrix_kernel of the near_field
+
+        near_mk.set_coeff(coeff);
+        near_mk.set_epsilon(reg);
+
+        far_mk.set_coeff(coeff);
+        far_mk.set_epsilon(reg);
+    }
     // build the approximation used in the near field
-    interpolator_type interpolator(order, tree_height, box.width(0));
+    interpolator_type interpolator(far_mk, order, tree_height, box.width(0));
     far_field_type far_field(interpolator);
     //  construct the fmm operator
     fmm_operator_type fmm_operator(near_field, far_field);
@@ -444,6 +491,9 @@ auto run(ParserType const& parser) -> void
     if(fmm_computation)
     {
         // FMM algorithm
+        std::cout << cpp_tools::colors::blue;
+        fmm_operator.settings(std::cout);
+        std::cout << cpp_tools::colors::reset;
         for(std::size_t i = 0; i < nb_runs; ++i)
         {
             // std::cout << "\tFMM run " << i << "\n";
@@ -908,11 +958,12 @@ auto get_parser(int argc, char* argv[], Args&&... args)
 auto main([[maybe_unused]] int argc, [[maybe_unused]] char* argv[]) -> int
 {
     // Parameter handling
-    auto parser = get_parser(
-      argc, argv, local_args::size{}, local_args::order{}, local_args::tree_height{}, local_args::group_size{},
-      local_args::nb_runs{}, local_args::direct_computation{}, local_args::fmm_computation{}, local_args::input_file{},
-      local_args::output_file{}, local_args::dimension{}, local_args::use_float{}, local_args::kernel{},
-      local_args::thread_count{}, local_args::interp_settings{}, local_args::operators_to_proceed{});
+    auto parser = get_parser(argc, argv, local_args::size{}, local_args::order{}, local_args::tree_height{},
+                             local_args::group_size{}, local_args::nb_runs{}, local_args::direct_computation{},
+                             local_args::fmm_computation{}, local_args::input_file{}, local_args::output_file{},
+                             local_args::dimension{}, local_args::use_float{}, local_args::kernel{},
+                             local_args::thread_count{}, local_args::interp_settings{},
+                             local_args::operators_to_proceed{}, local_args::gauss_coeff(), local_args::gauss_regul());
 
     select_precision(parser);
 
diff --git a/checks/check_smooth_kernel.cpp b/checks/check_smooth_kernel.cpp
index 5b5cc4f4ab63712844090722ecc4965ad2f1d752..1a3436b6c86abd22260a497a93399a616828b394 100644
--- a/checks/check_smooth_kernel.cpp
+++ b/checks/check_smooth_kernel.cpp
@@ -50,7 +50,22 @@ namespace local_args
         std::string input_hint = "int"; /*!< The input hint */
         type def = 4;
     };
-
+    struct gauss_coeff
+    {
+        cpp_tools::cl_parser::str_vec flags = {"--gauss-coeff", "-gc"};
+        std::string description = "Coefficient of the gaussian.";
+        using type = double;
+        std::string input_hint = "double"; /*!< The input hint */
+        type def = type(1.0);
+    };
+    struct gauss_regul
+    {
+        cpp_tools::cl_parser::str_vec flags = {"--gauss-regul", "-gr"};
+        std::string description = "Regularisation Coefficient of the gaussian on the diagonal term.";
+        using type = double;
+        std::string input_hint = "double"; /*!< The input hint */
+        type def = type(0.0);
+    };
     template<typename T>
     std::ostream& operator<<(std::ostream& os, const std::vector<T>& vec)
     {
@@ -108,6 +123,7 @@ auto run(ParserType const& parser) -> void
     static constexpr std::size_t nb_inputs{matrix_kernel_type::km};
     static constexpr std::size_t nb_outputs{matrix_kernel_type::kn};
 
+    std::cout << cpp_tools::colors::blue << "dimension: " << dimension << std::endl << cpp_tools::colors::reset;
 #ifdef PART_VAR
     using particle_type =
       scalfmm::container::particle<value_type, dimension, value_type, nb_inputs, value_type, nb_outputs, std::int64_t>;
@@ -123,10 +139,9 @@ auto run(ParserType const& parser) -> void
 
     // interpolation types
     using near_field_type = scalfmm::operators::near_field_operator<matrix_kernel_type>;
-    using interpolator_type =
-      scalfmm::interpolation::interpolator<value_type, dimension, matrix_kernel_type,
-                                           //    scalfmm::options::uniform_<scalfmm::options::fft_>>;
-                                           scalfmm::options::chebyshev_<scalfmm::options::low_rank_>>;
+    using interpolator_type = scalfmm::interpolation::interpolator<value_type, dimension, matrix_kernel_type,
+                                                                   scalfmm::options::uniform_<scalfmm::options::fft_>>;
+    //    scalfmm::options::chebyshev_<scalfmm::options::low_rank_>>;
     using far_field_type = scalfmm::operators::far_field_operator<interpolator_type>;
     using fmm_operator_type = scalfmm::operators::fmm_operators<near_field_type, far_field_type>;
 
@@ -184,12 +199,13 @@ auto run(ParserType const& parser) -> void
 
     // print name of the kernel
     std::cout << "- Kernel name: " << mk.name() << std::endl;
-
     // if we deal with a gaussian kernel, we set the different parameters
     if constexpr(local_utils::is_gaussian_kernel_v<matrix_kernel_type>)
     {
-        mk.set_coeff(5.);
-        mk.set_epsilon(1.5);
+        value_type coeff{parser.template get<local_args::gauss_coeff>()};
+        value_type reg{parser.template get<local_args::gauss_regul>()};
+        mk.set_coeff(coeff);
+        mk.set_epsilon(reg);
     }
 
     // test whether the kernel is smooth or not (in that case, just print message to the standard output)
@@ -201,14 +217,14 @@ auto run(ParserType const& parser) -> void
     {
         std::cout << "Kernel is not smooth" << std::endl;
     }
-
+    std::cout << " Good order " << int(nb_particles / std::pow(2, tree_height)) << std::endl;
     // build the approximation used in the near field
     interpolator_type interpolator(mk, order, tree_height, box.width(0));
     far_field_type far_field(interpolator);
     //  construct the fmm operator
     fmm_operator_type fmm_operator(near_field, far_field);
-    using settings = typename interpolator_type::settings;
-    settings s;
+    // using settings = typename interpolator_type::settings;
+    // settings s;
 
     std::cout << cpp_tools::colors::blue;
     fmm_operator.settings(std::cout);
@@ -218,6 +234,8 @@ auto run(ParserType const& parser) -> void
     group_tree_type tree(tree_height, order, box, group_size, group_size, container);
 
     tree.statistics("smmoth kernel", std::cout);
+
+    scalfmm::list::sequential::build_interaction_lists(tree, tree, fmm_operator);
     // now we have everything to call the fmm algorithm
     scalfmm::algorithms::fmm[scalfmm::options::_s(scalfmm::options::seq)](tree, fmm_operator);
 
@@ -280,25 +298,26 @@ auto main([[maybe_unused]] int argc, [[maybe_unused]] char* argv[]) -> int
     const std::size_t nb_particles{10000};
 
     // Parameter handling
-    auto parser = get_parser(argc, argv, local_args::nb_particles{}, local_args::tree_height{}, local_args::order{});
+    auto parser = get_parser(argc, argv, local_args::nb_particles{}, local_args::tree_height{}, local_args::order{},
+                             local_args::gauss_coeff(), local_args::gauss_regul());
 
     // Gaussian kernel
     {
         using matrix_kernel_type = scalfmm::matrix_kernels::gaussian<value_type>;
         run<value_type, dimension, matrix_kernel_type>(parser);
     }
-
+#ifdef ALL_KERNELS
     // One over r (modified)
     {
         using matrix_kernel_type = scalfmm::matrix_kernels::debug::one_over_r_modified;
         run<value_type, dimension, matrix_kernel_type>(parser);
     }
+#endif
 
     // One over r (classical)
     {
         using matrix_kernel_type = scalfmm::matrix_kernels::laplace::one_over_r;
         run<value_type, dimension, matrix_kernel_type>(parser);
     }
-
     return 0;
 }