sci_form/forcefield/etkdg_3d/mod.rs
1//! ETKDG 3D force field matching RDKit's construct3DForceField.
2//! Uses flat-bottom distance constraints (zero gradient within tolerance)
3//! instead of harmonic bond stretch + angle bend.
4
5pub mod builder;
6pub mod energy;
7pub mod gradient;
8pub mod optimizer;
9
10pub use builder::*;
11pub use energy::*;
12pub use gradient::*;
13pub use optimizer::*;
14
15const KNOWN_DIST_TOL: f64 = 0.01;
16const KNOWN_DIST_K: f64 = 100.0;
17
18/// A pre-computed flat-bottom distance constraint
19#[derive(Clone)]
20pub struct DistConstraint {
21 pub i: usize,
22 pub j: usize,
23 pub min_len: f64,
24 pub max_len: f64,
25 pub k: f64,
26}
27
28/// Pre-computed M6 torsion contribution (Fourier series)
29#[derive(Clone)]
30pub struct M6TorsionContrib {
31 pub i: usize,
32 pub j: usize,
33 pub k: usize,
34 pub l: usize,
35 pub signs: [f64; 6],
36 pub v: [f64; 6],
37}
38
39/// Pre-computed UFF inversion contribution (out-of-plane)
40/// Matches RDKit's InversionContrib exactly:
41/// E = K * (C0 + C1*sinY + C2*cos2W)
42/// where Y is the Wilson out-of-plane angle.
43/// For SP2 C/N/O: C0=1, C1=-1, C2=0 → E = K*(1-sinY)
44#[derive(Clone)]
45pub struct UFFInversionContrib {
46 /// Neighbor atom 1 (idx1 in RDKit)
47 pub at1: usize,
48 /// Central atom (idx2 in RDKit)
49 pub at2: usize,
50 /// Neighbor atom 2 (idx3 in RDKit)
51 pub at3: usize,
52 /// Neighbor atom 3 (idx4 in RDKit)
53 pub at4: usize,
54 pub force_constant: f64,
55 pub c0: f64,
56 pub c1: f64,
57 pub c2: f64,
58}
59
60/// Flat-bottom angle constraint (degrees) for SP (linear) centers.
61/// E = k * angleTerm^2, where angleTerm = clamp(angle - [min,max], 0)
62#[derive(Clone)]
63pub struct AngleConstraint {
64 pub i: usize, // outer atom 1
65 pub j: usize, // central atom
66 pub k: usize, // outer atom 2
67 pub min_deg: f64,
68 pub max_deg: f64,
69 pub force_k: f64,
70}
71
72/// Pre-computed 3D ETKDG force field terms.
73/// Distance constraints are split into three groups matching RDKit's construct3DForceField
74/// ordering: 1-2 bonds, 1-3 angles, and long-range. This ensures identical
75/// floating-point accumulation order in energy/gradient, which is critical
76/// for BFGS trajectory reproducibility.
77pub struct Etkdg3DFF {
78 pub dist_12: Vec<DistConstraint>,
79 pub dist_13: Vec<DistConstraint>,
80 pub dist_long: Vec<DistConstraint>,
81 pub angle_constraints: Vec<AngleConstraint>,
82 pub torsion_contribs: Vec<M6TorsionContrib>,
83 pub inversion_contribs: Vec<UFFInversionContrib>,
84 pub oop_k: f64,
85 pub torsion_k_omega: f64,
86 pub bounds_force_scaling: f64,
87 pub use_m6_torsions: bool,
88}