Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
124 changes: 61 additions & 63 deletions gufe/components/explicitmoleculecomponent.py
Original file line number Diff line number Diff line change
Expand Up @@ -35,68 +35,6 @@ def _ensure_ofe_name(mol: RDKitMol, name: str) -> str:
return name


def _check_partial_charges(mol: RDKitMol, logger=None) -> None:
"""
Checks for the presence of partial charges.

Note
----
We ensure the charges are set as atom properties
to ensure they are detected by OpenFF

Raises
------
ValueError
* If the partial charges are not of length atoms.
* If the sum of partial charges is not equal to the
formal charge.
UserWarning
* If partial charges are found.
* If the partial charges are near 0 for all atoms.
"""
if "atom.dprop.PartialCharge" not in mol.GetPropNames():
return

p_chgs = np.array(mol.GetProp("atom.dprop.PartialCharge").split(), dtype=float)

if len(p_chgs) != mol.GetNumAtoms():
errmsg = f"Incorrect number of partial charges: {len(p_chgs)} " f" were provided for {mol.GetNumAtoms()} atoms"
raise ValueError(errmsg)

if abs(sum(p_chgs) - Chem.GetFormalCharge(mol)) > 0.01:
errmsg = (
f"Sum of partial charges {sum(p_chgs)} differs from " f"RDKit formal charge {Chem.GetFormalCharge(mol)}"
)
raise ValueError(errmsg)

# set the charges on the atoms if not already set
for i, charge in enumerate(p_chgs):
atom = mol.GetAtomWithIdx(i)
if not atom.HasProp("PartialCharge"):
atom.SetDoubleProp("PartialCharge", charge)
else:
atom_charge = atom.GetDoubleProp("PartialCharge")
if not np.isclose(atom_charge, charge):
errmsg = (
f"non-equivalent partial charges between atom and " f"molecule properties: {atom_charge} {charge}"
)
raise ValueError(errmsg)

if np.all(np.isclose(p_chgs, 0.0)):
wmsg = "Partial charges provided all equal to zero. These may be ignored by some Protocols."
warnings.warn(wmsg)
else:
msg = "Partial charges have been provided"
if name := mol.GetProp("ofe-name"):
msg += f" for {name}"
msg += ", these will preferentially be used instead of generating new partial charges."

if logger is None:
logger = logging.getLogger(__name__)

logger.info(msg)


class ExplicitMoleculeComponent(Component):
"""Base class for explicit molecules.

Expand All @@ -110,7 +48,6 @@ class ExplicitMoleculeComponent(Component):

def __init__(self, rdkit: RDKitMol, name: str = ""):
name = _ensure_ofe_name(rdkit, name)
_check_partial_charges(rdkit, logger=self.logger)
conformers = list(rdkit.GetConformers())
if not conformers:
raise ValueError("Molecule was provided with no conformers.")
Expand All @@ -129,6 +66,8 @@ def __init__(self, rdkit: RDKitMol, name: str = ""):
self._smiles: Optional[str] = None
self._name = name

self._check_partial_charges()

def __getstate__(self):
# TODO: check that RDKit setting is set before issuing warning
if Chem.GetDefaultPickleProperties() != Chem.PropertyPickleOptions.AllProps:
Expand Down Expand Up @@ -186,3 +125,62 @@ def _from_dict(cls, d: dict):

def _to_dict(self) -> dict:
raise NotImplementedError()

def _check_partial_charges(self) -> None:
"""
Checks for the presence of partial charges.

Note
----
We ensure the charges are set as atom properties
to ensure they are detected by OpenFF

Raises
------
ValueError
* If the partial charges are not of length atoms.
* If the sum of partial charges is not equal to the
formal charge.
UserWarning
* If partial charges are found.
* If the partial charges are near 0 for all atoms.
"""
mol = self._rdkit

if "atom.dprop.PartialCharge" not in mol.GetPropNames():
return

p_chgs = np.array(mol.GetProp("atom.dprop.PartialCharge").split(), dtype=float)

if len(p_chgs) != mol.GetNumAtoms():
errmsg = (
f"Incorrect number of partial charges: {len(p_chgs)} " f" were provided for {mol.GetNumAtoms()} atoms"
)
raise ValueError(errmsg)

if abs(sum(p_chgs) - Chem.GetFormalCharge(mol)) > 0.01:
errmsg = (
f"Sum of partial charges {sum(p_chgs)} differs from " f"RDKit formal charge {Chem.GetFormalCharge(mol)}"
)
raise ValueError(errmsg)

# set the charges on the atoms if not already set
for i, charge in enumerate(p_chgs):
atom = mol.GetAtomWithIdx(i)
if not atom.HasProp("PartialCharge"):
atom.SetDoubleProp("PartialCharge", charge)
else:
atom_charge = atom.GetDoubleProp("PartialCharge")
if not np.isclose(atom_charge, charge):
errmsg = (
f"non-equivalent partial charges between atom and "
f"molecule properties: {atom_charge} {charge}"
)
raise ValueError(errmsg)

if np.all(np.isclose(p_chgs, 0.0)):
wmsg = "Partial charges provided all equal to zero. These may be ignored by some Protocols."
warnings.warn(wmsg)
else:
msg = f"Partial charges are present for {self.key} (name: '{self.name}')"
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

A note for others, we went with gufekey and the name of the component as writing out the smiles to the log might not be a good idea as this can be sensitive information that should not be shared.

self.logger.info(msg)
12 changes: 3 additions & 9 deletions gufe/tests/test_smallmoleculecomponent.py
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,7 @@

import gufe
from gufe import SmallMoleculeComponent
from gufe.components.explicitmoleculecomponent import _check_partial_charges, _ensure_ofe_name
from gufe.components.explicitmoleculecomponent import _ensure_ofe_name
from gufe.tokenization import TOKENIZABLE_REGISTRY

from .test_explicitmoleculecomponent import ExplicitMoleculeComponentMixin
Expand Down Expand Up @@ -248,17 +248,11 @@ def charged_off_ethane(self, named_ethane):
off_ethane.assign_partial_charges(partial_charge_method="am1bcc")
return off_ethane

def test_check_partial_charges_without_gufe_logger(self, charged_off_ethane, caplog):
rd_mol = charged_off_ethane.to_rdkit()
caplog.set_level(logging.INFO)
_check_partial_charges(rd_mol, logger=None)
assert "Partial charges have been provided for ethane, these" in caplog.text

def test_partial_charges_logging(self, charged_off_ethane, caplog):
caplog.set_level(logging.INFO)
SmallMoleculeComponent.from_openff(charged_off_ethane)

assert "Partial charges have been provided" in caplog.text
assert "Partial charges are present for SmallMoleculeComponent-" in caplog.text

def test_partial_charges_zero_warning(self, charged_off_ethane):
charged_off_ethane.partial_charges[:] = 0 * unit.elementary_charge
Expand Down Expand Up @@ -292,7 +286,7 @@ def test_partial_charges_applied_to_atoms(self, caplog):
caplog.set_level(logging.INFO)

ofe = SmallMoleculeComponent.from_rdkit(mol)
assert "Partial charges have been provided" in caplog.text
assert "Partial charges are present for SmallMoleculeComponent-" in caplog.text

# convert to openff and make sure the charges are set
off_mol = ofe.to_openff()
Expand Down
23 changes: 23 additions & 0 deletions news/partial_charge_warning_to_logger.rst
Original file line number Diff line number Diff line change
@@ -0,0 +1,23 @@
**Added:**

* <news item>

**Changed:**

* The message stating that partial charges are already present in an ``ExplicitMoleculeComponent`` is now included in ``logger.info``, rather than as a warning message. This should make output significantly less noisy for some users.

**Deprecated:**

* <news item>

**Removed:**

* <news item>

**Fixed:**

* <news item>

**Security:**

* <news item>