File fix-broken-fft-check.patch of Package python3-espressomd

Overview Repositories Revisions Requests Users Attributes Meta

File fix-broken-fft-check.patch of Package python3-espressomd (Revision 0b3503db5641f0f4d441be042d6df3ad)

Currently displaying revision 0b3503db5641f0f4d441be042d6df3ad , Show latest

From 12b1c3a77060a0d402984ce88b230744fa31c79d Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Jean-No=C3=ABl=20Grad?= <jgrad@icp.uni-stuttgart.de>
Date: Thu, 1 Sep 2022 13:37:40 +0200
Subject: [PATCH] core: Fix broken FFT check

Properly check the magnitude of complex residuals. Add a parameter
to disable the check for debugging purposes.
---
 src/core/electrostatics/p3m.cpp                    |  8 +++++---
 src/core/electrostatics/p3m.hpp                    |  3 ++-
 src/core/p3m/fft.cpp                               |  2 +-
 .../unit_tests/EspressoSystemStandAlone_test.cpp   |  2 +-
 src/python/espressomd/electrostatics.py            |  7 +++++++
 src/script_interface/electrostatics/CoulombP3M.hpp |  6 ++++--
 .../electrostatics/CoulombP3MGPU.hpp               |  6 ++++--
 testsuite/python/coulomb_interface.py              |  6 ++++--
 testsuite/python/icc_interface.py                  | 14 ++++++++++----
 testsuite/python/save_checkpoint.py                |  1 +
 testsuite/python/test_checkpoint.py                |  2 ++
 11 files changed, 41 insertions(+), 16 deletions(-)

diff --git a/src/core/electrostatics/p3m.cpp b/src/core/electrostatics/p3m.cpp
index d0d059157f1..b06300fed9e 100644
--- a/src/core/electrostatics/p3m.cpp
+++ b/src/core/electrostatics/p3m.cpp
@@ -270,9 +270,11 @@ void CoulombP3M::init() {
 }
 
 CoulombP3M::CoulombP3M(P3MParameters &&parameters, double prefactor,
-                       int tune_timings, bool tune_verbose)
+                       int tune_timings, bool tune_verbose,
+                       bool check_complex_residuals)
     : p3m{std::move(parameters)}, tune_timings{tune_timings},
-      tune_verbose{tune_verbose} {
+      tune_verbose{tune_verbose}, check_complex_residuals{
+                                      check_complex_residuals} {
 
   m_is_tuned = !p3m.params.tuning;
   p3m.params.tuning = false;
@@ -490,7 +492,7 @@ double CoulombP3M::long_range_kernel(bool force_flag, bool energy_flag,
     }
 
     /* Back FFT force component mesh */
-    auto const check_complex = !p3m.params.tuning;
+    auto const check_complex = !p3m.params.tuning and check_complex_residuals;
     for (int d = 0; d < 3; d++) {
       fft_perform_back(p3m.E_mesh[d].data(), check_complex, p3m.fft, comm_cart);
     }
diff --git a/src/core/electrostatics/p3m.hpp b/src/core/electrostatics/p3m.hpp
index 4f682bf7fa8..a8d02b055b4 100644
--- a/src/core/electrostatics/p3m.hpp
+++ b/src/core/electrostatics/p3m.hpp
@@ -89,13 +89,14 @@ struct CoulombP3M : public Coulomb::Actor<CoulombP3M> {
 
   int tune_timings;
   bool tune_verbose;
+  bool check_complex_residuals;
 
 private:
   bool m_is_tuned;
 
 public:
   CoulombP3M(P3MParameters &&parameters, double prefactor, int tune_timings,
-             bool tune_verbose);
+             bool tune_verbose, bool check_complex_residuals);
 
   bool is_tuned() const { return m_is_tuned; }
 
diff --git a/src/core/p3m/fft.cpp b/src/core/p3m/fft.cpp
index 55ff03c66d3..2c1d193155c 100644
--- a/src/core/p3m/fft.cpp
+++ b/src/core/p3m/fft.cpp
@@ -738,7 +738,7 @@ void fft_perform_back(double *data, bool check_complex, fft_data_struct &fft,
   for (int i = 0; i < fft.plan[1].new_size; i++) {
     fft.data_buf[i] = data[2 * i]; /* real value */
     // Vincent:
-    if (check_complex && (data[2 * i + 1] > 1e-5)) {
+    if (check_complex and std::abs(data[2 * i + 1]) > 1e-5) {
       printf("Complex value is not zero (i=%d,data=%g)!!!\n", i,
              data[2 * i + 1]);
       if (i > 100)
diff --git a/src/core/unit_tests/EspressoSystemStandAlone_test.cpp b/src/core/unit_tests/EspressoSystemStandAlone_test.cpp
index 94cbc51a24d..ddf19a283ed 100644
--- a/src/core/unit_tests/EspressoSystemStandAlone_test.cpp
+++ b/src/core/unit_tests/EspressoSystemStandAlone_test.cpp
@@ -115,7 +115,7 @@ static void mpi_set_tuned_p3m_local(double prefactor) {
                            0.654,
                            1e-3};
   auto solver =
-      std::make_shared<CoulombP3M>(std::move(p3m), prefactor, 1, false);
+      std::make_shared<CoulombP3M>(std::move(p3m), prefactor, 1, false, true);
   ::Coulomb::add_actor(solver);
 }
 
diff --git a/src/python/espressomd/electrostatics.py b/src/python/espressomd/electrostatics.py
index 6b91ec6d355..60a12cc4483 100644
--- a/src/python/espressomd/electrostatics.py
+++ b/src/python/espressomd/electrostatics.py
@@ -155,6 +155,7 @@ def default_params(self):
                 "mesh_off": [-1., -1., -1.],
                 "prefactor": 0.,
                 "check_neutrality": True,
+                "check_complex_residuals": True,
                 "tune": True,
                 "timings": 10,
                 "verbose": True}
@@ -231,6 +232,9 @@ class P3M(_P3MBase):
     check_neutrality : :obj:`bool`, optional
         Raise a warning if the system is not electrically neutral when
         set to ``True`` (default).
+    check_complex_residuals: :obj:`bool`, optional
+        Raise a warning if the backward Fourier transform has non-zero
+        complex residuals when set to ``True`` (default).
 
     """
     _so_name = "Coulomb::CoulombP3M"
@@ -281,6 +285,9 @@ class P3MGPU(_P3MBase):
     check_neutrality : :obj:`bool`, optional
         Raise a warning if the system is not electrically neutral when
         set to ``True`` (default).
+    check_complex_residuals: :obj:`bool`, optional
+        Raise a warning if the backward Fourier transform has non-zero
+        complex residuals when set to ``True`` (default).
 
     """
     _so_name = "Coulomb::CoulombP3MGPU"
diff --git a/src/script_interface/electrostatics/CoulombP3M.hpp b/src/script_interface/electrostatics/CoulombP3M.hpp
index 147ee7d0af7..0f167cfa442 100644
--- a/src/script_interface/electrostatics/CoulombP3M.hpp
+++ b/src/script_interface/electrostatics/CoulombP3M.hpp
@@ -69,6 +69,8 @@ class CoulombP3M : public Actor<CoulombP3M, ::CoulombP3M> {
         {"timings", AutoParameter::read_only,
          [this]() { return actor()->tune_timings; }},
         {"tune", AutoParameter::read_only, [this]() { return m_tune; }},
+        {"check_complex_residuals", AutoParameter::read_only,
+         [this]() { return actor()->check_complex_residuals; }},
     });
   }
 
@@ -85,8 +87,8 @@ class CoulombP3M : public Actor<CoulombP3M, ::CoulombP3M> {
                                get_value<double>(params, "accuracy")};
       m_actor = std::make_shared<CoreActorClass>(
           std::move(p3m), get_value<double>(params, "prefactor"),
-          get_value<int>(params, "timings"),
-          get_value<bool>(params, "verbose"));
+          get_value<int>(params, "timings"), get_value<bool>(params, "verbose"),
+          get_value<bool>(params, "check_complex_residuals"));
     });
     set_charge_neutrality_tolerance(params);
   }
diff --git a/src/script_interface/electrostatics/CoulombP3MGPU.hpp b/src/script_interface/electrostatics/CoulombP3MGPU.hpp
index d212a5e671f..3fe11a9de85 100644
--- a/src/script_interface/electrostatics/CoulombP3MGPU.hpp
+++ b/src/script_interface/electrostatics/CoulombP3MGPU.hpp
@@ -70,6 +70,8 @@ class CoulombP3MGPU : public Actor<CoulombP3MGPU, ::CoulombP3MGPU> {
         {"timings", AutoParameter::read_only,
          [this]() { return actor()->tune_timings; }},
         {"tune", AutoParameter::read_only, [this]() { return m_tune; }},
+        {"check_complex_residuals", AutoParameter::read_only,
+         [this]() { return actor()->check_complex_residuals; }},
     });
   }
 
@@ -86,8 +88,8 @@ class CoulombP3MGPU : public Actor<CoulombP3MGPU, ::CoulombP3MGPU> {
                                get_value<double>(params, "accuracy")};
       m_actor = std::make_shared<CoreActorClass>(
           std::move(p3m), get_value<double>(params, "prefactor"),
-          get_value<int>(params, "timings"),
-          get_value<bool>(params, "verbose"));
+          get_value<int>(params, "timings"), get_value<bool>(params, "verbose"),
+          get_value<bool>(params, "check_complex_residuals"));
     });
     m_actor->request_gpu();
     set_charge_neutrality_tolerance(params);
diff --git a/testsuite/python/coulomb_interface.py b/testsuite/python/coulomb_interface.py
index eff174e6e08..a5a1dc6b828 100644
--- a/testsuite/python/coulomb_interface.py
+++ b/testsuite/python/coulomb_interface.py
@@ -59,7 +59,8 @@ def tearDown(self):
             system, espressomd.electrostatics.P3M,
             dict(prefactor=2., epsilon=0., mesh_off=[0.6, 0.7, 0.8], r_cut=1.5,
                  cao=2, mesh=[8, 10, 8], alpha=12., accuracy=0.01, tune=False,
-                 check_neutrality=True, charge_neutrality_tolerance=7e-12))
+                 check_neutrality=True, charge_neutrality_tolerance=7e-12,
+                 check_complex_residuals=False))
         test_p3m_cpu_non_metallic = tests_common.generate_test_for_actor_class(
             system, espressomd.electrostatics.P3M,
             dict(prefactor=2., epsilon=3., mesh_off=[0.6, 0.7, 0.8], r_cut=1.5,
@@ -78,7 +79,8 @@ def tearDown(self):
             system, espressomd.electrostatics.P3MGPU,
             dict(prefactor=2., epsilon=0., mesh_off=[0.6, 0.7, 0.8], r_cut=1.5,
                  cao=2, mesh=[8, 10, 8], alpha=12., accuracy=0.01, tune=False,
-                 check_neutrality=True, charge_neutrality_tolerance=7e-12))
+                 check_neutrality=True, charge_neutrality_tolerance=7e-12,
+                 check_complex_residuals=False))
         test_p3m_gpu_non_metallic = tests_common.generate_test_for_actor_class(
             system, espressomd.electrostatics.P3MGPU,
             dict(prefactor=2., epsilon=3., mesh_off=[0.6, 0.7, 0.8], r_cut=1.5,
diff --git a/testsuite/python/icc_interface.py b/testsuite/python/icc_interface.py
index 2b968d2432b..6eadca2d0e7 100644
--- a/testsuite/python/icc_interface.py
+++ b/testsuite/python/icc_interface.py
@@ -121,15 +121,21 @@ def test_exceptions_small_r_cut(self):
 
     @utx.skipIfMissingFeatures(["P3M"])
     def test_exceptions_large_r_cut(self):
-        icc, (_, p) = self.setup_icc_particles_and_solver(max_iterations=1)
-        p3m = espressomd.electrostatics.P3M(**self.valid_p3m_parameters())
+        icc, (_, p) = self.setup_icc_particles_and_solver(
+            max_iterations=1, convergence=10.)
+        p3m = espressomd.electrostatics.P3M(
+            check_complex_residuals=False, **self.valid_p3m_parameters())
 
         self.system.actors.add(p3m)
         self.system.actors.add(icc)
-
         with self.assertRaisesRegex(Exception, f"Particle with id {p.id} has a charge .+ that is too large for the ICC algorithm"):
-            p.q = 1e9
+            p.q = 1e8
             self.system.integrator.run(0)
+
+        self.system.actors.remove(icc)
+        self.system.part.clear()
+        icc, (_, p) = self.setup_icc_particles_and_solver(max_iterations=1)
+        self.system.actors.add(icc)
         with self.assertRaisesRegex(Exception, "ICC failed to converge in the given number of maximal steps"):
             p.q = 0.
             self.system.integrator.run(0)
diff --git a/testsuite/python/save_checkpoint.py b/testsuite/python/save_checkpoint.py
index 09ebda6f002..55ab1b37b82 100644
--- a/testsuite/python/save_checkpoint.py
+++ b/testsuite/python/save_checkpoint.py
@@ -125,6 +125,7 @@
         cao=1,
         alpha=1.0,
         r_cut=1.0,
+        check_complex_residuals=False,
         timings=15,
         tune=False)
     if 'ELC' in modes:
diff --git a/testsuite/python/test_checkpoint.py b/testsuite/python/test_checkpoint.py
index 17363085c3f..98f9959f68f 100644
--- a/testsuite/python/test_checkpoint.py
+++ b/testsuite/python/test_checkpoint.py
@@ -504,6 +504,7 @@ def test_p3m(self):
         reference = {'prefactor': 1.0, 'accuracy': 0.1, 'mesh': 3 * [10],
                      'cao': 1, 'alpha': 1.0, 'r_cut': 1.0, 'tune': False,
                      'timings': 15, 'check_neutrality': True,
+                     'check_complex_residuals': False,
                      'charge_neutrality_tolerance': 1e-12}
         for key in reference:
             self.assertIn(key, state)
@@ -519,6 +520,7 @@ def test_elc(self):
         p3m_reference = {'prefactor': 1.0, 'accuracy': 0.1, 'mesh': 3 * [10],
                          'cao': 1, 'alpha': 1.0, 'r_cut': 1.0, 'tune': False,
                          'timings': 15, 'check_neutrality': True,
+                         'check_complex_residuals': False,
                          'charge_neutrality_tolerance': 7e-12}
         elc_reference = {'gap_size': 6.0, 'maxPWerror': 0.1,
                          'delta_mid_top': 0.9, 'delta_mid_bot': 0.1,

Places

File fix-broken-fft-check.patch of Package python3-espressomd (Revision 0b3503db5641f0f4d441be042d6df3ad)

Places