Switch the particle storage to struct of arrays

This uses the soa_derive crate to automatically derive all the needed
methods from the struct definition.

benches/argon.rs

@@ -45,7 +45,7 @@ fn cache_move_particle(bencher: &mut Bencher) {
     let mut rng = utils::get_rng(654646);
 
     let i: usize = rng.gen_range(0, system.size());
-    let mut delta = system.particle(i).position;
+    let mut delta = system.particles().position[i];
     delta += Vector3D::new(rng.gen(), rng.gen(), rng.gen());
 
     bencher.iter(||{

benches/nacl.rs

@@ -85,7 +85,7 @@ fn cache_move_particle_ewald(bencher: &mut Bencher) {
     let mut rng = utils::get_rng(41201154);
 
     let i: usize = rng.gen_range(0, system.size());
-    let mut delta = system.particle(i).position;
+    let mut delta = system.particles().position[i];
     delta += Vector3D::new(rng.gen(), rng.gen(), rng.gen());
 
     bencher.iter(||{
@@ -103,7 +103,7 @@ fn cache_move_particle_wolf(bencher: &mut Bencher) {
     let mut rng = utils::get_rng(474114);
 
     let i: usize = rng.gen_range(0, system.size());
-    let mut delta = system.particle(i).position;
+    let mut delta = system.particles().position[i];
     delta += Vector3D::new(rng.gen(), rng.gen(), rng.gen());
 
     bencher.iter(||{

benches/propane.rs

@@ -45,9 +45,9 @@ fn cache_move_particles(bencher: &mut Bencher) {
     let mut rng = utils::get_rng(84541545);
 
     let molecule = rng.choose(system.molecules()).unwrap();
+    let indexes = molecule.into_iter();
     let mut delta = vec![];
-    for i in molecule {
-        let position = system.particle(i).position;
+    for position in &system.particles().position[indexes] {
         delta.push(position + Vector3D::new(rng.gen(), rng.gen(), rng.gen()));
     }
 
@@ -64,8 +64,9 @@ fn cache_move_all_rigid_molecules(bencher: &mut Bencher) {
     let mut rng = utils::get_rng(7012121);
     for molecule in system.molecules().to_owned() {
         let delta = Vector3D::new(rng.gen(), rng.gen(), rng.gen());
-        for i in molecule {
-            system.particle_mut(i).position += delta;
+        let indexes = molecule.into_iter();
+        for position in &mut system.particles_mut().position[indexes] {
+            *position += delta;
         }
     }
 

benches/water.rs

@@ -95,9 +95,9 @@ fn cache_move_particles_wolf(bencher: &mut Bencher) {
     let mut rng = utils::get_rng(454548784);
 
     let molecule = rng.choose(system.molecules()).unwrap();
+    let indexes = molecule.into_iter();
     let mut delta = vec![];
-    for i in molecule {
-        let position = system.particle(i).position;
+    for position in &system.particles().position[indexes] {
         delta.push(position + Vector3D::new(rng.gen(), rng.gen(), rng.gen()));
     }
 
@@ -116,9 +116,9 @@ fn cache_move_particles_ewald(bencher: &mut Bencher) {
     let mut rng = utils::get_rng(9886565);
 
     let molecule = rng.choose(system.molecules()).unwrap();
+    let indexes = molecule.into_iter();
     let mut delta = vec![];
-    for i in molecule {
-        let position = system.particle(i).position;
+    for position in &system.particles().position[indexes] {
         delta.push(position + Vector3D::new(rng.gen(), rng.gen(), rng.gen()));
     }
 
@@ -137,8 +137,9 @@ fn cache_move_all_rigid_molecules_wolf(bencher: &mut Bencher) {
     let mut rng = utils::get_rng(3);
     for molecule in system.molecules().to_owned() {
         let delta = Vector3D::new(rng.gen(), rng.gen(), rng.gen());
-        for i in molecule {
-            system.particle_mut(i).position += delta;
+        let indexes = molecule.into_iter();
+        for position in &mut system.particles_mut().position[indexes] {
+            *position += delta;
         }
     }
 
@@ -157,8 +158,9 @@ fn cache_move_all_rigid_molecules_ewald(bencher: &mut Bencher) {
     let mut rng = utils::get_rng(2121);
     for molecule in system.molecules().to_owned() {
         let delta = Vector3D::new(rng.gen(), rng.gen(), rng.gen());
-        for i in molecule {
-            system.particle_mut(i).position += delta;
+        let indexes = molecule.into_iter();
+        for position in &mut system.particles_mut().position[indexes] {
+            *position += delta;
         }
     }
 

examples/binary.rs

@@ -5,7 +5,7 @@
 extern crate lumol;
 extern crate lumol_input as input;
 
-use lumol::sys::{Molecule, Particle, TrajectoryBuilder, UnitCell};
+use lumol::sys::{Molecule, Particle, ParticleVec, TrajectoryBuilder, UnitCell};
 use lumol::sys::{read_molecule, molecule_type};
 use lumol::sim::Simulation;
 use lumol::sim::mc::{MonteCarlo, Translate, Rotate};
@@ -30,7 +30,7 @@ fn main() {
     let co2 = {
         // We can read files to get molecule type
         let (molecule, atoms) = read_molecule("data/CO2.xyz").unwrap();
-        molecule_type(&molecule, &atoms)
+        molecule_type(&molecule, atoms.as_slice())
     };
     let h2o = {
         // Or define a new molecule by hand
@@ -41,7 +41,12 @@ fn main() {
         molecule.add_bond(0, 1);
         molecule.add_bond(1, 2);
 
-        molecule_type(&molecule, &[Particle::new("H"), Particle::new("O"), Particle::new("H")])
+        let mut atoms = ParticleVec::new();
+        atoms.push(Particle::new("H"));
+        atoms.push(Particle::new("O"));
+        atoms.push(Particle::new("H"));
+
+        molecule_type(&molecule, atoms.as_slice())
     };
 
     let mut mc = MonteCarlo::new(units::from(500.0, "K").unwrap());

examples/mc_npt_spce.rs

@@ -44,12 +44,10 @@ fn get_system() -> System {
     system.add_angle_potential("H", "O", "H", angle);
 
     // Set charges
-    let h = 0.42380;
-    let o = - 2.0 * h;
     for particle in system.particles_mut() {
         match particle.name.as_ref() {
-            "H" => particle.charge = h,
-            "O" => particle.charge = o,
+            "H" => *particle.charge = 0.42380,
+            "O" => *particle.charge = -2.0 * 0.42380,
             _ => panic!("Unknown particle name in charge setting."),
         }
     }

src/core/Cargo.toml

@@ -25,6 +25,7 @@ rayon = "0.8"
 thread_local = "0.3"
 caldyn = "0.4"
 itertools = "0.6"
+soa_derive = "0.3"
 
 [dev-dependencies]
 tempfile = "2.1"