Changed lengths function of UnitCell to yield vector. Small changes a…

…ccording to revision.

src/core/src/energy/global/ewald.rs

@@ -124,8 +124,8 @@ impl Ewald {
         // Because we do a spherical truncation in k space, we have to transform
         // kmax into a spherical cutoff 'radius'
         let lenghts = cell.lengths();
-        let max_lenght = f64::max(f64::max(lenghts.0, lenghts.1), lenghts.2);
-        let min_lenght = f64::min(f64::min(lenghts.0, lenghts.1), lenghts.2);
+        let max_lenght = f64::max(f64::max(lenghts[0], lenghts[1]), lenghts[2]);
+        let min_lenght = f64::min(f64::min(lenghts[0], lenghts[1]), lenghts[2]);
         let k_rc = self.kmax as f64 * (2.0 * PI / max_lenght);
         self.kmax2 = k_rc * k_rc;
 

src/core/src/sim/mc/moves/resize.rs

@@ -31,13 +31,16 @@ impl Resize {
     /// displacement of `delta`.
     pub fn new(pressure: f64, delta: f64) -> Resize {
         assert!(delta > 0.0, "delta must be positive in Resize move");
-        // TODO: set rc_max from the largest cut off in PairPotentials
+        // TODO: set rc_max from the largest cut off in PairPotentials.
+        // For the time being, we set rc_max to zero.
+        // That way, the simulation cell can get smaller than 2*rc_max
+        // which can lead to wrong energies.
         Resize {
             delta: 0.,
             range: Range::new(-delta, delta),
             old_system: System::new(),
             pressure: pressure,
-            rc_max: 100.0, // arbitrarily large value
+            rc_max: 0.0,
         }
     }
 }
@@ -58,8 +61,8 @@ impl MCMove for Resize {
         let new_cell = system.cell().clone();
         // check the radius of the smallest inscribed sphere and compare to the 
         // cut off distance.
-        // abort simulation when box gets too small
-        if new_cell.lengths_as_vec().iter().any(|&d| 0.5 * d < self.rc_max) {
+        // abort simulation when box gets smaller than twice the cutoff radius
+        if new_cell.lengths().iter().any(|&d| 0.5 * d <= self.rc_max) {
             fatal_error!(
                 "Tried to decrease the cell size but \
                 new size conflicts with the cut off radius. \
@@ -82,7 +85,7 @@ impl MCMove for Resize {
         cache.unused();
         // get system energy before change: this is stored in `self.old_system`
         let old_energy = self.old_system.potential_energy();
-        // get
+        // get system energy after preparation
         let new_energy = system.potential_energy();
         let delta_energy = new_energy - old_energy;
 

src/core/src/sys/cells.rs

@@ -158,7 +158,7 @@ impl UnitCell {
     }
 
     /// Get the distances between faces of the unit cell
-    pub fn lengths(&self) -> (f64, f64, f64) {
+    pub fn lengths(&self) -> [f64; 3] {
         assert!(self.shape != CellShape::Infinite);
 
         let (a, b, c) = (self.vect_a(), self.vect_b(), self.vect_c());
@@ -171,24 +171,7 @@ impl UnitCell {
         let y = f64::abs(nb[0]*b[0] + nb[1]*b[1] + nb[2]*b[2]);
         let z = f64::abs(nc[0]*c[0] + nc[1]*c[1] + nc[2]*c[2]);
 
-        return (x, y, z);
-    }
-
-    /// Get a vector of the distances between faces of the unit cell
-    pub fn lengths_as_vec(&self) -> Vec<f64> {
-        assert!(self.shape != CellShape::Infinite);
-
-        let (a, b, c) = (self.vect_a(), self.vect_b(), self.vect_c());
-        // Plans normal vectors
-        let na = (b ^ c).normalized();
-        let nb = (c ^ a).normalized();
-        let nc = (a ^ b).normalized();
-
-        let x = f64::abs(na[0]*a[0] + na[1]*a[1] + na[2]*a[2]);
-        let y = f64::abs(nb[0]*b[0] + nb[1]*b[1] + nb[2]*b[2]);
-        let z = f64::abs(nc[0]*c[0] + nc[1]*c[1] + nc[2]*c[2]);
-
-        vec![x, y, z]
+        [x, y, z]
     }
 
     /// Get the first angle of the cell
@@ -541,13 +524,13 @@ mod tests {
     #[test]
     fn lengths() {
         let ortho = UnitCell::ortho(3.0, 4.0, 5.0);
-        assert_eq!(ortho.lengths(), (3.0, 4.0, 5.0));
+        assert_eq!(ortho.lengths(), [3.0, 4.0, 5.0]);
 
         let triclinic = UnitCell::triclinic(3.0, 4.0, 5.0, 90.0, 90.0, 90.0);
-        assert_eq!(triclinic.lengths(), (3.0, 4.0, 5.0));
+        assert_eq!(triclinic.lengths(), [3.0, 4.0, 5.0]);
 
         let triclinic = UnitCell::triclinic(3.0, 4.0, 5.0, 90.0, 80.0, 100.0);
-        assert_eq!(triclinic.lengths(), (2.908132319388713, 3.9373265973230853, 4.921658246653857));
+        assert_eq!(triclinic.lengths(), [2.908132319388713, 3.9373265973230853, 4.921658246653857]);
     }
 
     #[test]