Cargo.toml

@@ -6,3 +6,5 @@ edition = "2021"
 # See more keys and their definitions at https://doc.rust-lang.org/cargo/reference/manifest.html
 
 [dependencies]
+nalgebra = "0.33.2"
+approx = "0.5"

src/context.rs

@@ -0,0 +1,20 @@
+use crate::integrator::Integrator;
+use crate::platform::Platform;
+use crate::state::{State, StateType};
+use crate::system::System;
+use nalgebra::Vector3;
+
+pub struct Context {
+    pub system: System,
+    pub integrator: Box<dyn Integrator>,
+    pub platform: Box<dyn Platform>,
+}
+impl Context {
+    pub fn set_positions(&mut self, positions: Vec<Vector3<f32>>) {
+        unimplemented!()
+    }
+
+    pub fn get_state(&self, state_type: StateType) -> State {
+        unimplemented!()
+    }
+}

src/force/harmonic_angle_force.rs

@@ -0,0 +1,22 @@
+use crate::context::Context;
+use crate::force::Force;
+
+#[derive(Debug, Clone, PartialEq)]
+pub struct HarmonicAngleForce {}
+impl Force for HarmonicAngleForce {
+    fn add_bond(&mut self, particle1: usize, particle2: usize, length: f64, k: f64) {
+        unimplemented!()
+    }
+
+    fn set_bond_parameters(&mut self, particle1: usize, particle2: usize, length: f64, k: f64) {
+        unimplemented!()
+    }
+
+    fn update_parameters_in_context(&mut self, context: &mut Context) {
+        unimplemented!()
+    }
+
+    fn set_uses_periodic_boundary_conditions(&mut self, status: bool) {
+        unimplemented!()
+    }
+}

src/force/harmonic_bond_force.rs

@@ -0,0 +1,127 @@
+use crate::context::Context;
+use crate::force::Force;
+
+#[derive(Debug, Clone, PartialEq)]
+pub struct HarmonicBondForce {}
+impl Force for HarmonicBondForce {
+    fn add_bond(&mut self, particle1: usize, particle2: usize, length: f64, k: f64) {
+        unimplemented!()
+    }
+
+    fn set_bond_parameters(&mut self, particle1: usize, particle2: usize, length: f64, k: f64) {
+        unimplemented!()
+    }
+
+    fn update_parameters_in_context(&mut self, context: &mut Context) {
+        unimplemented!()
+    }
+
+    fn set_uses_periodic_boundary_conditions(&mut self, status: bool) {
+        unimplemented!()
+    }
+}
+
+#[cfg(test)]
+mod test {
+    use super::HarmonicBondForce;
+    use crate::context::Context;
+    use crate::force::Force;
+    use crate::integrator::VerletIntegrator;
+    use crate::platform::CPU;
+    use crate::state::{State, StateType};
+    use crate::system::System;
+    use approx::assert_relative_eq;
+    use nalgebra::Vector3;
+
+    #[test]
+    fn test_harmonic_bond_force() {
+        let mut system = System {};
+        system.add_particle(1.0);
+        system.add_particle(1.0);
+        system.add_particle(1.0);
+        let integrator = VerletIntegrator { dt: 0.01 };
+
+        let mut forcefield = HarmonicBondForce {};
+        forcefield.add_bond(0, 1, 1.5, 0.8);
+        forcefield.add_bond(1, 2, 1.2, 0.7);
+
+        system.add_force(Box::new(forcefield.clone()));
+
+        let mut context = Context {
+            system,
+            integrator: Box::new(integrator),
+            platform: Box::new(CPU),
+        };
+        let mut positions = vec![Vector3::zeros(); 3];
+        positions[0] = Vector3::new(0.0, 2.0, 0.0);
+        positions[1] = Vector3::new(0.0, 0.0, 0.0);
+        positions[2] = Vector3::new(1.0, 0.0, 0.0);
+        context.set_positions(positions);
+
+        let state = context.get_state(StateType::Both);
+        let forces = state.get_forces();
+        assert_relative_eq!(forces[0], Vector3::new(0.0, -0.8 * 0.5, 0.0));
+        assert_relative_eq!(forces[2], Vector3::new(0.7 * 0.2, 0.0, 0.0));
+        assert_relative_eq!(forces[1], -forces[0] - forces[2]);
+
+        let energy = state.get_potential_energy();
+        assert_relative_eq!(energy, 0.5 * 0.8 * 0.5 * 0.5 + 0.5 * 0.7 * 0.2 * 0.2);
+
+        forcefield.set_bond_parameters(0, 0, 1.6, 0.9);
+        forcefield.set_bond_parameters(1, 1, 1.3, 0.8);
+        forcefield.update_parameters_in_context(&mut context);
+
+        let state = context.get_state(StateType::Both);
+        let forces = state.get_forces();
+        assert_relative_eq!(forces[0], Vector3::new(0.0, -0.9 * 0.4, 0.0));
+        assert_relative_eq!(forces[2], Vector3::new(0.8 * 0.3, 0.0, 0.0));
+        assert_relative_eq!(forces[1], -forces[0] - forces[2]);
+
+        let energy = state.get_potential_energy();
+        assert_relative_eq!(energy, 0.5 * 0.9 * 0.4 * 0.4 + 0.5 * 0.8 * 0.3 * 0.3);
+    }
+
+    #[test]
+    fn test_harmonic_bond_force_with_periodic_boundary_conditions() {
+        let mut system = System {};
+        system.add_particle(1.0);
+        system.add_particle(1.0);
+        system.set_default_periodic_box_vectors(
+            Vector3::new(3.0, 0.0, 0.0),
+            Vector3::new(0.0, 3.0, 0.0),
+            Vector3::new(0.0, 0.0, 3.0),
+        );
+
+        let integrator = VerletIntegrator { dt: 0.01 };
+        let mut forcefield = HarmonicBondForce {};
+        forcefield.add_bond(0, 1, 1.2, 0.8);
+        forcefield.set_uses_periodic_boundary_conditions(true);
+
+        system.add_force(Box::new(forcefield));
+
+        let mut context = Context {
+            system,
+            integrator: Box::new(integrator),
+            platform: Box::new(CPU),
+        };
+
+        let mut positions = vec![Vector3::zeros(); 2];
+        positions[0] = Vector3::new(0.0, 2.0, 0.0);
+        positions[1] = Vector3::new(0.0, 0.0, 0.0);
+        context.set_positions(positions);
+
+        let state = context.get_state(StateType::Both);
+        let forces = state.get_forces();
+        assert_relative_eq!(forces[0], Vector3::new(0.0, -0.8 * 0.2, 0.0));
+        assert_relative_eq!(forces[1], Vector3::new(0.0, 0.8 * 0.2, 0.0));
+
+        let energy = state.get_potential_energy();
+        assert_relative_eq!(energy, 0.5 * 0.8 * 0.2 * 0.2);
+    }
+
+    #[test]
+    fn test_harmonic_bond_force_parallel_computation() {
+        // TODO: Rust-OpenMM
+        unimplemented!()
+    }
+}

src/force/mod.rs

@@ -0,0 +1,16 @@
+use std::fmt::Debug;
+
+use crate::context::Context;
+
+pub trait Force: Debug {
+    fn add_bond(&mut self, particle1: usize, particle2: usize, length: f64, k: f64);
+    fn set_bond_parameters(&mut self, particle1: usize, particle2: usize, length: f64, k: f64);
+    fn update_parameters_in_context(&mut self, context: &mut Context);
+    fn set_uses_periodic_boundary_conditions(&mut self, status: bool);
+}
+
+pub mod harmonic_angle_force;
+pub use harmonic_angle_force::HarmonicAngleForce;
+
+pub mod harmonic_bond_force;
+pub use harmonic_bond_force::HarmonicBondForce;

src/integrator/mod.rs

@@ -0,0 +1,6 @@
+pub trait Integrator {}
+
+pub struct VerletIntegrator {
+    pub dt: f64,
+}
+impl Integrator for VerletIntegrator {}

src/lib.rs

@@ -1,228 +1,8 @@
-pub mod system {
-    use crate::force::Force;
-    use crate::virtual_site::VirtualSite;
+pub mod context;
+pub mod force;
+pub mod platform;
+pub mod state;
+pub mod system;
 
-    pub struct System {}
-    impl System {
-        pub fn add_particle(&mut self, _mass: f32) {
-            // Todo: Rust-OpenMM
-            unimplemented!()
-        }
-
-        pub fn get_num_particles(&self) -> usize {
-            // Todo: Rust-OpenMM
-            unimplemented!()
-        }
-
-        pub fn get_particle_mass(&self, _index: usize) -> f32 {
-            // Todo: Rust-OpenMM
-            unimplemented!()
-        }
-
-        pub fn set_particle_mass(&mut self, _index: usize, _mass: f32) {
-            // Todo: Rust-OpenMM
-            unimplemented!()
-        }
-
-        pub fn is_virtual_site(&self, _index: usize) -> bool {
-            // Todo: Rust-OpenMM
-            unimplemented!()
-        }
-
-        pub fn get_virtual_site(&self, _index: usize) -> Option<Box<dyn VirtualSite>> {
-            // Todo: Rust-OpenMM
-            unimplemented!()
-        }
-
-        pub fn set_virtual_site(&mut self, _index: usize, _vsite: Option<Box<dyn VirtualSite>>) {
-            // Todo: Rust-OpenMM
-            unimplemented!()
-        }
-
-        pub fn add_constraint(&mut self, _p1: usize, _p2: usize, _dist: f32) {
-            // Todo: Rust-OpenMM
-            unimplemented!()
-        }
-
-        pub fn remove_constraint(&mut self, _index: usize) {
-            // Todo: Rust-OpenMM
-            unimplemented!()
-        }
-
-        pub fn get_num_constraints(&self) -> usize {
-            // Todo: Rust-OpenMM
-            unimplemented!()
-        }
-
-        pub fn get_constraint_parameters(&self, _index: usize) -> (usize, usize, f32) {
-            // Todo: Rust-OpenMM
-            unimplemented!()
-        }
-
-        pub fn set_constraint_parameters(
-            &mut self,
-            _index: usize,
-            _p1: usize,
-            _p2: usize,
-            _dist: f32,
-        ) {
-            // Todo: Rust-OpenMM
-            unimplemented!()
-        }
-
-        pub fn add_force(&mut self, _force: Box<dyn Force>) {
-            // Todo: Rust-OpenMM
-            unimplemented!()
-        }
-
-        pub fn remove_force(&mut self, _index: usize) {
-            // Todo: Rust-OpenMM
-            unimplemented!()
-        }
-
-        pub fn get_num_forces(&self) -> usize {
-            // Todo: Rust-OpenMM
-            unimplemented!()
-        }
-
-        pub fn get_force(&self, _index: usize) -> Box<dyn Force> {
-            // Todo: Rust-OpenMM
-            unimplemented!()
-        }
-    }
-}
-
-pub mod force {
-    use std::fmt::Debug;
-
-    pub trait Force: Debug {}
-
-    #[derive(Debug, Clone, PartialEq)]
-    pub struct HarmonicAngleForce {}
-    impl Force for HarmonicAngleForce {}
-
-    #[derive(Debug, Clone, PartialEq)]
-    pub struct HarmonicBondForce {}
-    impl Force for HarmonicBondForce {}
-}
-
-pub mod virtual_site {
-    use std::fmt::Debug;
-
-    pub trait VirtualSite: Debug {}
-
-    #[derive(Debug, Clone, PartialEq)]
-    pub struct TwoParticleAverageSite {
-        pub p1: usize,
-        pub p2: usize,
-        pub weight1: f32,
-        pub weight2: f32,
-    }
-    impl TwoParticleAverageSite {
-        pub fn new(_p1: usize, _p2: usize, _weight1: f32, _weight2: f32) -> Self {
-            unimplemented!()
-        }
-    }
-    impl VirtualSite for TwoParticleAverageSite {}
-}
-
-#[cfg(test)]
-mod test {
-    use crate::force::{HarmonicAngleForce, HarmonicBondForce};
-    use crate::system::System;
-    use crate::virtual_site::TwoParticleAverageSite;
-
-    // Todo: Rust-OpenMM
-    #[ignore]
-    #[test]
-    fn test_system_particles() {
-        let num_particles = 10;
-        let mut system = System {};
-
-        for i in 0..num_particles {
-            system.add_particle(1.0 + 0.1 * i as f32);
-        }
-        system.set_particle_mass(5, 100.0);
-        assert_eq!(system.get_num_particles(), num_particles);
-
-        for i in 0..num_particles {
-            let expected_mass = if i == 5 { 100.0 } else { 1.0 + 0.1 * i as f32 };
-            assert_eq!(system.get_particle_mass(i), expected_mass);
-        }
-    }
-
-    // Todo: Rust-OpenMM
-    #[ignore]
-    #[test]
-    fn test_system_constraints() {
-        let num_particles = 10;
-        let mut system = System {};
-
-        for i in 0..num_particles - 1 {
-            system.add_constraint(i, i + 1, 0.2 * i as f32);
-        }
-        system.remove_constraint(5);
-        system.set_constraint_parameters(3, 0, 5, 99.0);
-        assert_eq!(system.get_num_constraints(), num_particles - 2);
-
-        for i in 0..num_particles - 2 {
-            let (p1, p2, dist) = system.get_constraint_parameters(i);
-            if i == 3 {
-                assert_eq!(p1, 0);
-                assert_eq!(p2, 5);
-                assert_eq!(dist, 99.0);
-            } else {
-                let j = if i < 5 { i } else { i + 1 };
-                assert_eq!(p1, j);
-                assert_eq!(p2, j + 1);
-                assert_eq!(dist, 0.2 * j as f32);
-            }
-        }
-    }
-
-    // Todo: Rust-OpenMM
-    #[ignore]
-    #[test]
-    fn test_system_force() {
-        let mut system = System {};
-        let bonds = HarmonicBondForce {};
-        system.add_force(Box::new(bonds));
-
-        let angles = HarmonicAngleForce {};
-        system.add_force(Box::new(angles));
-
-        assert_eq!(system.get_num_forces(), 2);
-        // Todo: Rust-OpenMM
-        // assert_eq!(system.get_force(0), Box::new(bonds));
-        // assert_eq!(system.get_force(1), Box::new(angles));
-
-        system.remove_force(0);
-        assert_eq!(system.get_num_forces(), 1);
-        // Todo: Rust-OpenMM
-        // assert_eq!(system.get_force(0), Box::new(angles));
-    }
-
-    // Todo: Rust-OpenMM
-    #[ignore]
-    #[test]
-    fn test_system_virtual_site() {
-        let mut system = System {};
-        let num_particles = 10;
-        for i in 0..num_particles {
-            system.add_particle(1.0);
-            assert!(!system.is_virtual_site(i));
-        }
-
-        let vsite = TwoParticleAverageSite::new(2, 3, 0.4, 0.6);
-        system.set_virtual_site(4, Some(Box::new(vsite)));
-        for i in 0..num_particles {
-            assert_eq!(system.is_virtual_site(i), i == 4);
-        }
-        // Todo: Rust-OpenMM
-        // assert_eq!(system.get_virtual_site(4), Some(Box::new(vsite)));
-        system.set_virtual_site(4, None);
-        for i in 0..num_particles {
-            assert!(!system.is_virtual_site(i));
-        }
-    }
-}
+pub mod integrator;
+pub mod virtual_site;

src/platform/mod.rs

@@ -0,0 +1,7 @@
+pub trait Platform {}
+
+pub struct CPU;
+pub struct OpenCL;
+
+impl Platform for CPU {}
+impl Platform for OpenCL {}

src/state.rs

@@ -0,0 +1,23 @@
+use nalgebra::Vector3;
+
+pub enum StateType {
+    Both,
+    Forces,
+    Energy,
+}
+
+pub struct State {
+    pub state_type: StateType,
+    pub forces: Vec<Vector3<f32>>,
+    pub energy: f32,
+}
+
+impl State {
+    pub fn get_forces(&self) -> &Vec<Vector3<f32>> {
+        &self.forces
+    }
+
+    pub fn get_potential_energy(&self) -> f32 {
+        self.energy
+    }
+}

src/system.rs

@@ -0,0 +1,195 @@
+use crate::force::Force;
+use crate::virtual_site::VirtualSite;
+use nalgebra::Vector3;
+pub struct System {}
+impl System {
+    pub fn add_particle(&mut self, _mass: f32) {
+        // Todo: Rust-OpenMM
+        unimplemented!()
+    }
+
+    pub fn get_num_particles(&self) -> usize {
+        // Todo: Rust-OpenMM
+        unimplemented!()
+    }
+
+    pub fn get_particle_mass(&self, _index: usize) -> f32 {
+        // Todo: Rust-OpenMM
+        unimplemented!()
+    }
+
+    pub fn set_particle_mass(&mut self, _index: usize, _mass: f32) {
+        // Todo: Rust-OpenMM
+        unimplemented!()
+    }
+
+    pub fn is_virtual_site(&self, _index: usize) -> bool {
+        // Todo: Rust-OpenMM
+        unimplemented!()
+    }
+
+    pub fn get_virtual_site(&self, _index: usize) -> Option<Box<dyn VirtualSite>> {
+        // Todo: Rust-OpenMM
+        unimplemented!()
+    }
+
+    pub fn set_virtual_site(&mut self, _index: usize, _vsite: Option<Box<dyn VirtualSite>>) {
+        // Todo: Rust-OpenMM
+        unimplemented!()
+    }
+
+    pub fn add_constraint(&mut self, _p1: usize, _p2: usize, _dist: f32) {
+        // Todo: Rust-OpenMM
+        unimplemented!()
+    }
+
+    pub fn remove_constraint(&mut self, _index: usize) {
+        // Todo: Rust-OpenMM
+        unimplemented!()
+    }
+
+    pub fn get_num_constraints(&self) -> usize {
+        // Todo: Rust-OpenMM
+        unimplemented!()
+    }
+
+    pub fn get_constraint_parameters(&self, _index: usize) -> (usize, usize, f32) {
+        // Todo: Rust-OpenMM
+        unimplemented!()
+    }
+
+    pub fn set_constraint_parameters(&mut self, _index: usize, _p1: usize, _p2: usize, _dist: f32) {
+        // Todo: Rust-OpenMM
+        unimplemented!()
+    }
+
+    pub fn add_force(&mut self, _force: Box<dyn Force>) {
+        // Todo: Rust-OpenMM
+        unimplemented!()
+    }
+
+    pub fn remove_force(&mut self, _index: usize) {
+        // Todo: Rust-OpenMM
+        unimplemented!()
+    }
+
+    pub fn get_num_forces(&self) -> usize {
+        // Todo: Rust-OpenMM
+        unimplemented!()
+    }
+
+    pub fn get_force(&self, _index: usize) -> Box<dyn Force> {
+        // Todo: Rust-OpenMM
+        unimplemented!()
+    }
+
+    pub fn set_default_periodic_box_vectors(
+        &mut self,
+        _a: Vector3<f32>,
+        _b: Vector3<f32>,
+        _c: Vector3<f32>,
+    ) {
+        // Todo: Rust-OpenMM
+        unimplemented!()
+    }
+}
+#[cfg(test)]
+mod test {
+    use crate::force::{HarmonicAngleForce, HarmonicBondForce};
+    use crate::system::System;
+    use crate::virtual_site::TwoParticleAverageSite;
+
+    // Todo: Rust-OpenMM
+    #[ignore]
+    #[test]
+    fn test_system_particles() {
+        let num_particles = 10;
+        let mut system = System {};
+
+        for i in 0..num_particles {
+            system.add_particle(1.0 + 0.1 * i as f32);
+        }
+        system.set_particle_mass(5, 100.0);
+        assert_eq!(system.get_num_particles(), num_particles);
+
+        for i in 0..num_particles {
+            let expected_mass = if i == 5 { 100.0 } else { 1.0 + 0.1 * i as f32 };
+            assert_eq!(system.get_particle_mass(i), expected_mass);
+        }
+    }
+
+    // Todo: Rust-OpenMM
+    #[ignore]
+    #[test]
+    fn test_system_constraints() {
+        let num_particles = 10;
+        let mut system = System {};
+
+        for i in 0..num_particles - 1 {
+            system.add_constraint(i, i + 1, 0.2 * i as f32);
+        }
+        system.remove_constraint(5);
+        system.set_constraint_parameters(3, 0, 5, 99.0);
+        assert_eq!(system.get_num_constraints(), num_particles - 2);
+
+        for i in 0..num_particles - 2 {
+            let (p1, p2, dist) = system.get_constraint_parameters(i);
+            if i == 3 {
+                assert_eq!(p1, 0);
+                assert_eq!(p2, 5);
+                assert_eq!(dist, 99.0);
+            } else {
+                let j = if i < 5 { i } else { i + 1 };
+                assert_eq!(p1, j);
+                assert_eq!(p2, j + 1);
+                assert_eq!(dist, 0.2 * j as f32);
+            }
+        }
+    }
+
+    // Todo: Rust-OpenMM
+    #[ignore]
+    #[test]
+    fn test_system_force() {
+        let mut system = System {};
+        let bonds = HarmonicBondForce {};
+        system.add_force(Box::new(bonds));
+
+        let angles = HarmonicAngleForce {};
+        system.add_force(Box::new(angles));
+
+        assert_eq!(system.get_num_forces(), 2);
+        // Todo: Rust-OpenMM
+        // assert_eq!(system.get_force(0), Box::new(bonds));
+        // assert_eq!(system.get_force(1), Box::new(angles));
+
+        system.remove_force(0);
+        assert_eq!(system.get_num_forces(), 1);
+        // Todo: Rust-OpenMM
+        // assert_eq!(system.get_force(0), Box::new(angles));
+    }
+
+    // Todo: Rust-OpenMM
+    #[ignore]
+    #[test]
+    fn test_system_virtual_site() {
+        let mut system = System {};
+        let num_particles = 10;
+        for i in 0..num_particles {
+            system.add_particle(1.0);
+            assert!(!system.is_virtual_site(i));
+        }
+
+        let vsite = TwoParticleAverageSite::new(2, 3, 0.4, 0.6);
+        system.set_virtual_site(4, Some(Box::new(vsite)));
+        for i in 0..num_particles {
+            assert_eq!(system.is_virtual_site(i), i == 4);
+        }
+        // Todo: Rust-OpenMM
+        // assert_eq!(system.get_virtual_site(4), Some(Box::new(vsite)));
+        system.set_virtual_site(4, None);
+        for i in 0..num_particles {
+            assert!(!system.is_virtual_site(i));
+        }
+    }
+}

src/virtual_site.rs

@@ -0,0 +1,17 @@
+use std::fmt::Debug;
+
+pub trait VirtualSite: Debug {}
+
+#[derive(Debug, Clone, PartialEq)]
+pub struct TwoParticleAverageSite {
+    pub p1: usize,
+    pub p2: usize,
+    pub weight1: f32,
+    pub weight2: f32,
+}
+impl TwoParticleAverageSite {
+    pub fn new(_p1: usize, _p2: usize, _weight1: f32, _weight2: f32) -> Self {
+        unimplemented!()
+    }
+}
+impl VirtualSite for TwoParticleAverageSite {}