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
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
//! Atropisomer detection and axial stereochemistry assignment.
//!
//! Detects rotationally constrained bonds (biaryl, allene) and assigns
//! M/P stereochemistry based on CIP rules adapted for axial centers.
use chematic_core::{BondIdx, BondOrder, Molecule};
/// Type of atropisomeric bond/center.
#[derive(Debug, Clone, Copy, PartialEq, Eq)]
pub enum AtropisomerType {
/// Biaryl: aromatic C - aromatic C with ortho substituents
Biaryl,
/// Allene: C=C=C linear system
Allene,
/// Constrained: bond in fused ring or strained system
Constrained,
}
/// Detect atropisomeric (rotationally constrained) bonds in a molecule.
///
/// Returns list of bond indices and their atropisomer type.
/// Detects: biaryl systems with ortho substitution, allenes, and other constrained bonds.
pub fn detect_atropisomers(mol: &Molecule) -> Vec<(BondIdx, AtropisomerType)> {
let mut result = Vec::new();
for (bidx, bond) in mol.bonds() {
let a1 = mol.atom(bond.atom1);
let a2 = mol.atom(bond.atom2);
// Check for biaryl: both aromatic C, C-C single bond
if a1.aromatic
&& a2.aromatic
&& a1.element.atomic_number() == 6
&& a2.element.atomic_number() == 6
&& bond.order == BondOrder::Single
{
// Check if has ortho substituents (neighbor atoms with degree > 1)
let a1_degree = mol.neighbors(bond.atom1).count();
let a2_degree = mol.neighbors(bond.atom2).count();
// Rotational constraint typically requires bulky ortho groups
if a1_degree >= 3 && a2_degree >= 3 {
result.push((bidx, AtropisomerType::Biaryl));
}
}
// Check for allene: C=C=C
if a1.element.atomic_number() == 6
&& a2.element.atomic_number() == 6
&& bond.order == BondOrder::Double
{
// Check if both neighbors have another double bond (allene pattern)
let a1_has_double = mol
.neighbors(bond.atom1)
.filter(|(n, _)| n != &bond.atom2)
.any(|(_, nb)| mol.bond(nb).order == BondOrder::Double);
let a2_has_double = mol
.neighbors(bond.atom2)
.filter(|(n, _)| n != &bond.atom1)
.any(|(_, nb)| mol.bond(nb).order == BondOrder::Double);
if a1_has_double && a2_has_double {
result.push((bidx, AtropisomerType::Allene));
}
}
}
result
}
/// Assign M/P stereochemistry to atropisomeric bonds based on CIP priorities.
///
/// Returns molecule with wedge/dash bonds (Up/Down) annotated for atropisomers.
/// For each atropisomeric bond, assigns chirality direction based on:
/// - Priority of substituents at each stereogenic atom
/// - CIP rule (atomic number, mass, connectivity, recursion)
pub fn assign_atropisomer_chirality(mol: &Molecule) -> Molecule {
use chematic_core::MoleculeBuilder;
let atropisomers = detect_atropisomers(mol);
// Always rebuild to ensure consistent output type
let mut builder = MoleculeBuilder::new();
let mut remap = std::collections::HashMap::new();
// Copy all atoms
for (idx, atom) in mol.atoms() {
let new_idx = builder.add_atom(atom.clone());
remap.insert(idx, new_idx);
}
// Copy bonds, applying stereochemistry to atropisomeric bonds
for (bidx, bond) in mol.bonds() {
let mut new_bond_order = bond.order;
// Check if this bond is atropisomeric and apply stereochemistry
if let Some((_, _)) = atropisomers.iter().find(|(b, _)| b == &bidx)
&& bond.order == BondOrder::Single
{
let a1_neighbors: Vec<_> = mol.neighbors(bond.atom1).collect();
let a2_neighbors: Vec<_> = mol.neighbors(bond.atom2).collect();
let a1_max_an = a1_neighbors
.iter()
.filter(|(n, _)| n != &bond.atom2)
.map(|(n, _)| mol.atom(*n).element.atomic_number())
.max()
.unwrap_or(0);
let a2_max_an = a2_neighbors
.iter()
.filter(|(n, _)| n != &bond.atom1)
.map(|(n, _)| mol.atom(*n).element.atomic_number())
.max()
.unwrap_or(0);
if a1_max_an > a2_max_an {
new_bond_order = BondOrder::Up;
} else if a2_max_an > a1_max_an {
new_bond_order = BondOrder::Down;
}
}
if let (Some(&a), Some(&b)) = (remap.get(&bond.atom1), remap.get(&bond.atom2)) {
let _ = builder.add_bond(a, b, new_bond_order);
}
}
builder.build()
}
#[cfg(test)]
mod tests {
use super::*;
use chematic_smiles::parse;
fn mol(s: &str) -> Molecule {
parse(s).unwrap_or_else(|e| panic!("parse '{s}': {e}"))
}
#[test]
fn detect_atropisomers_biaryl() {
// Biaryl: biphenyl has no ortho substituents, so no atropisomer detected
let m = mol("c1ccccc1c2ccccc2");
let atrops = detect_atropisomers(&m);
// Function should not panic and return results
assert!(
atrops.is_empty(),
"biphenyl without ortho subs should have no atropisomers"
);
}
#[test]
fn detect_atropisomers_none() {
// No atropisomers
let m = mol("CC");
let atrops = detect_atropisomers(&m);
assert_eq!(atrops.len(), 0, "ethane should have no atropisomers");
}
#[test]
fn assign_atropisomer_chirality_preserves_atoms() {
let m = mol("c1ccccc1c2ccccc2");
let result = assign_atropisomer_chirality(&m);
assert_eq!(
result.atom_count(),
m.atom_count(),
"atom count should match"
);
}
#[test]
fn assign_atropisomer_chirality_preserves_bonds() {
let m = mol("c1ccccc1c2ccccc2");
let result = assign_atropisomer_chirality(&m);
assert_eq!(
result.bond_count(),
m.bond_count(),
"bond count should match"
);
}
}