Allow to set alpha from the input files

doc/src/input/simulations/index.rst

@@ -1,6 +1,6 @@
-*****************************
-The ``[simulations]`` section
-*****************************
+*******************************
+The ``[[simulations]]`` section
+*******************************
 
 The way to propagate the system is defined in the ``[[simulations]]`` section of
 the input file. This section always contains at least two keys: ``nsteps``

doc/src/input/systems/interactions/electrostatic.rst

@@ -59,6 +59,17 @@ Usually 7-8 is a good value for pure water, for a very periodic charges
 distribution (like a crystal) a lower value, such as 5 is sufficient, and for
 more heterogeneous system, higher values of ``kmax`` are needed.
 
+It is also possible (but not required) to set the ``alpha`` parameter of Ewald
+solver directly with the corresponding keyword. A good value of ``alpha`` one
+that satisfies :math:`\exp \left(-\alpha \frac L 2 \right) << 1`. If no value is
+provided in the input file, the default value of :math:`\pi / \text{cutoff}` is
+used.
+
+.. code::
+
+    [coulomb]
+    ewald = {cutoff = "9 A", kmax = 7, alpha = "0.33451"}
+
 Wolf solver
 -----------
 

src/input/src/interactions/toml.rs

@@ -50,10 +50,11 @@ impl FromToml for Mie {
         let epsilon = extract::str("epsilon", table, "Mie potential")?;
         let m = extract::number("m", table, "Mie potential")?;
         let n = extract::number("n", table, "Mie potential")?;
-        
+
         if m < 3.0 {
             warn!("'m' is smaller than 3. Tail corrections for Mie potential are set to zero.");
         };
+
         Ok(Mie::new(
                 units::from_str(sigma)?,
                 units::from_str(epsilon)?,
@@ -166,6 +167,11 @@ impl FromToml for Ewald {
     fn from_toml(table: &Table) -> Result<Ewald> {
         let cutoff = extract::str("cutoff", table, "Ewald coulombic potential")?;
         let kmax = extract::uint("kmax", table, "Ewald coulombic potential")?;
-        Ok(Ewald::new(units::from_str(cutoff)?, kmax as usize, None))
+        let alpha = if table.contains_key("alpha") {
+            Some(extract::number("alpha", table, "Ewald coulombic potential")?)
+        } else {
+            None
+        };
+        Ok(Ewald::new(units::from_str(cutoff)?, kmax as usize, alpha))
     }
 }

src/input/tests/interactions/bad/coulomb.toml

@@ -47,6 +47,14 @@ ewald = {cutoff = "6 A", kmax = -7}
 [input]
 version = 1
 
+[coulomb]
+ewald = {cutoff = "6 A", kmax = 7, alpha = false}
+#^ 'alpha' must be a number in Ewald coulombic potential
+
++++
+[input]
+version = 1
+
 [coulomb]
 ewald = true
 #^ Coulombic solver 'ewald' must be a table

src/input/tests/interactions/good/ewald.toml

@@ -8,3 +8,15 @@ restriction = "exclude13"
 [charges]
 A = -8
 B = 3
+
++++
+
+[input]
+version = 1
+
+[coulomb]
+ewald = {cutoff = "189 A", kmax = 112, alpha = 0.28}
+
+[charges]
+A = -8
+B = 3

tests/data/md-nacl/ewald.toml

@@ -17,7 +17,7 @@ atoms = ["Cl", "Cl"]
 lj = {sigma = "4.612 A", epsilon = "0.02502 kcal/mol"}
 
 [coulomb]
-ewald = {cutoff = "5.5 A", kmax = 10, alpha = ""}
+ewald = {cutoff = "5.5 A", kmax = 10}
 
 [charges]
 Na = 1.0