Skip to content

Commit

Permalink
Added fix_mol_electronic_configuration() to species
Browse files Browse the repository at this point in the history
Also added get_atom_theoretical_charge()
  • Loading branch information
alongd committed Oct 2, 2023
1 parent 438af0c commit 793fdec
Showing 1 changed file with 62 additions and 5 deletions.
67 changes: 62 additions & 5 deletions arc/species/species.py
Original file line number Diff line number Diff line change
Expand Up @@ -1935,11 +1935,7 @@ def _scissors(self,
added_radical = list()
for mol, label in zip([mol1, mol2], [label1, label2]):
for atom in mol.atoms:
theoretical_charge = elements.PeriodicSystem.valence_electrons[atom.symbol] \
- atom.get_total_bond_order() \
- atom.radical_electrons - \
2 * atom.lone_pairs
if theoretical_charge == atom.charge + 1:
if get_atom_theoretical_charge(atom) == atom.charge + 1:
# we're missing a radical electron on this atom
if label not in added_radical or label == 'H':
atom.radical_electrons += 1
Expand Down Expand Up @@ -2708,3 +2704,64 @@ def split_mol(mol: Molecule) -> Tuple[List[Molecule], List[List[int]]]:
molecules.append(Molecule(atoms=[mol.atoms[index] for index in frag_indices]))
fragments.append(frag_indices)
return molecules, fragments


def get_atom_theoretical_charge(atom: Atom) -> float:
"""
Get the theoretical charge of an atom based on its electronic properties.
Args:
atom (Atom): The atom to check.
Returns:
float: The theoretical charge of the atom.
"""
return elements.PeriodicSystem.valence_electrons[atom.symbol] \
- atom.get_total_bond_order() \
- atom.radical_electrons - \
2 * atom.lone_pairs


def fix_mol_electronic_configuration(mol: Optional[Molecule] = None) -> Optional[Molecule]:
"""
Fix the electronic configuration (number of radicals, lone pairs, and formal charges) of a Molecule object
that was perceived from xyz. Assumes the multiplicity attribute of the given molecule is correct.
Args:
mol (Molecule, optional): The molecule to fix.
Returns:
Molecule: The fixed molecule.
"""
if mol is None:
return None
multiplicity = mol.multiplicity
n_radicals = mol.get_radical_count()
if n_radicals + 1 > mol.multiplicity:
# Find adjacent atoms with opposite-spin radicals, and pair the radicals up into a bond.
mols = mol.generate_resonance_structures(keep_isomorphic=False, filter_structures=True, save_order=True)
for structure in mols + [mol]:
structure.multiplicity = multiplicity
visited = list()
updated = False
for mol_1 in mols:
mol_1.multiplicity = mol.multiplicity
for atom_1 in mol_1.atoms:
if atom_1 not in visited and atom_1.radical_electrons >= 1:
for atom_2, bond in atom_1.edges.items():
if atom_2.radical_electrons >= 1:
atom_1.radical_electrons -= 1
atom_2.radical_electrons -= 1
bond.order += 1
visited.append(atom_2)
updated = True
break
visited.append(atom_1)
# Update formal charges.
for atom in mol_1.atoms:
theoretical_charge = get_atom_theoretical_charge(atom)
if theoretical_charge != atom.charge:
atom.charge = theoretical_charge
if updated:
return mol_1
return mols[0]

0 comments on commit 793fdec

Please sign in to comment.