sci-form 0.14.2

High-performance 3D molecular conformer generation using ETKDG distance geometry
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61

//! Triangle inequality smoothing for distance bounds matrices.
//!
//! Enforces the triangle inequality on all (i, j, k) triples:
//!   upper(i,j) ≤ upper(i,k) + upper(k,j)
//!   lower(i,j) ≥ lower(i,k) − upper(k,j)
//!
//! Returns `false` if the bounds are infeasible (lower > upper after tightening).

use nalgebra::DMatrix;

pub fn triangle_smooth(b: &mut DMatrix<f64>) -> bool {
    triangle_smooth_tol(b, 0.0)
}

pub fn triangle_smooth_tol(b: &mut DMatrix<f64>, tol: f64) -> bool {
    let n = b.nrows();
    for k in 0..n {
        for i in 0..n {
            if i == k {
                continue;
            }
            let (uk, lk) = if i < k {
                (b[(i, k)], b[(k, i)])
            } else {
                (b[(k, i)], b[(i, k)])
            };
            for j in (i + 1)..n {
                if j == k {
                    continue;
                }
                let (ukj, lkj) = if j < k {
                    (b[(j, k)], b[(k, j)])
                } else {
                    (b[(k, j)], b[(j, k)])
                };
                let (su, d1, d2) = (uk + ukj, lk - ukj, lkj - uk);
                if b[(i, j)] > su {
                    b[(i, j)] = su;
                }
                // Lower bound must satisfy BOTH triangle constraints
                let li = b[(j, i)].max(d1).max(d2);
                b[(j, i)] = li;
                let lower = b[(j, i)];
                let upper = b[(i, j)];
                if lower > upper {
                    if tol > 0.0 && lower > 1e-6 && (lower - upper) / lower < tol {
                        b[(i, j)] = lower; // set upper = lower
                    } else {
                        return false;
                    }
                }
            }
        }
    }
    true
}

pub fn smooth_bounds_matrix(mut b: DMatrix<f64>) -> DMatrix<f64> {
    triangle_smooth(&mut b);
    b
}