Skip to main content

conspire/physics/molecular/single_chain/frc/
mod.rs

1#[cfg(test)]
2mod test;
3
4use crate::math::Current;
5use crate::{
6    math::{CrossProduct, Quantity, Scalar, random::random_uniform},
7    mechanics::Vector,
8    physics::molecular::single_chain::{
9        Configuration, Ensemble, Inextensible, Isometric, Isotensional, Legendre, MonteCarlo,
10        SingleChain, SingleChainError, Thermodynamics,
11    },
12    units::Length,
13};
14use std::f64::consts::TAU;
15
16/// The freely-rotating chain model.
17#[derive(Clone, Debug)]
18pub struct FreelyRotatingChain {
19    /// The link angle $`\theta_b`$.
20    pub link_angle: Scalar,
21    /// The link length $`\ell_b`$.
22    pub link_length: Scalar,
23    /// The number of links $`N_b`$.
24    pub number_of_links: u8,
25    /// The thermodynamic ensemble.
26    pub ensemble: Ensemble,
27}
28
29impl SingleChain for FreelyRotatingChain {
30    fn link_length(&self) -> Quantity<Length> {
31        Quantity::new(self.link_length)
32    }
33    fn number_of_links(&self) -> u8 {
34        self.number_of_links
35    }
36}
37
38impl Inextensible for FreelyRotatingChain {
39    fn maximum_nondimensional_extension(&self) -> Scalar {
40        1.0
41    }
42}
43
44impl Thermodynamics for FreelyRotatingChain {
45    fn ensemble(&self) -> Ensemble {
46        self.ensemble
47    }
48}
49
50impl Isometric for FreelyRotatingChain {
51    fn nondimensional_helmholtz_free_energy(
52        &self,
53        _nondimensional_force: Scalar,
54    ) -> Result<Scalar, SingleChainError> {
55        unimplemented!()
56    }
57    fn nondimensional_force(
58        &self,
59        _nondimensional_force: Scalar,
60    ) -> Result<Scalar, SingleChainError> {
61        unimplemented!()
62    }
63    fn nondimensional_stiffness(
64        &self,
65        _nondimensional_force: Scalar,
66    ) -> Result<Scalar, SingleChainError> {
67        unimplemented!()
68    }
69    fn nondimensional_spherical_distribution(
70        &self,
71        _nondimensional_force: Scalar,
72    ) -> Result<Scalar, SingleChainError> {
73        unimplemented!()
74    }
75}
76
77impl Isotensional for FreelyRotatingChain {
78    fn nondimensional_gibbs_free_energy_per_link(
79        &self,
80        _nondimensional_force: Scalar,
81    ) -> Result<Scalar, SingleChainError> {
82        unimplemented!()
83    }
84    fn nondimensional_extension(
85        &self,
86        _nondimensional_force: Scalar,
87    ) -> Result<Scalar, SingleChainError> {
88        unimplemented!()
89    }
90    fn nondimensional_compliance(
91        &self,
92        _nondimensional_force: Scalar,
93    ) -> Result<Scalar, SingleChainError> {
94        unimplemented!()
95    }
96}
97
98impl Legendre for FreelyRotatingChain {}
99
100impl MonteCarlo for FreelyRotatingChain {
101    fn random_nondimensional_link_vectors(&self, nondimensional_force: Scalar) -> Configuration {
102        if nondimensional_force != 0.0 {
103            unimplemented!()
104        }
105        let cos_theta = 2.0 * random_uniform() - 1.0;
106        let sin_theta = (1.0 - cos_theta * cos_theta).sqrt();
107        let phi = TAU * random_uniform();
108        let (sin_phi, cos_phi) = phi.sin_cos();
109        const AY: Vector<Current> = Vector::<Current>::const_from([0.0, 1.0, 0.0]);
110        const AZ: Vector<Current> = Vector::<Current>::const_from([0.0, 0.0, 1.0]);
111        let mut a = AY;
112        let mut b =
113            Vector::<Current>::const_from([sin_theta * cos_phi, sin_theta * sin_phi, cos_theta]);
114        let (sin_theta, cos_theta) = self.link_angle.sin_cos();
115        (0..self.number_of_links())
116            .map(|link| {
117                if link > 0 {
118                    a = if b[1].abs() < 0.9 { AY } else { AZ };
119                    let u = a.cross(&b).normalized();
120                    let v = b.cross(&u);
121                    let phi = TAU * random_uniform();
122                    let (sin_phi, cos_phi) = phi.sin_cos();
123                    b = &b * cos_theta + (&u * cos_phi + &v * sin_phi) * sin_theta;
124                }
125                b.clone()
126            })
127            .collect()
128    }
129}