chematic_ff/
params.rs

1use chematic_core::BondOrder;
2use std::f64::consts::PI;
3
4use crate::dreiding::DREIDINGType;
5
6pub fn dreiding_vdw(t: DREIDINGType) -> (f64, f64) {
7    let (r0, depth) = match t {
8        DREIDINGType::H_ => (2.646, 0.044),
9        DREIDINGType::C_3 => (3.473, 0.066),
10        DREIDINGType::C_2 => (3.473, 0.066),
11        DREIDINGType::C_1 => (3.473, 0.066),
12        DREIDINGType::C_R => (3.473, 0.066),
13        DREIDINGType::N_3 => (3.243, 0.068),
14        DREIDINGType::N_2 => (3.243, 0.068),
15        DREIDINGType::N_1 => (3.243, 0.068),
16        DREIDINGType::N_R => (3.243, 0.068),
17        DREIDINGType::O_3 => (3.144, 0.052),
18        DREIDINGType::O_2 => (3.144, 0.052),
19        DREIDINGType::O_R => (3.144, 0.052),
20        DREIDINGType::S_3 => (3.952, 0.274),
21        DREIDINGType::S_R => (3.952, 0.274),
22        DREIDINGType::P_3 => (3.741, 0.210),
23        DREIDINGType::F_ => (3.048, 0.050),
24        DREIDINGType::Cl => (3.516, 0.227),
25        DREIDINGType::Br => (3.641, 0.389),
26        DREIDINGType::I_ => (3.984, 0.680),
27        DREIDINGType::X_ => (3.473, 0.066),
28    };
29    (r0, depth)
30}
31
32pub fn dreiding_bond_len(t1: DREIDINGType, t2: DREIDINGType, order: BondOrder) -> f64 {
33    let key = (
34        dreiding_type_order(t1),
35        dreiding_type_order(t2),
36        bond_order_value(order),
37    );
38
39    match key {
40        (1, 1, 1) => 1.090,
41        (1, 6, 1) => 1.090,
42        (1, 7, 1) => 1.010,
43        (1, 8, 1) => 0.960,
44        (1, 16, 1) => 1.336,
45        (1, 17, 1) => 0.980,
46        (1, 35, 1) => 1.084,
47        (1, 53, 1) => 1.008,
48
49        (6, 6, 1) => 1.540,
50        (6, 6, 2) => 1.340,
51        (6, 6, 3) => 1.204,
52        (6, 6, 4) => 1.395,
53        (6, 7, 1) => 1.469,
54        (6, 7, 2) => 1.279,
55        (6, 7, 3) => 1.158,
56        (6, 7, 4) => 1.340,
57        (6, 8, 1) => 1.427,
58        (6, 8, 2) => 1.217,
59        (6, 8, 4) => 1.370,
60        (6, 16, 1) => 1.820,
61        (6, 17, 1) => 1.770,
62        (6, 35, 1) => 1.940,
63        (6, 53, 1) => 2.140,
64
65        (7, 7, 1) => 1.450,
66        (7, 7, 2) => 1.225,
67        (7, 7, 3) => 1.098,
68        (7, 8, 1) => 1.400,
69        (7, 8, 2) => 1.218,
70        (7, 16, 1) => 1.754,
71        (7, 17, 1) => 1.690,
72
73        (8, 8, 1) => 1.480,
74        (8, 8, 2) => 1.208,
75        (8, 16, 1) => 1.680,
76        (8, 17, 1) => 1.611,
77
78        (16, 16, 1) => 2.048,
79        (16, 16, 2) => 1.549,
80
81        _ => match order {
82            BondOrder::Triple => 1.20,
83            BondOrder::Double => 1.34,
84            BondOrder::Aromatic => 1.40,
85            _ => 1.54,
86        },
87    }
88}
89
90pub fn dreiding_angle(t: DREIDINGType) -> f64 {
91    let deg = match t {
92        DREIDINGType::C_3 => 109.47,
93        DREIDINGType::C_2 => 120.0,
94        DREIDINGType::C_1 => 180.0,
95        DREIDINGType::C_R => 120.0,
96        DREIDINGType::N_3 => 106.7,
97        DREIDINGType::N_2 => 120.0,
98        DREIDINGType::N_1 => 180.0,
99        DREIDINGType::N_R => 120.0,
100        DREIDINGType::O_3 => 104.51,
101        DREIDINGType::O_2 => 120.0,
102        DREIDINGType::O_R => 120.0,
103        DREIDINGType::S_3 => 99.0,
104        DREIDINGType::S_R => 120.0,
105        DREIDINGType::P_3 => 93.0,
106        DREIDINGType::H_ => 109.47,
107        DREIDINGType::F_ => 109.47,
108        DREIDINGType::Cl => 109.47,
109        DREIDINGType::Br => 109.47,
110        DREIDINGType::I_ => 109.47,
111        DREIDINGType::X_ => 109.47,
112    };
113    deg * PI / 180.0
114}
115
116pub fn dreiding_torsion_barrier(_t1: DREIDINGType, _t2: DREIDINGType) -> f64 {
117    2.0
118}
119
120fn dreiding_type_order(t: DREIDINGType) -> u8 {
121    match t {
122        DREIDINGType::H_ => 1,
123        DREIDINGType::C_3 | DREIDINGType::C_2 | DREIDINGType::C_1 | DREIDINGType::C_R => 6,
124        DREIDINGType::N_3 | DREIDINGType::N_2 | DREIDINGType::N_1 | DREIDINGType::N_R => 7,
125        DREIDINGType::O_3 | DREIDINGType::O_2 | DREIDINGType::O_R => 8,
126        DREIDINGType::S_3 | DREIDINGType::S_R => 16,
127        DREIDINGType::P_3 => 15,
128        DREIDINGType::F_ => 9,
129        DREIDINGType::Cl => 17,
130        DREIDINGType::Br => 35,
131        DREIDINGType::I_ => 53,
132        DREIDINGType::X_ => 0,
133    }
134}
135
136fn bond_order_value(order: BondOrder) -> u8 {
137    match order {
138        BondOrder::Single => 1,
139        BondOrder::Double => 2,
140        BondOrder::Triple => 3,
141        BondOrder::Quadruple => 4,
142        BondOrder::Aromatic => 4,
143        _ => 1,
144    }
145}
146
147#[cfg(test)]
148mod tests {
149    use super::*;
150
151    #[test]
152    fn test_vdw_carbon() {
153        let (r0, depth) = dreiding_vdw(DREIDINGType::C_3);
154        assert!((r0 - 3.473).abs() < 0.001);
155        assert!((depth - 0.066).abs() < 0.001);
156    }
157
158    #[test]
159    fn test_bond_len_cc_single() {
160        let r = dreiding_bond_len(DREIDINGType::C_3, DREIDINGType::C_3, BondOrder::Single);
161        assert!((r - 1.540).abs() < 0.001);
162    }
163
164    #[test]
165    fn test_bond_len_cc_double() {
166        let r = dreiding_bond_len(DREIDINGType::C_2, DREIDINGType::C_2, BondOrder::Double);
167        assert!((r - 1.340).abs() < 0.001);
168    }
169
170    #[test]
171    fn test_angle_sp3() {
172        let theta = dreiding_angle(DREIDINGType::C_3);
173        let expected = 109.47 * PI / 180.0;
174        assert!((theta - expected).abs() < 0.01);
175    }
176
177    #[test]
178    fn test_angle_sp2() {
179        let theta = dreiding_angle(DREIDINGType::C_2);
180        let expected = 120.0 * PI / 180.0;
181        assert!((theta - expected).abs() < 0.01);
182    }
183}
chematic_ff/params.rs

chematic_ff/
params.rs
