Skip to main content

sci_form/forcefield/
mmff94.rs

1use super::traits::ForceFieldContribution;
2
3/// Dispersión estérica repulsiva/atractiva regida por el Potencial Amortiguado 14-7 (Buffered 14-7) de Halgren.
4pub struct Mmff94BufferedVanDerWaals {
5    pub atom_i_idx: usize,
6    pub atom_j_idx: usize,
7    pub radius_star: f64,   // Parámetro dimensional empírico cruzado R*ij
8    pub epsilon_depth: f64, // Factor de profundidad termodinámica eps_ij
9}
10
11impl ForceFieldContribution for Mmff94BufferedVanDerWaals {
12    fn evaluate_energy_and_inject_gradient(&self, coords: &[f64], grad: &mut [f64]) -> f64 {
13        let root_i = self.atom_i_idx * 3;
14        let root_j = self.atom_j_idx * 3;
15
16        let delta_x = coords[root_i] - coords[root_j];
17        let delta_y = coords[root_i + 1] - coords[root_j + 1];
18        let delta_z = coords[root_i + 2] - coords[root_j + 2];
19
20        let dist_squared = delta_x * delta_x + delta_y * delta_y + delta_z * delta_z;
21        let mut dist_r = dist_squared.sqrt();
22
23        // Tope asintótico absoluto inferior para colisiones
24        if dist_r < 1e-8 {
25            dist_r = 1e-8;
26        }
27
28        // Algebra fraccionaria amortiguadora: E_vdW = eps * (1.07 R* / (R + 0.07 R*))^7 * ((1.12 R*^7 / (R^7 + 0.12 R*^7)) - 2)
29        let r_star_powered_7 = self.radius_star.powi(7);
30        let dist_r_powered_7 = dist_r.powi(7);
31
32        let repulsive_denominator = dist_r + 0.07 * self.radius_star;
33        let repulsive_term = (1.07 * self.radius_star / repulsive_denominator).powi(7);
34
35        let attractive_denominator = dist_r_powered_7 + 0.12 * r_star_powered_7;
36        let attractive_term = (1.12 * r_star_powered_7 / attractive_denominator) - 2.0;
37
38        let vdw_total_energy = self.epsilon_depth * repulsive_term * attractive_term;
39
40        // Derivación espacial analítica
41        let gradient_rep_term = -7.0 * repulsive_term / repulsive_denominator;
42        let gradient_attr_term = -7.0 * dist_r.powi(6) * (1.12 * r_star_powered_7)
43            / (attractive_denominator * attractive_denominator);
44
45        let force_scalar_magnitude = self.epsilon_depth
46            * (gradient_rep_term * attractive_term + repulsive_term * gradient_attr_term);
47
48        // Factorización cartesiana
49        let vector_prefactor = force_scalar_magnitude / dist_r;
50        let grad_x = vector_prefactor * delta_x;
51        let grad_y = vector_prefactor * delta_y;
52        let grad_z = vector_prefactor * delta_z;
53
54        grad[root_i] += grad_x;
55        grad[root_i + 1] += grad_y;
56        grad[root_i + 2] += grad_z;
57
58        grad[root_j] -= grad_x;
59        grad[root_j + 1] -= grad_y;
60        grad[root_j + 2] -= grad_z;
61
62        vdw_total_energy
63    }
64}