Skip to main content

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

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