diff --git a/.github/workflows/beta.yaml b/.github/workflows/beta.yaml index 6b48d625a..46b87ef5e 100644 --- a/.github/workflows/beta.yaml +++ b/.github/workflows/beta.yaml @@ -16,7 +16,7 @@ jobs: fail-fast: false matrix: os: - - macos-12 + - macos-latest - ubuntu-latest python-version: - "3.10" @@ -33,6 +33,8 @@ jobs: steps: - uses: actions/checkout@v4 + with: + fetch-depth: 0 - name: Install conda environment uses: mamba-org/setup-micromamba@v1 diff --git a/.github/workflows/ci.yaml b/.github/workflows/ci.yaml index 1ae6a5579..97bb05500 100644 --- a/.github/workflows/ci.yaml +++ b/.github/workflows/ci.yaml @@ -24,7 +24,7 @@ jobs: strategy: matrix: os: - - macos-12 + - macos-latest - ubuntu-latest python-version: - "3.10" @@ -44,6 +44,8 @@ jobs: steps: - uses: actions/checkout@v4 + with: + fetch-depth: 0 - name: Install conda environment uses: mamba-org/setup-micromamba@v1 @@ -124,7 +126,7 @@ jobs: python devtools/scripts/molecule-regressions.py - name: Run mypy - if: ${{ matrix.python-version == '3.10' }} + if: ${{ matrix.python-version == '3.11' }} run: | # As of 01/23, JAX with mypy is too slow to use without a pre-built cache # https://github.com/openforcefield/openff-interchange/pull/578#issuecomment-1369979875 diff --git a/.github/workflows/examples.yaml b/.github/workflows/examples.yaml index 2f13fa39f..ac174f8d0 100644 --- a/.github/workflows/examples.yaml +++ b/.github/workflows/examples.yaml @@ -25,7 +25,7 @@ jobs: matrix: os: - ubuntu-latest - - macos-12 + - macos-latest python-version: - "3.10" @@ -35,6 +35,8 @@ jobs: steps: - uses: actions/checkout@v4 + with: + fetch-depth: 0 - name: Install conda environment uses: mamba-org/setup-micromamba@v1 @@ -46,9 +48,6 @@ jobs: - name: Install package run: | - # These packages are brought in by conda (via the toolkit) and must be removed manually - # since pip doesn't know about the -base split and does not uninstall the -base package - micromamba remove --force openff-interchange openff-interchange-base python -m pip install . - name: Environment Information @@ -65,7 +64,7 @@ jobs: - name: Run docexamples run: | - pytest --doctest-modules openff/interchange/ --ignore=openff/interchange/_tests + python -m pytest --doctest-modules openff/interchange/ --ignore=openff/interchange/_tests - name: Run example notebooks if: always() diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index 18b9feba1..80dd27779 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -1,6 +1,6 @@ repos: - repo: https://github.com/pre-commit/pre-commit-hooks - rev: v4.6.0 + rev: v5.0.0 hooks: - id: check-yaml - id: end-of-file-fixer @@ -12,7 +12,7 @@ repos: hooks: - id: add-trailing-comma - repo: https://github.com/astral-sh/ruff-pre-commit - rev: v0.6.2 + rev: v0.7.3 hooks: - id: ruff - id: ruff-format @@ -23,18 +23,18 @@ repos: args: ["openff/interchange/", "-e", "openff/interchange/_tests/"] pass_filenames: false - repo: https://github.com/adamchainz/blacken-docs - rev: 1.18.0 + rev: 1.19.1 hooks: - id: blacken-docs files: ^docs/ - repo: https://github.com/igorshubovych/markdownlint-cli - rev: v0.41.0 + rev: v0.42.0 hooks: - id: markdownlint exclude: .github args: ["--disable", "MD013", "MD033", "MD024", "MD046", "--ignore", "docs/using/experimental.md", "--ignore", "docs/using/status.md", "--"] - repo: https://github.com/nbQA-dev/nbQA - rev: 1.8.7 + rev: 1.9.1 hooks: - id: nbqa-pyupgrade files: ^examples @@ -45,12 +45,12 @@ repos: args: - --fix - repo: https://github.com/kynan/nbstripout - rev: 0.7.1 + rev: 0.8.0 hooks: - id: nbstripout files: ^examples - repo: https://github.com/asottile/pyupgrade - rev: v3.17.0 + rev: v3.19.0 hooks: - id: pyupgrade args: diff --git a/.readthedocs.yml b/.readthedocs.yml index 5757a3610..6f88bb5ed 100644 --- a/.readthedocs.yml +++ b/.readthedocs.yml @@ -1,16 +1,23 @@ -# .readthedocs.yml - version: 2 build: - os: ubuntu-22.04 + os: ubuntu-24.04 tools: - # https://docs.readthedocs.io/en/stable/config-file/v2.html#build-tools-python - python: "mambaforge-4.10" + python: "mambaforge-latest" + jobs: + pre_build: + - mamba list sphinx: configuration: docs/conf.py + # Interchange's docs should not produce warnings; if they start, + # something is broken (most likely automatic API doc generation) fail_on_warning: true conda: environment: devtools/conda-envs/docs_env.yaml + +python: + install: + - method: pip + path: . diff --git a/MANIFEST.in b/MANIFEST.in index 5e0270c12..e5dc6aa83 100644 --- a/MANIFEST.in +++ b/MANIFEST.in @@ -1,6 +1,5 @@ include LICENSE include MANIFEST.in -include versioneer.py graft openff global-exclude *.py[cod] __pycache__ *.so diff --git a/devtools/conda-envs/beta_env.yaml b/devtools/conda-envs/beta_env.yaml index f6c952343..a8d8a2bb4 100644 --- a/devtools/conda-envs/beta_env.yaml +++ b/devtools/conda-envs/beta_env.yaml @@ -8,7 +8,7 @@ dependencies: - python - numpy >=1.21 - pydantic >=1.10.17,<3 - - openmm >=7.6 + - openmm # OpenFF stack - openff-toolkit ~=0.16.4 - openff-nagl ~=0.3.7 diff --git a/devtools/conda-envs/dev_env.yaml b/devtools/conda-envs/dev_env.yaml index a0be5ce0a..27ce7bcbd 100644 --- a/devtools/conda-envs/dev_env.yaml +++ b/devtools/conda-envs/dev_env.yaml @@ -6,7 +6,7 @@ dependencies: # Core - python =3.10 - pip - - versioneer-518 + - versioningit - pip - numpy - pydantic =2 @@ -15,7 +15,7 @@ dependencies: - openff-units - ambertools =23 # Optional features - - openmm + - openmm =8.1.2 # smirnoff-plugins =2024 # de-forcefields # add back after smirnoff-plugins update - openff-nagl diff --git a/devtools/conda-envs/docs_env.yaml b/devtools/conda-envs/docs_env.yaml index 260cfc5c3..da1d6e9a2 100644 --- a/devtools/conda-envs/docs_env.yaml +++ b/devtools/conda-envs/docs_env.yaml @@ -5,11 +5,12 @@ dependencies: # Base depends # https://docs.readthedocs.io/en/stable/config-file/v2.html#build-tools-python - python =3.10 + - versioningit - pip - numpy =1 - pydantic =2 - openff-toolkit-base ~=0.16.4 - - openmm =8 + - openmm =8.1.2 - mbuild - foyer =1 - nglview @@ -20,8 +21,8 @@ dependencies: - myst-parser - numpydoc - autodoc-pydantic =2 - - sphinx ~=4.4 + - sphinx ==6.1.3 - sphinxcontrib-mermaid - sphinx-notfound-page - pip: - - git+https://github.com/openforcefield/openff-sphinx-theme.git@main + - git+https://github.com/openforcefield/openff-sphinx-theme.git@main diff --git a/devtools/conda-envs/examples_env.yaml b/devtools/conda-envs/examples_env.yaml index 482e52c34..f07c74c96 100644 --- a/devtools/conda-envs/examples_env.yaml +++ b/devtools/conda-envs/examples_env.yaml @@ -5,7 +5,7 @@ channels: dependencies: # Core - python - - versioneer-518 + - versioningit - numpy - pydantic =2 # OpenFF stack diff --git a/devtools/conda-envs/test_env.yaml b/devtools/conda-envs/test_env.yaml index a62789c9d..581ee7876 100644 --- a/devtools/conda-envs/test_env.yaml +++ b/devtools/conda-envs/test_env.yaml @@ -4,7 +4,7 @@ channels: dependencies: # Core - python - - versioneer-518 + - versioningit - numpy - pydantic =2 # OpenFF stack @@ -14,7 +14,7 @@ dependencies: # Needs to be explicitly listed to not be dropped when AmberTools is removed - rdkit # Optional features - - openmm + - openmm =8.1.2 # smirnoff-plugins =2024 # de-forcefields # add back after smirnoff-plugins update - openff-nagl diff --git a/docs/conf.py b/docs/conf.py index 69e2a9560..b026e5fba 100644 --- a/docs/conf.py +++ b/docs/conf.py @@ -14,10 +14,8 @@ # Incase the project was not installed import importlib import os -import sys - -sys.path.insert(0, os.path.abspath("..")) +from openff.interchange import __version__ # -- Project information ----------------------------------------------------- @@ -27,18 +25,12 @@ "Computational Molecular Science Python Cookiecutter version 1.2" ) author = "Open Force Field Initiative" - -# The short X.Y version -version = "" -# The full version, including alpha/beta/rc tags -release = "" +version = __version__ # -- General configuration --------------------------------------------------- # If your documentation needs a minimal Sphinx version, state it here. -# -# needs_sphinx = '1.0' # Add any Sphinx extension module names here, as strings. They can be # extensions coming with Sphinx (named 'sphinx.ext.*') or your custom @@ -89,7 +81,6 @@ autosummary_context = { # Modules to exclude from API docs "exclude_modules": [ - "openff.interchange.conftest", "openff.interchange._tests", # Maybe this is excluded by default? ], } @@ -102,7 +93,7 @@ autodoc_default_options = { "member-order": "bysource", "undoc-members": True, - "inherited-members": [], + "inherited-members": "", "show-inheritance": True, } autodoc_preserve_defaults = True @@ -169,7 +160,7 @@ # # This is also used if you do content translation via gettext catalogs. # Usually you set "language" from the command line for these cases. -language = None +language = "en" # List of patterns, relative to source directory, that match files and # directories to ignore when looking for source files. @@ -258,7 +249,7 @@ ( master_doc, "interchange.tex", - "openff-interchange Documentation", + "Interchange Documentation", "interchange", "manual", ), @@ -270,7 +261,7 @@ # One entry per manual page. List of tuples # (source start file, name, description, authors, manual section). man_pages = [ - (master_doc, "interchange", "openff-interchange Documentation", [author], 1), + (master_doc, "interchange", "Interchange Documentation", [author], 1), ] @@ -283,7 +274,7 @@ ( master_doc, "interchange", - "openff-interchange Documentation", + "Interchange Documentation", author, "interchange", "A molecular interchange object from the Open Force Field Initiative", @@ -292,4 +283,13 @@ ] -# -- Extension configuration ------------------------------------------------- +# -- ReadTheDocs configuration ------------------------------------------------- + +if os.environ.get("READTHEDOCS", "") == "True": + # Define the canonical URL if you are using a custom domain on Read the Docs + html_baseurl = os.environ.get("READTHEDOCS_CANONICAL_URL", "") + + # Tell Jinja2 templates the build is running on Read the Docs + if "html_context" not in globals(): + html_context = {} + html_context["READTHEDOCS"] = True diff --git a/docs/docutils.conf b/docs/docutils.conf new file mode 100644 index 000000000..d22dcc571 --- /dev/null +++ b/docs/docutils.conf @@ -0,0 +1,2 @@ +[html writers] +table-style: colwidths-grid diff --git a/docs/index.md b/docs/index.md index 15bc0ffb4..c9b3249e0 100644 --- a/docs/index.md +++ b/docs/index.md @@ -31,6 +31,7 @@ using/plugins.md using/status.md using/edges.md using/experimental.md +using/advanced.md ``` diff --git a/docs/installation.md b/docs/installation.md index e5ee9cceb..aad720e5a 100644 --- a/docs/installation.md +++ b/docs/installation.md @@ -1,6 +1,6 @@ # Installation -These instructions assume that the `mamba` package manager is installed. If you do not have Conda/Mamba or a drop-in replacement installed, see the [OpenFF installation documentation](openff.docs:install). +These instructions assume that the `mamba` package manager is installed. If you do not have Conda/Mamba or a drop-in replacement installed, see the [OpenFF installation documentation](inv:openff.docs#install). ## Quick Installation diff --git a/docs/releasehistory.md b/docs/releasehistory.md index 4d51f1ea8..e65f075f1 100644 --- a/docs/releasehistory.md +++ b/docs/releasehistory.md @@ -15,9 +15,23 @@ Please note that all releases prior to a version 1.0.0 are considered pre-releas ### New features -* `Interchange.from_openmm` now processes virtual sites, but only `openmm.ThreeParticleAverageSite`s. +* #1081 `Interchange.from_openmm` now processes virtual sites, but only `openmm.ThreeParticleAverageSite`s. +* #1053 Logs, at the level of `logging.INFO`, how charges are assigned by SMIRNOFF force fields to each atom and virtual site. -## 0.4.0 - 2024 +### Bug fixes + +* The `charge_from_molecules` argument must include only molecules that contain partial charges and are non-isomorphic with each other. +* The `charge_from_molecules` argument as used by the OpenFF Toolkit is handled internally as `molecules_with_preset_charges`. + +### Performance improvements + +* #1097 Migrates version handling to `versioningit`, which should result in shorter import times. + +### Documentation improvements + +* Documents charge assignment hierarchy in the user guide. + +## 0.4.0 - 2024-11-04 ### Breaking changes and behavior changes @@ -58,6 +72,7 @@ Please note that all releases prior to a version 1.0.0 are considered pre-releas ### Bug fixes +* Removes an internal use of the `@experimental` decorator which prevented `Interchange.from_openmm` from being fully removed from an experimental state. * Fixes a regression in which some `ElectrostaticsCollection.charges` properties did not return cached values. * Better process atom names in `Interchange.from_openmm` * Fixes regression tests. diff --git a/docs/using/advanced.md b/docs/using/advanced.md new file mode 100644 index 000000000..8623d02f7 --- /dev/null +++ b/docs/using/advanced.md @@ -0,0 +1,44 @@ +# Advanced features + +## Charge assignment logging + +SMIRNOFF force fields support several different partial charge assignment methods. These are applied, in the following order + +1. Look for preset charges from the `charge_from_molecules` argument +1. Look for chemical environment matches within the `` section +1. Look for chemical environment matches within the `` section +1. Try to run AM1-BCC according to the `` section or some variant + +If a molecule gets charges from one method, attempts to match charges for later methods are skipped. Note that preset charges override the force field and are not checked for consistency; any charges provided to the `charge_from_molecules` argument technically modify the force field. For more on how SMIRNOFF defines this behavior, see [this issue](https://github.com/openforcefield/standards/issues/68) and linked discussions. + +After all (mass-bearing) atoms have partial charges assigned, virtual sites are given charges by transferring charge from the atoms they are associated with ("orientation" atoms) according to the parameters of the force field. (The value of these parameters can be 0.0.) + +Given this complexity, it may be useful to track how each atom actually got charges assigned. Interchange has opt-in logging to track this behavior. This uses the [standard library `logging` module](https://docs.python.org/3/library/logging.html) at the `INFO` level. The easiest way to get started is by adding something like `logging.basicConfig(level=logging.INFO)` to the beginning of a script or program. For example, this script: + +```python +import logging + +from openff.toolkit import ForceField, Molecule + +logging.basicConfig(level=logging.INFO) + +ForceField("openff-2.2.0.offxml").create_interchange( + Molecule.from_smiles("CCO").to_topology() +) +``` + +will produce output including something like + +```shell +INFO:openff.interchange.smirnoff._nonbonded:Charge section ToolkitAM1BCC, using charge method am1bcc, applied to (topology) atom index 0 +INFO:openff.interchange.smirnoff._nonbonded:Charge section ToolkitAM1BCC, using charge method am1bcc, applied to (topology) atom index 1 +INFO:openff.interchange.smirnoff._nonbonded:Charge section ToolkitAM1BCC, using charge method am1bcc, applied to (topology) atom index 2 +INFO:openff.interchange.smirnoff._nonbonded:Charge section ToolkitAM1BCC, using charge method am1bcc, applied to (topology) atom index 3 +INFO:openff.interchange.smirnoff._nonbonded:Charge section ToolkitAM1BCC, using charge method am1bcc, applied to (topology) atom index 4 +INFO:openff.interchange.smirnoff._nonbonded:Charge section ToolkitAM1BCC, using charge method am1bcc, applied to (topology) atom index 5 +INFO:openff.interchange.smirnoff._nonbonded:Charge section ToolkitAM1BCC, using charge method am1bcc, applied to (topology) atom index 6 +INFO:openff.interchange.smirnoff._nonbonded:Charge section ToolkitAM1BCC, using charge method am1bcc, applied to (topology) atom index 7 +INFO:openff.interchange.smirnoff._nonbonded:Charge section ToolkitAM1BCC, using charge method am1bcc, applied to (topology) atom index 8 +``` + +This functionality is only available with SMIRNOFF force fields. diff --git a/docs/using/edges.md b/docs/using/edges.md index 3e2749e1e..236397abf 100644 --- a/docs/using/edges.md +++ b/docs/using/edges.md @@ -1,5 +1,35 @@ # Sharp edges +## Quirks of charge assignment + +### Charge assignment hierarchy + +Interchange, following the [SMIRNOFF specification](https://openforcefield.github.io/standards/standards/smirnoff/#partial-charge-and-electrostatics-models), assigns charges to (heavy) atoms with the following priority: + +1. **Preset charges**: Look for molecule matches in the `charge_from_molecules` argument +2. **Library charges**: Look for chemical environment matches in the `` section of the force field +3. **Charge increment models**: Look for chemical environment matches in the `` section of the force field +4. **AM1-BCC**: Try to run some variant of AM1-BCC as described by the `` section of the force field + +If charges are successfully assigned using a method, no lower-priority methods in this list are attempted. For example: + +* A force field with library charge parameters for peptides (i.e. a biopolymer force field, covering all appropriate residues) will NOT try to call AM1-BCC on a biopolymers. +* If a ligand is successfully assigned preset charges, chemical environment matching of library charges and charge increments will be skipped, as will AM1-BCC. +* If a variant of AM1-BCC (i.e. using something other than AM1 and/or using custom BCCs) is encoded in a `` section, other AM1-BCC implementations will not be called. +* If preset charges are not provided, a force field like OpenFF's Parsley or Sage lines without (many) library charges or charge increments will attempt AM1-BCC on all molecules (except water and monoatomic ions). + +After all of these steps are complete and all heavy atoms given partial charges, virtual sites are assigned charges using the values of `charge_increment`s in the virtual site parameters. Because virtual site charge are only described by the force field, using preset charges with virtual sites is discouraged. + +### Preset charges + +The charges specified by the force field can be overridden by providing molecules with partial charges to the `charge_from_molecules` argument. This may be used to make use of alternate implementations of the appropriate charge generation method, or to provide different charges to the force field. Charges provided via `charge_from_molecules` are called "preset charges" because they are pre-set by the user, rather than computed by the force field. The following restrictions are in place when using preset charges: + +* All molecules in the the `charge_from_molecules` list must be non-isomorphic with each other. +* All molecules in the the `charge_from_molecules` list must have partial charges. +* All copies of a molecule in the topology will be parametrized with the charges from an isomorphic molecule from the `charge_from_molecules` list. + +Using preset charges with virtual sites is discouraged as it can provide surprising results. + ## Quirks of core OpenFF objects Future refactors may remove the side effects of these quirks, but currently there are some diff --git a/docs/using/plugins.md b/docs/using/plugins.md index 2ace5953d..41ed8a6e4 100644 --- a/docs/using/plugins.md +++ b/docs/using/plugins.md @@ -8,7 +8,7 @@ Currently, systems using custom interactions can only be exported to OpenMM. The [`ParameterHandler`]: openff.toolkit.typing.engines.smirnoff.parameters.ParameterHandler [`SMIRNOFFCollection`]: openff.interchange.plugins.SMIRNOFFCollection -[SetupTools entry points]: setuptools:userguide/entry_point +[SetupTools entry points]: inv:setuptools#userguide/entry_point [`openmm.Force`]: openmm.openmm.Force ## Creating a custom `ParameterHandler` @@ -46,9 +46,9 @@ The high-level objectives of a parameter handler are to: * define optional attributes such as `id` [`ParameterType`]: openff.toolkit.typing.engines.smirnoff.parameters.ParameterType -[entry point group]: setuptools:userguide/entry_point +[entry point group]: inv:setuptools#userguide/entry_point [`check_handler_compatibility`]: openff.toolkit.typing.engines.smirnoff.parameters.ParameterHandler.check_handler_compatibility -[Toolkit documentation]: openff.toolkit:developing +[Toolkit documentation]: inv:openff.toolkit#developing ## Creating a custom `SMIRNOFFCollection` diff --git a/openff/interchange/__init__.py b/openff/interchange/__init__.py index 517a18a1f..3e7020b6e 100644 --- a/openff/interchange/__init__.py +++ b/openff/interchange/__init__.py @@ -1,16 +1,16 @@ """A project (and object) for storing, manipulating, and converting molecular mechanics data.""" import importlib +from importlib.metadata import version from types import ModuleType from typing import TYPE_CHECKING -from openff.interchange._version import get_versions - if TYPE_CHECKING: # Type checkers can't see lazy-imported objects from openff.interchange.components.interchange import Interchange -__version__ = get_versions()["version"] + +__version__ = version("openff.interchange") __all__: list[str] = [ "Interchange", diff --git a/openff/interchange/_tests/conftest.py b/openff/interchange/_tests/conftest.py index ce566e68c..061341aa8 100644 --- a/openff/interchange/_tests/conftest.py +++ b/openff/interchange/_tests/conftest.py @@ -91,7 +91,7 @@ def sage_with_sigma_hole(sage): smirks="[#6:1]-[#17:2]", distance=1.4 * unit.angstrom, type="BondCharge", - match="once", + match="all_permutations", charge_increment1=0.1 * unit.elementary_charge, charge_increment2=0.2 * unit.elementary_charge, ), diff --git a/openff/interchange/_tests/unit_tests/smirnoff/test_charge_assignment_logging.py b/openff/interchange/_tests/unit_tests/smirnoff/test_charge_assignment_logging.py new file mode 100644 index 000000000..3e8bedda3 --- /dev/null +++ b/openff/interchange/_tests/unit_tests/smirnoff/test_charge_assignment_logging.py @@ -0,0 +1,447 @@ +"""Test charge assignment logging.""" + +import logging +import re +from collections import defaultdict + +import pytest +from openff.toolkit import ForceField, Molecule, Topology +from openff.toolkit.utils import OPENEYE_AVAILABLE + +from openff.interchange._tests import get_protein + +""" +Hierarchy is +1. Match molceules with preset charges +2. Match chemical environment against library charges +3. Match chemical environment against charge increments +4. Run AM1-BCC (or a variant) on remaining molecules + +Test cases +---------- + +0. Sage with a ligand in vacuum +* Ligand gets AM1-BCC + +1. Sage with a ligand in water/ions +* Water gets library charges +* Ions get library charges +* Non-water solvent gets AM1-BCC +* Ligands get AM1-BCC + +2. Sage with a ligand in non-water solvent +* Non-water solvent gets AM1-BCC +* Ligands get AM1-BCC + +2. Sage with a ligand in mixed solvent +* Water gets library charges +* Ions get library charges +* Non-water solvent gets AM1-BCC +* Ligands get AM1-BCC + +4. ff14SB with Sage +* Protein gets library charges +* Water gets library charges +* Ions get library charges +* Non-water solvent gets AM1-BCC +* Ligands get AM1-BCC + +5. ff14SB with Sage and preset charges on Protein A +* Protein A gets preset charges +* Other proteins get library charges +* Water gets library charges +* Ions get library charges +* Non-water solvent gets AM1-BCC +* Ligands get AM1-BCC + +6. Sage with ligand and OPC water +* Water gets library charges +* Water virtual sites get ?? +* Ions get library charges +* Non-water solvent gets AM1-BCC +* Ligands get AM1-BCC + +7. Sage with preset charges on ligand A +* Water gets library charges +* Ions get library charges +* Non-water solvent gets AM1-BCC +* Ligand A gets preset charges +* Other ligands get AM1-BCC + +8. Sage with preset charges on water +* Water gets preset charges +* Ions get library charges +* Non-water solvent gets AM1-BCC +* Ligands get AM1-BCC + +9. Sage with (ligand) virtual site parameters +* Water gets preset charges +* Ions get library charges +* Non-water solvent gets AM1-BCC +* Ligand heavy atoms get AM1-BCC and charge increments +* Virtual sites get charge increments + +10. AM1-with-custom-BCCs Sage with ligand and ions water +* Water gets library charges +* Ions get library charges +* Ligand gets charge increments + +Other details +* Specifics of charge method (AM1-BCC, AM1-BCC-ELF10, AM1-BCC via NAGL) +* Molecules with preset charges can be similar but not exact enough +* Preset charges and virtual sites is underspecified and cannot be tested +""" + +# TODO: Handle whether or not NAGL/graph charges are applied +AM1BCC_KEY = "am1bccelf10" if OPENEYE_AVAILABLE else "am1bcc" +NAGL_KEY = "openff-gnn-am1bcc-0.1.0-rc.2.pt" + + +def map_methods_to_atom_indices(caplog: pytest.LogCaptureFixture) -> dict[str, list[int]]: + """ + Map partial charge assignment methods to (sorted) atom indices. + """ + info = defaultdict(list) + + for record in caplog.records: + # skip logged warnings from upstreams/other packages + if record.name.startswith("openff.interchange"): + assert record.levelname == "INFO", "Only INFO logs are expected." + else: + continue + + message = record.msg + + if message.startswith("Charge section LibraryCharges"): + info["library"].append(int(message.split("atom index ")[-1])) + + elif message.startswith("Charge section ToolkitAM1BCC"): + info[message.split(", ")[1].split(" ")[-1]].append(int(message.split("atom index ")[-1])) + + # without also pulling the virtual site - particle mapping (which is different for each engine) + # it's hard to store more information than the orientation atoms that are affected by each + # virtual site's charges + elif message.startswith("Charge section VirtualSites"): + orientation_atoms: list[int] = [ + int(value.strip()) for value in re.findall(r"\((.*?)\)", message)[0].split(",") + ] + + for atom in orientation_atoms: + info["orientation"].append(atom) + + elif message.startswith("Preset charges"): + info["preset"].append(int(message.split("atom index")[-1])) + + elif message.startswith("Charge section ChargeIncrementModel"): + if "using charge method" in message: + info[f"chargeincrements_{message.split(',')[1].split(' ')[-1]}"].append( + int(message.split("atom index ")[-1]), + ) + + elif "applying charge increment" in message: + # TODO: Store the "other" atoms? + info["chargeincrements"].append(int(message.split("atom ")[1].split(" ")[0])) + + else: + raise ValueError(f"Unexpected log message {message}") + + return {key: sorted(val) for key, val in info.items()} + + +@pytest.fixture +def sage_with_nagl(sage): + from openff.toolkit.typing.engines.smirnoff.parameters import ChargeIncrementModelHandler, ChargeIncrementType + + sage.register_parameter_handler( + parameter_handler=ChargeIncrementModelHandler( + version=0.3, + partial_charge_method=NAGL_KEY, + ), + ) + + # Add dummy "BCCs" for testing, even though this model has BCCs built into + # the partial charge assignment method + sage["ChargeIncrementModel"].add_parameter( + parameter=ChargeIncrementType( + smirks="[#6:1]-[#1:2]", + charge_increment=[ + "0.1 elementary_charge", + "-0.1 elementary_charge", + ], + ), + ) + + return sage + + +@pytest.fixture +def ions() -> Topology: + return Topology.from_molecules( + [ + Molecule.from_smiles(smiles) + for smiles in [ + "[Na+]", + "[Cl-]", + "[Na+]", + "[Cl-]", + ] + ], + ) + + +@pytest.fixture +def ligand() -> Molecule: + return Molecule.from_mapped_smiles("[C:1]([C:2]([O:3][H:9])([H:7])[H:8])([H:4])([H:5])[H:6]") + + +@pytest.fixture +def solvent() -> Molecule: + return Molecule.from_mapped_smiles("[H:3][C:1]([H:4])([H:5])[N:2]([H:6])[H:7]") + + +def ligand_in_solvent(ligand, solvent) -> Topology: + return Topology.from_molecules( + [ + ligand, + solvent, + solvent, + solvent, + ], + ) + + +@pytest.fixture +def water_and_ions(water, ions) -> Topology: + topology = Topology.from_molecules(3 * [water]) + for molecule in ions.molecules: + topology.add_molecule(molecule) + + return topology + + +@pytest.fixture +def ligand_and_water_and_ions(ligand, water_and_ions) -> Topology: + topology = ligand.to_topology() + + for molecule in water_and_ions.molecules: + topology.add_molecule(molecule) + + return topology + + +""" +0.xSage with just a ligand +1.xSage with a ligand in water/ions +2.xSage with a ligand in non-water solvent +3.xSage with a ligand in mixed solvent +4.xff14SB with Sage +5.xff14SB with Sage and preset charges on Protein A +6.xSage with ligand and OPC water +7.xSage with preset charges on ligand A +8.xSage with preset charges on water +9.xSage with (ligand) virtual site parameters +10.xAM1-with-custom-BCCs Sage with ligand and ions water +""" + + +def test_case0(caplog, sage, ligand): + with caplog.at_level(logging.INFO): + sage.create_interchange(ligand.to_topology()) + + info = map_methods_to_atom_indices(caplog) + + # atoms 0 through 8 are ethanol, getting AM1-BCC + assert info[AM1BCC_KEY] == [*range(0, 9)] + + +def test_case1(caplog, sage, ligand_and_water_and_ions): + with caplog.at_level(logging.INFO): + sage.create_interchange(ligand_and_water_and_ions) + + info = map_methods_to_atom_indices(caplog) + + # atoms 0 through 8 are ethanol, getting AM1-BCC + assert info[AM1BCC_KEY] == [*range(0, 9)] + + # atoms 9 through 21 are water + ions, getting library charges + assert info["library"] == [*range(9, 22)] + + +def test_case2(caplog, sage, ligand, solvent): + topology = Topology.from_molecules([ligand, solvent, solvent, solvent]) + + with caplog.at_level(logging.INFO): + sage.create_interchange(topology) + + info = map_methods_to_atom_indices(caplog) + + # everything should get AM1-BCC charges + assert info[AM1BCC_KEY] == [*range(0, topology.n_atoms)] + + +def test_case3(caplog, sage, ligand_and_water_and_ions, solvent): + for index in range(3): + ligand_and_water_and_ions.add_molecule(solvent) + + ligand_and_water_and_ions.molecule(0).assign_partial_charges(partial_charge_method="gasteiger") + + with caplog.at_level(logging.INFO): + sage.create_interchange( + ligand_and_water_and_ions, + ) + + info = map_methods_to_atom_indices(caplog) + + # atoms 0 through 8 are ethanol, getting AM1-BCC, + # and also solvent molecules (starting at index 22) + assert info[AM1BCC_KEY] == [*range(0, 9), *range(22, 22 + 3 * 7)] + + # atoms 9 through 21 are water + ions, getting library charges + assert info["library"] == [*range(9, 22)] + + +@pytest.mark.slow +@pytest.mark.parametrize("preset_on_protein", [True, False]) +def test_cases4_5(caplog, ligand_and_water_and_ions, preset_on_protein): + ff = ForceField("ff14sb_off_impropers_0.0.3.offxml", "openff-2.0.0.offxml") + + complex = get_protein("MainChain_ALA_ALA").to_topology() + + for molecule in ligand_and_water_and_ions.molecules: + complex.add_molecule(molecule) + + if preset_on_protein: + complex.molecule(0).assign_partial_charges(partial_charge_method="zeros") + + with caplog.at_level(logging.INFO): + if preset_on_protein: + ff.create_interchange(complex, charge_from_molecules=[complex.molecule(0)]) + else: + ff.create_interchange(complex) + + info = map_methods_to_atom_indices(caplog) + + assert info[AM1BCC_KEY] == [*range(complex.molecule(0).n_atoms, complex.molecule(0).n_atoms + 9)] + + if preset_on_protein: + # protein gets preset charges + assert info["preset"] == [*range(0, complex.molecule(0).n_atoms)] + + # everything after the protein and ligand should get library charges + assert info["library"] == [ + *range(complex.molecule(0).n_atoms + 9, complex.n_atoms), + ] + else: + # the protein and everything after the ligand should get library charges + assert info["library"] == [ + *range(0, complex.molecule(0).n_atoms), + *range(complex.molecule(0).n_atoms + 9, complex.n_atoms), + ] + + +def test_case6(caplog, ligand, water): + force_field = ForceField("openff-2.0.0.offxml", "opc.offxml") + + topology = Topology.from_molecules([ligand, water, water, water]) + + with caplog.at_level(logging.INFO): + force_field.create_interchange(topology) + + info = map_methods_to_atom_indices(caplog) + + # atoms 0 through 8 are ethanol, getting AM1-BCC charges + assert info[AM1BCC_KEY] == [*range(0, 9)] + + # atoms 9 through 17 are water atoms, getting library charges + assert info["library"] == [*range(9, 18)] + + # particles 18 through 20 are water virtual sites, but the current logging strategy + # makes it hard to match these up (and the particle indices are different OpenMM/GROMACS/etc) + + # can still check that orientation atoms are subject to virtual site + # charge increments (even if the increment is +0.0 e) + assert info["orientation"] == [*range(9, 18)] + + +def test_case7(caplog, sage, ligand_and_water_and_ions): + ligand_and_water_and_ions.molecule(0).assign_partial_charges(partial_charge_method="gasteiger") + + with caplog.at_level(logging.INFO): + sage.create_interchange( + ligand_and_water_and_ions, + charge_from_molecules=[ligand_and_water_and_ions.molecule(0)], + ) + + info = map_methods_to_atom_indices(caplog) + + # atoms 0 through 8 are ethanol, getting preset charges + assert info["preset"] == [*range(0, 9)] + + # atoms 9 through 21 are water + ions, getting library charges + assert info["library"] == [*range(9, 22)] + + +def test_case8(caplog, sage, water_and_ions): + water_and_ions.molecule(0).assign_partial_charges(partial_charge_method="gasteiger") + + with caplog.at_level(logging.INFO): + sage.create_interchange( + water_and_ions, + charge_from_molecules=[water_and_ions.molecule(0)], + ) + + info = map_methods_to_atom_indices(caplog) + + # atoms 0 through 8 are water, getting preset charges + assert info["preset"] == [*range(0, 9)] + + # atoms 9 through 12 are ions, getting library charges + assert info["library"] == [*range(9, 13)] + + +def test_case9(caplog, sage_with_bond_charge): + with caplog.at_level(logging.INFO): + sage_with_bond_charge.create_interchange( + Molecule.from_mapped_smiles( + "[H:3][C:1]([H:4])([H:5])[Cl:2]", + ).to_topology(), + ) + + info = map_methods_to_atom_indices(caplog) + + # atoms 0 through 5 are ligand, getting AM1-BCC charges + assert info[AM1BCC_KEY] == [*range(0, 5)] + + # atoms 0 and 1 are the orientation atoms of the sigma hole virtual site + assert info["orientation"] == [0, 1] + + +def test_case10(caplog, sage_with_nagl, ligand): + from openff.toolkit.utils.nagl_wrapper import NAGLToolkitWrapper + from openff.toolkit.utils.rdkit_wrapper import RDKitToolkitWrapper + from openff.toolkit.utils.toolkit_registry import ToolkitRegistry, toolkit_registry_manager + + pytest.importorskip("openff.nagl") + pytest.importorskip("rdkit") + + with caplog.at_level(logging.INFO): + with toolkit_registry_manager( + toolkit_registry=ToolkitRegistry( + toolkit_precedence=[NAGLToolkitWrapper, RDKitToolkitWrapper], + ), + ): + sage_with_nagl.create_interchange(ligand.to_topology()) + + info = map_methods_to_atom_indices(caplog) + + # atoms 0 through 5 are ligand, getting NAGL charges + assert info[f"chargeincrements_{NAGL_KEY}"] == [*range(0, 9)] + + # TODO: These are logged symmetrically (i.e. hydrogens are listed) + # even though the charges appear to be correct, assert should + # simply by == [0, 1] since the hydrogens shouldn't be caught + assert 0 in info["chargeincrements"] + assert 1 in info["chargeincrements"] + + # the standard AM1-BCC should not have ran + assert AM1BCC_KEY not in info diff --git a/openff/interchange/_tests/unit_tests/smirnoff/test_create.py b/openff/interchange/_tests/unit_tests/smirnoff/test_create.py index 58e70516c..c2cf1fcc4 100644 --- a/openff/interchange/_tests/unit_tests/smirnoff/test_create.py +++ b/openff/interchange/_tests/unit_tests/smirnoff/test_create.py @@ -7,6 +7,7 @@ from openff.interchange._tests import get_test_file_path from openff.interchange.constants import kcal_mol, kcal_mol_a2 from openff.interchange.exceptions import ( + PresetChargesError, UnassignedAngleError, UnassignedBondError, UnassignedTorsionError, @@ -21,6 +22,7 @@ _upconvert_vdw_handler, library_charge_from_molecule, ) +from openff.interchange.warnings import PresetChargesAndVirtualSitesWarning def _get_interpolated_bond_k(bond_handler) -> float: @@ -160,7 +162,8 @@ def test_library_charges_from_molecule(): assert library_charges.charge == [*mol.partial_charges] -class TestChargeFromMolecules: +class TestPresetCharges: + @pytest.mark.slow def test_charge_from_molecules_basic(self, sage): molecule = Molecule.from_smiles("CCO") molecule.assign_partial_charges(partial_charge_method="am1bcc") @@ -231,6 +234,102 @@ def test_charges_from_molecule_reordered( assert numpy.allclose(expected_charges, found_charges) + def test_error_when_any_missing_partial_charges(self, sage): + ethanol = Molecule.from_smiles("CCO") + + ethanol.assign_partial_charges(partial_charge_method="am1bcc") + + topology = Topology.from_molecules( + [ + Molecule.from_smiles("C"), + ethanol, + ], + ) + + with pytest.raises( + PresetChargesError, + match="All.*must have partial charges", + ): + sage.create_interchange( + topology, + charge_from_molecules=[Molecule.from_smiles("C")], + ) + + def test_error_multiple_isomorphic_molecules(self, sage, ethanol, reversed_ethanol): + # these ethanol fixtures *do* have charges, if they didn't have charges it would be + # ambiguous which error should be raised + topology = Topology.from_molecules([ethanol, reversed_ethanol, ethanol]) + + with pytest.raises( + PresetChargesError, + match="All molecules.*isomorphic.*unique", + ): + sage.create_interchange( + topology, + charge_from_molecules=[ethanol, reversed_ethanol], + ) + + def test_virtual_site_charge_increments_applied_after_preset_charges_water( + self, + tip4p, + water, + ): + """ + Test that charge increments to/from virtual sites are still applied after preset charges are set. + + This is funky user input, since preset charges override library charges, which are important for water. + """ + water.partial_charges = Quantity( + [-2.0, 1.0, 1.0], + "elementary_charge", + ) + + with pytest.warns(PresetChargesAndVirtualSitesWarning): + # This dict is probably ordered O H H VS + charges = tip4p.create_interchange( + water.to_topology(), + charge_from_molecules=[water], + )["Electrostatics"].charges + + # tip4p_fb.offxml is meant to result in charges of -1.0517362213526 (VS) 0 (O) and 0.5258681106763 (H) + # the force field has 0.0 for all library charges (skipping AM1-BCC), so preset charges screw this up. + # the resulting charges, if using preset charges of [-2, 1, 1], should be the result of adding together + # [-2, 1, 1, 0] (preset charges) and + # [0, 1.5258681106763001, 1.5258681106763001, -1.0517362213526] (charge increments) + + numpy.testing.assert_allclose( + [charge.m for charge in charges.values()], + [-2.0, 1.5258681106763001, 1.5258681106763001, -1.0517362213526], + ) + + def test_virtual_site_charge_increments_applied_after_preset_charges_ligand( + self, + sage_with_sigma_hole, + ): + ligand = Molecule.from_mapped_smiles("[C:1]([Cl:2])([Cl:3])([Cl:4])[Cl:5]") + + ligand.partial_charges = Quantity( + [1.2, -0.3, -0.3, -0.3, -0.3], + "elementary_charge", + ) + + with pytest.warns(PresetChargesAndVirtualSitesWarning): + # This dict is probably ordered C Cl Cl Cl Cl VS VS VS VS + charges = sage_with_sigma_hole.create_interchange( + ligand.to_topology(), + charge_from_molecules=[ligand], + )["Electrostatics"].charges + + # this fake sigma hole parameter shifts 0.1 e from the carbon and 0.2 e from the chlorine, so the + # resulting charges should be like + # [1.2, -0.3, -0.3, -0.3, -0.3] + + # [0.4, 0.2, 0.2, 0.2, 0.2, -0.3, -0.3, -0.3, -0.3] + + numpy.testing.assert_allclose( + [charge.m for charge in charges.values()], + [1.6, -0.1, -0.1, -0.1, -0.1, -0.3, -0.3, -0.3, -0.3], + ) + @skip_if_missing("openmm") class TestPartialBondOrdersFromMolecules: diff --git a/openff/interchange/_version.py b/openff/interchange/_version.py deleted file mode 100644 index 72757a8bf..000000000 --- a/openff/interchange/_version.py +++ /dev/null @@ -1,735 +0,0 @@ -# This file helps to compute a version number in source trees obtained from -# git-archive tarball (such as those provided by githubs download-from-tag -# feature). Distribution tarballs (built by setup.py sdist) and build -# directories (produced by setup.py build) will contain a much shorter file -# that just contains the computed version number. - -# This file is released into the public domain. -# Generated by versioneer-0.29 -# https://github.com/python-versioneer/python-versioneer - -"""Git implementation of _version.py.""" - -import errno -import functools -import os -import re -import subprocess -import sys -from collections.abc import Callable -from typing import Any - - -def get_keywords() -> dict[str, str]: - """Get the keywords needed to look up the version information.""" - # these strings will be replaced by git during git-archive. - # setup.py/versioneer.py will grep for the variable names, so they must - # each be defined on a line of their own. _version.py will just call - # get_keywords(). - git_refnames = "$Format:%d$" - git_full = "$Format:%H$" - git_date = "$Format:%ci$" - keywords = {"refnames": git_refnames, "full": git_full, "date": git_date} - return keywords - - -class VersioneerConfig: - """Container for Versioneer configuration parameters.""" - - VCS: str - style: str - tag_prefix: str - parentdir_prefix: str - versionfile_source: str - verbose: bool - - -def get_config() -> VersioneerConfig: - """Create, populate and return the VersioneerConfig() object.""" - # these strings are filled in when 'setup.py versioneer' creates - # _version.py - cfg = VersioneerConfig() - cfg.VCS = "git" - cfg.style = "pep440" - cfg.tag_prefix = "v" - cfg.parentdir_prefix = "None" - cfg.versionfile_source = "openff/interchange/_version.py" - cfg.verbose = False - return cfg - - -class NotThisMethod(Exception): - """Exception raised if a method is not valid for the current scenario.""" - - -LONG_VERSION_PY: dict[str, str] = {} -HANDLERS: dict[str, dict[str, Callable]] = {} - - -def register_vcs_handler(vcs: str, method: str) -> Callable: # decorator - """Create decorator to mark a method as the handler of a VCS.""" - - def decorate(f: Callable) -> Callable: - """Store f in HANDLERS[vcs][method].""" - if vcs not in HANDLERS: - HANDLERS[vcs] = {} - HANDLERS[vcs][method] = f - return f - - return decorate - - -def run_command( - commands: list[str], - args: list[str], - cwd: str | None = None, - verbose: bool = False, - hide_stderr: bool = False, - env: dict[str, str] | None = None, -) -> tuple[str | None, int | None]: - """Call the given command(s).""" - assert isinstance(commands, list) - process = None - - popen_kwargs: dict[str, Any] = {} - if sys.platform == "win32": - # This hides the console window if pythonw.exe is used - startupinfo = subprocess.STARTUPINFO() - startupinfo.dwFlags |= subprocess.STARTF_USESHOWWINDOW - popen_kwargs["startupinfo"] = startupinfo - - for command in commands: - try: - dispcmd = str([command, *args]) - # remember shell=False, so use git.cmd on windows, not just git - process = subprocess.Popen( - [command, *args], - cwd=cwd, - env=env, - stdout=subprocess.PIPE, - stderr=(subprocess.PIPE if hide_stderr else None), - **popen_kwargs, - ) - break - except OSError as e: - if e.errno == errno.ENOENT: - continue - if verbose: - print(f"unable to run {dispcmd}") - print(e) - return None, None - else: - if verbose: - print(f"unable to find command, tried {commands}") - return None, None - stdout = process.communicate()[0].strip().decode() - if process.returncode != 0: - if verbose: - print(f"unable to run {dispcmd} (error)") - print(f"stdout was {stdout}") - return None, process.returncode - return stdout, process.returncode - - -def versions_from_parentdir( - parentdir_prefix: str, - root: str, - verbose: bool, -) -> dict[str, Any]: - """Try to determine the version from the parent directory name. - - Source tarballs conventionally unpack into a directory that includes both - the project name and a version string. We will also support searching up - two directory levels for an appropriately named parent directory - """ - rootdirs = [] - - for _ in range(3): - dirname = os.path.basename(root) - if dirname.startswith(parentdir_prefix): - return { - "version": dirname[len(parentdir_prefix) :], - "full-revisionid": None, - "dirty": False, - "error": None, - "date": None, - } - rootdirs.append(root) - root = os.path.dirname(root) # up a level - - if verbose: - print( - f"Tried directories {rootdirs!s} but none started with prefix {parentdir_prefix}", - ) - raise NotThisMethod("rootdir doesn't start with parentdir_prefix") - - -@register_vcs_handler("git", "get_keywords") -def git_get_keywords(versionfile_abs: str) -> dict[str, str]: - """Extract version information from the given file.""" - # the code embedded in _version.py can just fetch the value of these - # keywords. When used from setup.py, we don't want to import _version.py, - # so we do it with a regexp instead. This function is not used from - # _version.py. - keywords: dict[str, str] = {} - try: - with open(versionfile_abs) as fobj: - for line in fobj: - if line.strip().startswith("git_refnames ="): - mo = re.search(r'=\s*"(.*)"', line) - if mo: - keywords["refnames"] = mo.group(1) - if line.strip().startswith("git_full ="): - mo = re.search(r'=\s*"(.*)"', line) - if mo: - keywords["full"] = mo.group(1) - if line.strip().startswith("git_date ="): - mo = re.search(r'=\s*"(.*)"', line) - if mo: - keywords["date"] = mo.group(1) - except OSError: - pass - return keywords - - -@register_vcs_handler("git", "keywords") -def git_versions_from_keywords( - keywords: dict[str, str], - tag_prefix: str, - verbose: bool, -) -> dict[str, Any]: - """Get version information from git keywords.""" - if "refnames" not in keywords: - raise NotThisMethod("Short version file found") - date = keywords.get("date") - if date is not None: - # Use only the last line. Previous lines may contain GPG signature - # information. - date = date.splitlines()[-1] - - # git-2.2.0 added "%cI", which expands to an ISO-8601 -compliant - # datestamp. However we prefer "%ci" (which expands to an "ISO-8601 - # -like" string, which we must then edit to make compliant), because - # it's been around since git-1.5.3, and it's too difficult to - # discover which version we're using, or to work around using an - # older one. - date = date.strip().replace(" ", "T", 1).replace(" ", "", 1) - refnames = keywords["refnames"].strip() - if refnames.startswith("$Format"): - if verbose: - print("keywords are unexpanded, not using") - raise NotThisMethod("unexpanded keywords, not a git-archive tarball") - refs = {r.strip() for r in refnames.strip("()").split(",")} - # starting in git-1.8.3, tags are listed as "tag: foo-1.0" instead of - # just "foo-1.0". If we see a "tag: " prefix, prefer those. - TAG = "tag: " - tags = {r[len(TAG) :] for r in refs if r.startswith(TAG)} - if not tags: - # Either we're using git < 1.8.3, or there really are no tags. We use - # a heuristic: assume all version tags have a digit. The old git %d - # expansion behaves like git log --decorate=short and strips out the - # refs/heads/ and refs/tags/ prefixes that would let us distinguish - # between branches and tags. By ignoring refnames without digits, we - # filter out many common branch names like "release" and - # "stabilization", as well as "HEAD" and "master". - tags = {r for r in refs if re.search(r"\d", r)} - if verbose: - print("discarding '{}', no digits".format(",".join(refs - tags))) - if verbose: - print("likely tags: {}".format(",".join(sorted(tags)))) - for ref in sorted(tags): - # sorting will prefer e.g. "2.0" over "2.0rc1" - if ref.startswith(tag_prefix): - r = ref[len(tag_prefix) :] - # Filter out refs that exactly match prefix or that don't start - # with a number once the prefix is stripped (mostly a concern - # when prefix is '') - if not re.match(r"\d", r): - continue - if verbose: - print(f"picking {r}") - return { - "version": r, - "full-revisionid": keywords["full"].strip(), - "dirty": False, - "error": None, - "date": date, - } - # no suitable tags, so version is "0+unknown", but full hex is still there - if verbose: - print("no suitable tags, using unknown + full revision id") - return { - "version": "0+unknown", - "full-revisionid": keywords["full"].strip(), - "dirty": False, - "error": "no suitable tags", - "date": None, - } - - -@register_vcs_handler("git", "pieces_from_vcs") -def git_pieces_from_vcs( - tag_prefix: str, - root: str, - verbose: bool, - runner: Callable = run_command, -) -> dict[str, Any]: - """Get version from 'git describe' in the root of the source tree. - - This only gets called if the git-archive 'subst' keywords were *not* - expanded, and _version.py hasn't already been rewritten with a short - version string, meaning we're inside a checked out source tree. - """ - GITS = ["git"] - if sys.platform == "win32": - GITS = ["git.cmd", "git.exe"] - - # GIT_DIR can interfere with correct operation of Versioneer. - # It may be intended to be passed to the Versioneer-versioned project, - # but that should not change where we get our version from. - env = os.environ.copy() - env.pop("GIT_DIR", None) - runner = functools.partial(runner, env=env) - - _, rc = runner( - GITS, - ["rev-parse", "--git-dir"], - cwd=root, - hide_stderr=not verbose, - ) - if rc != 0: - if verbose: - print(f"Directory {root} not under git control") - raise NotThisMethod("'git rev-parse --git-dir' returned error") - - # if there is a tag matching tag_prefix, this yields TAG-NUM-gHEX[-dirty] - # if there isn't one, this yields HEX[-dirty] (no NUM) - describe_out, rc = runner( - GITS, - [ - "describe", - "--tags", - "--dirty", - "--always", - "--long", - "--match", - f"{tag_prefix}[[:digit:]]*", - ], - cwd=root, - ) - # --long was added in git-1.5.5 - if describe_out is None: - raise NotThisMethod("'git describe' failed") - describe_out = describe_out.strip() - full_out, rc = runner(GITS, ["rev-parse", "HEAD"], cwd=root) - if full_out is None: - raise NotThisMethod("'git rev-parse' failed") - full_out = full_out.strip() - - pieces: dict[str, Any] = {} - pieces["long"] = full_out - pieces["short"] = full_out[:7] # maybe improved later - pieces["error"] = None - - branch_name, rc = runner( - GITS, - ["rev-parse", "--abbrev-ref", "HEAD"], - cwd=root, - ) - # --abbrev-ref was added in git-1.6.3 - if rc != 0 or branch_name is None: - raise NotThisMethod("'git rev-parse --abbrev-ref' returned error") - branch_name = branch_name.strip() - - if branch_name == "HEAD": - # If we aren't exactly on a branch, pick a branch which represents - # the current commit. If all else fails, we are on a branchless - # commit. - branches, rc = runner(GITS, ["branch", "--contains"], cwd=root) - # --contains was added in git-1.5.4 - if rc != 0 or branches is None: - raise NotThisMethod("'git branch --contains' returned error") - branches = branches.split("\n") - - # Remove the first line if we're running detached - if "(" in branches[0]: - branches.pop(0) - - # Strip off the leading "* " from the list of branches. - branches = [branch[2:] for branch in branches] - if "master" in branches: - branch_name = "master" - elif not branches: - branch_name = None - else: - # Pick the first branch that is returned. Good or bad. - branch_name = branches[0] - - pieces["branch"] = branch_name - - # parse describe_out. It will be like TAG-NUM-gHEX[-dirty] or HEX[-dirty] - # TAG might have hyphens. - git_describe = describe_out - - # look for -dirty suffix - dirty = git_describe.endswith("-dirty") - pieces["dirty"] = dirty - if dirty: - git_describe = git_describe[: git_describe.rindex("-dirty")] - - # now we have TAG-NUM-gHEX or HEX - - if "-" in git_describe: - # TAG-NUM-gHEX - mo = re.search(r"^(.+)-(\d+)-g([0-9a-f]+)$", git_describe) - if not mo: - # unparsable. Maybe git-describe is misbehaving? - pieces["error"] = f"unable to parse git-describe output: '{describe_out}'" - return pieces - - # tag - full_tag = mo.group(1) - if not full_tag.startswith(tag_prefix): - if verbose: - fmt = "tag '%s' doesn't start with prefix '%s'" - print(fmt % (full_tag, tag_prefix)) - pieces["error"] = f"tag '{full_tag}' doesn't start with prefix '{tag_prefix}'" - return pieces - pieces["closest-tag"] = full_tag[len(tag_prefix) :] - - # distance: number of commits since tag - pieces["distance"] = int(mo.group(2)) - - # commit: short hex revision ID - pieces["short"] = mo.group(3) - - else: - # HEX: no tags - pieces["closest-tag"] = None - out, rc = runner(GITS, ["rev-list", "HEAD", "--left-right"], cwd=root) - pieces["distance"] = len(out.split()) # total number of commits - - # commit date: see ISO-8601 comment in git_versions_from_keywords() - date = runner(GITS, ["show", "-s", "--format=%ci", "HEAD"], cwd=root)[0].strip() - # Use only the last line. Previous lines may contain GPG signature - # information. - date = date.splitlines()[-1] - pieces["date"] = date.strip().replace(" ", "T", 1).replace(" ", "", 1) - - return pieces - - -def plus_or_dot(pieces: dict[str, Any]) -> str: - """Return a + if we don't already have one, else return a .""" - if "+" in pieces.get("closest-tag", ""): - return "." - return "+" - - -def render_pep440(pieces: dict[str, Any]) -> str: - """Build up version string, with post-release "local version identifier". - - Our goal: TAG[+DISTANCE.gHEX[.dirty]] . Note that if you - get a tagged build and then dirty it, you'll get TAG+0.gHEX.dirty - - Exceptions: - 1: no tags. git_describe was just HEX. 0+untagged.DISTANCE.gHEX[.dirty] - """ - if pieces["closest-tag"]: - rendered = pieces["closest-tag"] - if pieces["distance"] or pieces["dirty"]: - rendered += plus_or_dot(pieces) - rendered += "%d.g%s" % (pieces["distance"], pieces["short"]) - if pieces["dirty"]: - rendered += ".dirty" - else: - # exception #1 - rendered = "0+untagged.%d.g%s" % ( - pieces["distance"], - pieces["short"], - ) - if pieces["dirty"]: - rendered += ".dirty" - return rendered - - -def render_pep440_branch(pieces: dict[str, Any]) -> str: - """TAG[[.dev0]+DISTANCE.gHEX[.dirty]] . - - The ".dev0" means not master branch. Note that .dev0 sorts backwards - (a feature branch will appear "older" than the master branch). - - Exceptions: - 1: no tags. 0[.dev0]+untagged.DISTANCE.gHEX[.dirty] - """ - if pieces["closest-tag"]: - rendered = pieces["closest-tag"] - if pieces["distance"] or pieces["dirty"]: - if pieces["branch"] != "master": - rendered += ".dev0" - rendered += plus_or_dot(pieces) - rendered += "%d.g%s" % (pieces["distance"], pieces["short"]) - if pieces["dirty"]: - rendered += ".dirty" - else: - # exception #1 - rendered = "0" - if pieces["branch"] != "master": - rendered += ".dev0" - rendered += "+untagged.%d.g%s" % ( - pieces["distance"], - pieces["short"], - ) - if pieces["dirty"]: - rendered += ".dirty" - return rendered - - -def pep440_split_post(ver: str) -> tuple[str, int | None]: - """Split pep440 version string at the post-release segment. - - Returns the release segments before the post-release and the - post-release version number (or -1 if no post-release segment is present). - """ - vc = str.split(ver, ".post") - return vc[0], int(vc[1] or 0) if len(vc) == 2 else None - - -def render_pep440_pre(pieces: dict[str, Any]) -> str: - """TAG[.postN.devDISTANCE] -- No -dirty. - - Exceptions: - 1: no tags. 0.post0.devDISTANCE - """ - if pieces["closest-tag"]: - if pieces["distance"]: - # update the post release segment - tag_version, post_version = pep440_split_post(pieces["closest-tag"]) - rendered = tag_version - if post_version is not None: - rendered += ".post%d.dev%d" % (post_version + 1, pieces["distance"]) - else: - rendered += ".post0.dev%d" % (pieces["distance"]) - else: - # no commits, use the tag as the version - rendered = pieces["closest-tag"] - else: - # exception #1 - rendered = "0.post0.dev%d" % pieces["distance"] - return rendered - - -def render_pep440_post(pieces: dict[str, Any]) -> str: - """TAG[.postDISTANCE[.dev0]+gHEX] . - - The ".dev0" means dirty. Note that .dev0 sorts backwards - (a dirty tree will appear "older" than the corresponding clean one), - but you shouldn't be releasing software with -dirty anyways. - - Exceptions: - 1: no tags. 0.postDISTANCE[.dev0] - """ - if pieces["closest-tag"]: - rendered = pieces["closest-tag"] - if pieces["distance"] or pieces["dirty"]: - rendered += ".post%d" % pieces["distance"] - if pieces["dirty"]: - rendered += ".dev0" - rendered += plus_or_dot(pieces) - rendered += "g{}".format(pieces["short"]) - else: - # exception #1 - rendered = "0.post%d" % pieces["distance"] - if pieces["dirty"]: - rendered += ".dev0" - rendered += "+g{}".format(pieces["short"]) - return rendered - - -def render_pep440_post_branch(pieces: dict[str, Any]) -> str: - """TAG[.postDISTANCE[.dev0]+gHEX[.dirty]] . - - The ".dev0" means not master branch. - - Exceptions: - 1: no tags. 0.postDISTANCE[.dev0]+gHEX[.dirty] - """ - if pieces["closest-tag"]: - rendered = pieces["closest-tag"] - if pieces["distance"] or pieces["dirty"]: - rendered += ".post%d" % pieces["distance"] - if pieces["branch"] != "master": - rendered += ".dev0" - rendered += plus_or_dot(pieces) - rendered += "g{}".format(pieces["short"]) - if pieces["dirty"]: - rendered += ".dirty" - else: - # exception #1 - rendered = "0.post%d" % pieces["distance"] - if pieces["branch"] != "master": - rendered += ".dev0" - rendered += "+g{}".format(pieces["short"]) - if pieces["dirty"]: - rendered += ".dirty" - return rendered - - -def render_pep440_old(pieces: dict[str, Any]) -> str: - """TAG[.postDISTANCE[.dev0]] . - - The ".dev0" means dirty. - - Exceptions: - 1: no tags. 0.postDISTANCE[.dev0] - """ - if pieces["closest-tag"]: - rendered = pieces["closest-tag"] - if pieces["distance"] or pieces["dirty"]: - rendered += ".post%d" % pieces["distance"] - if pieces["dirty"]: - rendered += ".dev0" - else: - # exception #1 - rendered = "0.post%d" % pieces["distance"] - if pieces["dirty"]: - rendered += ".dev0" - return rendered - - -def render_git_describe(pieces: dict[str, Any]) -> str: - """TAG[-DISTANCE-gHEX][-dirty]. - - Like 'git describe --tags --dirty --always'. - - Exceptions: - 1: no tags. HEX[-dirty] (note: no 'g' prefix) - """ - if pieces["closest-tag"]: - rendered = pieces["closest-tag"] - if pieces["distance"]: - rendered += "-%d-g%s" % (pieces["distance"], pieces["short"]) - else: - # exception #1 - rendered = pieces["short"] - if pieces["dirty"]: - rendered += "-dirty" - return rendered - - -def render_git_describe_long(pieces: dict[str, Any]) -> str: - """TAG-DISTANCE-gHEX[-dirty]. - - Like 'git describe --tags --dirty --always -long'. - The distance/hash is unconditional. - - Exceptions: - 1: no tags. HEX[-dirty] (note: no 'g' prefix) - """ - if pieces["closest-tag"]: - rendered = pieces["closest-tag"] - rendered += "-%d-g%s" % (pieces["distance"], pieces["short"]) - else: - # exception #1 - rendered = pieces["short"] - if pieces["dirty"]: - rendered += "-dirty" - return rendered - - -def render(pieces: dict[str, Any], style: str) -> dict[str, Any]: - """Render the given version pieces into the requested style.""" - if pieces["error"]: - return { - "version": "unknown", - "full-revisionid": pieces.get("long"), - "dirty": None, - "error": pieces["error"], - "date": None, - } - - if not style or style == "default": - style = "pep440" # the default - - if style == "pep440": - rendered = render_pep440(pieces) - elif style == "pep440-branch": - rendered = render_pep440_branch(pieces) - elif style == "pep440-pre": - rendered = render_pep440_pre(pieces) - elif style == "pep440-post": - rendered = render_pep440_post(pieces) - elif style == "pep440-post-branch": - rendered = render_pep440_post_branch(pieces) - elif style == "pep440-old": - rendered = render_pep440_old(pieces) - elif style == "git-describe": - rendered = render_git_describe(pieces) - elif style == "git-describe-long": - rendered = render_git_describe_long(pieces) - else: - raise ValueError(f"unknown style '{style}'") - - return { - "version": rendered, - "full-revisionid": pieces["long"], - "dirty": pieces["dirty"], - "error": None, - "date": pieces.get("date"), - } - - -def get_versions() -> dict[str, Any]: - """Get version information or return default if unable to do so.""" - # I am in _version.py, which lives at ROOT/VERSIONFILE_SOURCE. If we have - # __file__, we can work backwards from there to the root. Some - # py2exe/bbfreeze/non-CPython implementations don't do __file__, in which - # case we can only use expanded keywords. - - cfg = get_config() - verbose = cfg.verbose - - try: - return git_versions_from_keywords( - get_keywords(), - cfg.tag_prefix, - verbose, - ) - except NotThisMethod: - pass - - try: - root = os.path.realpath(__file__) - # versionfile_source is the relative path from the top of the source - # tree (where the .git directory might live) to this file. Invert - # this to find the root from __file__. - for _ in cfg.versionfile_source.split("/"): - root = os.path.dirname(root) - except NameError: - return { - "version": "0+unknown", - "full-revisionid": None, - "dirty": None, - "error": "unable to find root of source tree", - "date": None, - } - - try: - pieces = git_pieces_from_vcs(cfg.tag_prefix, root, verbose) - return render(pieces, cfg.style) - except NotThisMethod: - pass - - try: - if cfg.parentdir_prefix: - return versions_from_parentdir(cfg.parentdir_prefix, root, verbose) - except NotThisMethod: - pass - - return { - "version": "0+unknown", - "full-revisionid": None, - "dirty": None, - "error": "unable to compute version", - "date": None, - } diff --git a/openff/interchange/common/_nonbonded.py b/openff/interchange/common/_nonbonded.py index 9643107b4..551c157c8 100644 --- a/openff/interchange/common/_nonbonded.py +++ b/openff/interchange/common/_nonbonded.py @@ -101,7 +101,7 @@ class ElectrostaticsCollection(_NonbondedCollection): nonperiodic_potential: Literal["Coulomb", "cutoff", "no-cutoff"] = Field("Coulomb") exception_potential: Literal["Coulomb"] = Field("Coulomb") - _charges: dict[Any, _ElementaryChargeQuantity] = PrivateAttr(default_factory=dict) + _charges: dict[Any, _ElementaryChargeQuantity] = PrivateAttr() _charges_cached: bool = PrivateAttr(default=False) @computed_field # type: ignore[prop-decorator] diff --git a/openff/interchange/components/interchange.py b/openff/interchange/components/interchange.py index cc092df20..b71f9f7f9 100644 --- a/openff/interchange/components/interchange.py +++ b/openff/interchange/components/interchange.py @@ -5,7 +5,7 @@ from pathlib import Path from typing import TYPE_CHECKING, Literal, Union, overload -from openff.toolkit import ForceField, Molecule, Quantity, Topology, unit +from openff.toolkit import Molecule, Quantity, Topology, unit from openff.utilities.utilities import has_package, requires_package from pydantic import Field @@ -42,14 +42,15 @@ from openff.interchange.smirnoff._gbsa import SMIRNOFFGBSACollection from openff.interchange.warnings import InterchangeDeprecationWarning -if has_package("foyer"): - from foyer.forcefield import Forcefield as FoyerForcefield -if has_package("nglview"): - import nglview - if TYPE_CHECKING: import openmm import openmm.app + from openff.toolkit import ForceField + + if has_package("foyer"): + from foyer import Forcefield as FoyerForcefield + if has_package("nglview"): + import nglview class Interchange(_BaseModel): @@ -84,7 +85,7 @@ class Interchange(_BaseModel): @classmethod def from_smirnoff( cls, - force_field: ForceField, + force_field: "ForceField", topology: Topology | list[Molecule], box=None, positions=None, @@ -109,9 +110,12 @@ def from_smirnoff( The positions associated with atoms in the input topology. If ``None``, positions are taken from the molecules in topology, if present on all molecules. charge_from_molecules : `List[openff.toolkit.molecule.Molecule]`, optional - If specified, partial charges will be taken from the given molecules - instead of being determined by the force field. In either case, charges - on the input topology are ignored. + If specified, partial charges for any molecules isomorphic to those + given will be taken from the given molecules' `partial_charges` + attribute instead of being determined by the force field. All + molecules in this list must have partial charges assigned and must + not be isomorphic with any other molecules in the list. For all values + of this argument, charges on the input topology are ignored. partial_bond_orders_from_molecules : List[openff.toolkit.molecule.Molecule], optional If specified, partial bond orders will be taken from the given molecules instead of being determined by the force field. @@ -120,12 +124,14 @@ def from_smirnoff( Notes ----- - If the `Molecule` objects in the `topology` argument each contain conformers, the returned `Interchange` object - will have its positions set via concatenating the 0th conformer of each `Molecule`. + If the ``Molecule`` objects in the ``topology`` argument each contain + conformers, the returned ``Interchange`` object will have its positions + set via concatenating the 0th conformer of each ``Molecule``. - If the `Molecule` objects in the `topology` argument have stored partial charges, these are ignored and charges - are assigned according to the contents of the force field. To override the force field and use preset charges, - use the `charge_from_molecules` argument. + If the ``Molecule`` objects in the ``topology`` argument have stored + partial charges, these are ignored and charges are assigned according to + the contents of the force field. To override the force field and use + preset charges, use the ``charge_from_molecules`` argument. Examples -------- @@ -150,7 +156,7 @@ def from_smirnoff( topology=topology, box=box, positions=positions, - charge_from_molecules=charge_from_molecules, + molecules_with_preset_charges=charge_from_molecules, partial_bond_orders_from_molecules=partial_bond_orders_from_molecules, allow_nonintegral_charges=allow_nonintegral_charges, ) @@ -342,7 +348,7 @@ def to_gromacs( Molecule names in written files are not guaranteed to match the `Moleclue.name` attribute of the molecules in the topology, especially if they are empty strings or not unique. - See `to_gro` and `to_top` for more details. + See :py:meth:`to_gro ` and :py:meth:`to_top ` for more details. """ from openff.interchange.interop.gromacs.export._export import GROMACSWriter @@ -430,6 +436,9 @@ def to_gro(self, file_path: Path | str, decimal: int = 3): decimal: int, default=3 The number of decimal places to use when writing the GROMACS coordinate file. + Notes + ----- + Residue IDs must be positive integers (or string representations thereof). Residue IDs greater than 99,999 are reduced modulo 100,000 in line with common GROMACS practice. @@ -569,9 +578,10 @@ def to_openmm_topology( Parameters ---------- collate - If False, the default, all virtual sites will be added to a single residue at the end of the topology. - If True, virtual sites will be collated with their associated molecule and added to the residue of the last - atom in the molecule they belong to. + If ``False``, the default, all virtual sites will be added to a + single residue at the end of the topology. If ``True``, virtual + sites will be collated with their associated molecule and added to + the residue of the last atom in the molecule they belong to. """ from openff.interchange.interop.openmm._topology import to_openmm_topology @@ -593,9 +603,14 @@ def to_openmm_simulation( **kwargs, ) -> "openmm.app.simulation.Simulation": """ - Export this Interchange to an OpenMM `Simulation` object. + Export this Interchange to an OpenMM ``Simulation`` object. + + A :py:class:`Simulation ` encapsulates + all the information needed for a typical OpenMM simulation into a single + object with a simple API. - Positions are set on the `Simulation` if present on the `Interchange`. + Positions are set on the ``Simulation`` if present on the + ``Interchange``. Additional forces, such as a barostat, should be added with the ``additional_forces`` argument to avoid having to re-initialize diff --git a/openff/interchange/components/potentials.py b/openff/interchange/components/potentials.py index 5161f2cb8..123ec4161 100644 --- a/openff/interchange/components/potentials.py +++ b/openff/interchange/components/potentials.py @@ -2,9 +2,10 @@ import ast import json -from typing import Annotated, Any +from typing import TYPE_CHECKING, Annotated, Any, TypeAlias import numpy +from numpy.typing import ArrayLike from openff.toolkit import Quantity from openff.utilities.utilities import has_package, requires_package from pydantic import ( @@ -25,15 +26,13 @@ ) from openff.interchange.pydantic import _BaseModel -if has_package("jax"): - from jax import numpy as jax_numpy - -from numpy.typing import ArrayLike - -if has_package("jax"): - from jax import Array +if TYPE_CHECKING: + if has_package("jax"): + from jax import Array + else: + Array: TypeAlias = Any # type: ignore[no-redef] else: - Array = Any # type: ignore + Array: TypeAlias = ArrayLike class Potential(_BaseModel): @@ -269,6 +268,8 @@ def get_force_field_parameters( raise NotImplementedError if use_jax: + from jax import numpy as jax_numpy + return jax_numpy.array( [[v.m for v in p.parameters.values()] for p in self.potentials.values()], ) @@ -318,6 +319,8 @@ def get_system_parameters( q.append(p[index]) if use_jax: + from jax import numpy as jax_numpy + return jax_numpy.array(q) else: return numpy.array(q) diff --git a/openff/interchange/exceptions.py b/openff/interchange/exceptions.py index 57413a220..1124bbb9e 100644 --- a/openff/interchange/exceptions.py +++ b/openff/interchange/exceptions.py @@ -318,6 +318,12 @@ class MissingPartialChargesError(InterchangeException): """ +class PresetChargesError(InterchangeException): + """ + The `charge_from_molecules` is structured in an undefined way. + """ + + class UnassignedValenceError(InterchangeException): """Exception raised when there are valence terms for which a ParameterHandler did not find matches.""" diff --git a/openff/interchange/foyer/_base.py b/openff/interchange/foyer/_base.py index 52fe4841c..c4b8fdd94 100644 --- a/openff/interchange/foyer/_base.py +++ b/openff/interchange/foyer/_base.py @@ -1,15 +1,14 @@ from abc import abstractmethod from copy import copy +from typing import TYPE_CHECKING from openff.toolkit import Topology -from openff.utilities.utilities import has_package from openff.interchange.components.potentials import Collection, Potential from openff.interchange.models import PotentialKey, TopologyKey -if has_package("foyer"): - from foyer.exceptions import MissingForceError, MissingParametersError - from foyer.forcefield import Forcefield +if TYPE_CHECKING: + from foyer import Forcefield # Is this the safest way to achieve PotentialKey id separation? @@ -66,6 +65,8 @@ def store_matches( def store_potentials(self, force_field: "Forcefield") -> None: """Populate self.potentials with key-val pairs of [PotentialKey, Potential].""" + from foyer.exceptions import MissingForceError, MissingParametersError + for pot_key in self.key_map.values(): try: params = force_field.get_parameters( diff --git a/openff/interchange/foyer/_nonbonded.py b/openff/interchange/foyer/_nonbonded.py index 75a7e2fbd..947eb58ff 100644 --- a/openff/interchange/foyer/_nonbonded.py +++ b/openff/interchange/foyer/_nonbonded.py @@ -11,7 +11,7 @@ from openff.interchange.models import PotentialKey, TopologyKey if has_package("foyer"): - from foyer.forcefield import Forcefield + from foyer import Forcefield class FoyerVDWHandler(vdWCollection): @@ -60,7 +60,7 @@ class FoyerElectrostaticsHandler(ElectrostaticsCollection): force_field_key: str = "atoms" cutoff: _DistanceQuantity = 9.0 * unit.angstrom - _charges: dict[TopologyKey, Quantity] = PrivateAttr(default_factory=dict) + _charges: dict[TopologyKey, Quantity] = PrivateAttr(dict()) def store_charges( self, diff --git a/openff/interchange/interop/openmm/_import/_import.py b/openff/interchange/interop/openmm/_import/_import.py index 5485616d9..5c8cf28a7 100644 --- a/openff/interchange/interop/openmm/_import/_import.py +++ b/openff/interchange/interop/openmm/_import/_import.py @@ -7,7 +7,6 @@ from openff.utilities.utilities import has_package, requires_package from pydantic import ValidationError -from openff.interchange._experimental import experimental from openff.interchange.common._nonbonded import vdWCollection from openff.interchange.common._valence import ( AngleCollection, @@ -38,7 +37,6 @@ @requires_package("openmm") -@experimental def from_openmm( *, system: "openmm.System", diff --git a/openff/interchange/models.py b/openff/interchange/models.py index 95a159c62..7b84c062b 100644 --- a/openff/interchange/models.py +++ b/openff/interchange/models.py @@ -298,6 +298,8 @@ class SingleAtomChargeTopologyKey(LibraryChargeTopologyKey): Shim class for storing the result of charge_from_molecules. """ + extras: dict = dict() # noqa: RUF012 + class ChargeModelTopologyKey(_BaseModel): """Subclass of `TopologyKey` for use with charge models only.""" diff --git a/openff/interchange/smirnoff/_create.py b/openff/interchange/smirnoff/_create.py index 56dfab7d1..5a6da6d65 100644 --- a/openff/interchange/smirnoff/_create.py +++ b/openff/interchange/smirnoff/_create.py @@ -1,3 +1,5 @@ +import warnings + from openff.toolkit import ForceField, Molecule, Quantity, Topology from openff.toolkit.typing.engines.smirnoff import ParameterHandler from openff.toolkit.typing.engines.smirnoff.plugins import load_handler_plugins @@ -8,6 +10,7 @@ from openff.interchange.components.toolkit import _check_electrostatics_handlers from openff.interchange.exceptions import ( MissingParameterHandlerError, + PresetChargesError, SMIRNOFFHandlersNotImplementedError, ) from openff.interchange.plugins import load_smirnoff_plugins @@ -25,6 +28,7 @@ SMIRNOFFProperTorsionCollection, ) from openff.interchange.smirnoff._virtual_sites import SMIRNOFFVirtualSiteCollection +from openff.interchange.warnings import PresetChargesAndVirtualSitesWarning _SUPPORTED_PARAMETER_HANDLERS: set[str] = { "Constraints", @@ -93,17 +97,65 @@ def validate_topology(value): ) +def _preprocess_preset_charges( + molecules_with_preset_charges: list[Molecule] | None, +) -> list[Molecule] | None: + """ + Pre-process the molecules_with_preset_charges argument. + + If molecules_with_preset_charges is None, return None. + + If molecules_with_preset_charges is list[Molecule], ensure that + + 1. The input is a list of Molecules + 2. Each molecule has assign partial charges + 3. Ensure no molecules are isomorphic with another in the list + + """ + if molecules_with_preset_charges is None: + return None + + # This relies on Molecule.__eq__(), which may change and/or not have the same equality criteria + # as we want here. + # See https://github.com/openforcefield/openff-interchange/pull/1070#discussion_r1792728179 + molecule_set = {molecule for molecule in molecules_with_preset_charges} + + if len(molecule_set) != len(molecules_with_preset_charges): + raise PresetChargesError( + "All molecules in the molecules_with_preset_charges list must be isomorphically unique from each other", + ) + + for molecule in molecules_with_preset_charges: + if molecule.partial_charges is None: + raise PresetChargesError( + "All molecules in the molecules_with_preset_charges list must have partial charges assigned.", + ) + + return molecules_with_preset_charges + + def _create_interchange( force_field: ForceField, topology: Topology | list[Molecule], box: Quantity | None = None, positions: Quantity | None = None, - charge_from_molecules: list[Molecule] | None = None, + molecules_with_preset_charges: list[Molecule] | None = None, partial_bond_orders_from_molecules: list[Molecule] | None = None, allow_nonintegral_charges: bool = False, ) -> Interchange: + molecules_with_preset_charges = _preprocess_preset_charges(molecules_with_preset_charges) + _check_supported_handlers(force_field) + if molecules_with_preset_charges is not None and "VirtualSites" in force_field.registered_parameter_handlers: + warnings.warn( + "Preset charges were provided (via `charge_from_molecules`) alongside a force field that includes " + "virtual site parameters. Note that virtual sites will be applied charges from the force field and " + "cannot be given preset charges. Virtual sites may also affect the charges of their orientation " + "atoms, even if those atoms are given preset charges.", + PresetChargesAndVirtualSitesWarning, + ) + # interchange = Interchange(topology=topology) # or maybe interchange = Interchange(topology=validate_topology(topology)) @@ -138,7 +190,7 @@ def _create_interchange( interchange, force_field, interchange.topology, - charge_from_molecules, + molecules_with_preset_charges, allow_nonintegral_charges, ) _plugins(interchange, force_field, interchange.topology) @@ -274,7 +326,7 @@ def _electrostatics( interchange: Interchange, force_field: ForceField, topology: Topology, - charge_from_molecules: list[Molecule] | None = None, + molecules_with_preset_charges: list[Molecule] | None = None, allow_nonintegral_charges: bool = False, ): if "Electrostatics" not in force_field.registered_parameter_handlers: @@ -304,7 +356,7 @@ def _electrostatics( if handler is not None ], topology=topology, - charge_from_molecules=charge_from_molecules, + molecules_with_preset_charges=molecules_with_preset_charges, allow_nonintegral_charges=allow_nonintegral_charges, ), }, diff --git a/openff/interchange/smirnoff/_nonbonded.py b/openff/interchange/smirnoff/_nonbonded.py index 4a0286a7b..de4b800aa 100644 --- a/openff/interchange/smirnoff/_nonbonded.py +++ b/openff/interchange/smirnoff/_nonbonded.py @@ -1,5 +1,6 @@ import copy import functools +import logging import warnings from collections.abc import Iterable from typing import Any, Literal, Optional, Union @@ -42,6 +43,8 @@ from openff.interchange.smirnoff._base import SMIRNOFFCollection from openff.interchange.warnings import ForceFieldModificationWarning +logger = logging.getLogger(__name__) + ElectrostaticsHandlerType = Union[ ElectrostaticsHandler, ToolkitAM1BCCHandler, @@ -273,7 +276,7 @@ class SMIRNOFFElectrostaticsCollection(ElectrostaticsCollection, SMIRNOFFCollect ) # type: ignore[assignment] exception_potential: Literal["Coulomb"] = Field("Coulomb") - _charges: dict[Any, _ElementaryChargeQuantity] = PrivateAttr(default_factory=dict) + _charges: dict[Any, _ElementaryChargeQuantity] = PrivateAttr(dict()) _charges_cached: bool = PrivateAttr(default=False) @classmethod @@ -304,6 +307,19 @@ def charges( ) -> dict[TopologyKey | LibraryChargeTopologyKey | VirtualSiteKey, _ElementaryChargeQuantity]: """Get the total partial charge on each atom, including virtual sites.""" if len(self._charges) == 0 or self._charges_cached is False: + # TODO: Clearing this dict **should not be necessary** but in some cases in appears + # that this class attribute persists in some cases, possibly only on a single + # thread. Ideas to try include + # * Drop @computed_field ? + # * Don't handle caching ourselves and let Pydantic do it (i.e. have this + # function simply retrun ._get_charges(...) + # + # Hopefully this isn't a major issue - caches for large systems should still + # only be build once + # + # Some more context: https://github.com/openforcefield/openff-interchange/issues/842#issuecomment-2394211357 + self._charges.clear() + self._charges.update(self._get_charges(include_virtual_sites=True)) self._charges_cached = True @@ -330,6 +346,7 @@ def _get_charges( total_charge: Quantity = numpy.sum(parameter_value) # assumes virtual sites can only have charges determined in one step + charges[topology_key] = -1.0 * total_charge # Apply increments to "orientation" atoms @@ -349,7 +366,7 @@ def _get_charges( if potential_key.associated_handler in ( "LibraryCharges", "ToolkitAM1BCCHandler", - "charge_from_molecules", + "molecules_with_preset_charges", "ExternalSource", ): charges[atom_index] = _add_charges( @@ -384,6 +401,11 @@ def _get_charges( parameter_value, ) + logger.info( + "Charge section ChargeIncrementModel, applying charge increment from atom " + f"{topology_key.this_atom_index} to atoms {topology_key.other_atom_indices}", + ) + else: raise NotImplementedError() @@ -416,7 +438,7 @@ def create( cls, parameter_handler: Any, topology: Topology, - charge_from_molecules=None, + molecules_with_preset_charges=None, allow_nonintegral_charges: bool = False, ) -> Self: """ @@ -455,7 +477,7 @@ def create( handler.store_matches( parameter_handlers, topology, - charge_from_molecules=charge_from_molecules, + molecules_with_preset_charges=molecules_with_preset_charges, allow_nonintegral_charges=allow_nonintegral_charges, ) handler._charges = dict() @@ -678,7 +700,15 @@ def _find_charge_model_matches( ) potentials[potential_key] = Potential(parameters={"charge": partial_charge}) - matches[SingleAtomChargeTopologyKey(this_atom_index=atom_index)] = potential_key + matches[ + SingleAtomChargeTopologyKey( + this_atom_index=atom_index, + extras={ + "handler": handler_name, + "partial_charge_method": partial_charge_method, + }, + ) + ] = potential_key return partial_charge_method, matches, potentials @@ -777,12 +807,12 @@ def _assign_charges_from_molecules( cls, topology: Topology, unique_molecule: Molecule, - charge_from_molecules=Optional[list[Molecule]], + molecules_with_preset_charges=Optional[list[Molecule]], ) -> tuple[bool, dict, dict]: - if charge_from_molecules is None: + if molecules_with_preset_charges is None: return False, dict(), dict() - for molecule_with_charges in charge_from_molecules: + for molecule_with_charges in molecules_with_preset_charges: if molecule_with_charges.is_isomorphic_with(unique_molecule): break else: @@ -806,11 +836,12 @@ def _assign_charges_from_molecules( index_in_topology = atom_map[index_in_molecule_with_charges] topology_key = SingleAtomChargeTopologyKey( this_atom_index=index_in_topology, + extras={"handler": "preset"}, ) potential_key = PotentialKey( id=mapped_smiles, mult=index_in_molecule_with_charges, # Not sure this prevents clashes in some corner cases - associated_handler="charge_from_molecules", + associated_handler="molecules_with_preset_charges", bond_order=None, ) potential = Potential(parameters={"charge": partial_charge}) @@ -823,7 +854,7 @@ def store_matches( self, parameter_handler: ElectrostaticsHandlerType | list[ElectrostaticsHandlerType], topology: Topology, - charge_from_molecules=None, + molecules_with_preset_charges=None, allow_nonintegral_charges: bool = False, ) -> None: """ @@ -846,10 +877,10 @@ def store_matches( flag, matches, potentials = self._assign_charges_from_molecules( topology, unique_molecule, - charge_from_molecules, + molecules_with_preset_charges, ) # TODO: Here is where the toolkit calls self.check_charges_assigned(). Do we skip this - # entirely given that we are not accepting `charge_from_molecules`? + # entirely given that we are not accepting `molecules_with_preset_charges`? if not flag: # TODO: Rename this method to something like `_find_matches` @@ -881,6 +912,42 @@ def store_matches( # Have this new key (on a duplicate molecule) point to the same potential # as the old key (on a unique/reference molecule) + if type(new_key) is LibraryChargeTopologyKey: + logger.info( + "Charge section LibraryCharges applied to topology atom index " + f"{topology_atom_index}", + ) + + elif type(new_key) is SingleAtomChargeTopologyKey: + if new_key.extras["handler"] == "ToolkitAM1BCCHandler": + logger.info( + "Charge section ToolkitAM1BCC, using charge method " + f"{new_key.extras['partial_charge_method']}, " + f"applied to topology atom index {topology_atom_index}", + ) + + elif new_key.extras["handler"] == "preset": + logger.info( + f"Preset charges applied to atom index {topology_atom_index}", + ) + + else: + raise ValueError(f"Unhandled handler {new_key.extras['handler']}") + + elif type(new_key) is ChargeModelTopologyKey: + logger.info( + "Charge section ChargeIncrementModel, using charge method " + f"{new_key.partial_charge_method}, " + f"applied to topology atom index {new_key.this_atom_index}", + ) + + elif type(new_key) is ChargeIncrementTopologyKey: + # here is where the actual increments could be logged + pass + + else: + raise ValueError(f"Unhandled key type {type(new_key)}") + self.key_map[new_key] = matches[key] topology_charges = [0.0] * topology.n_atoms diff --git a/openff/interchange/smirnoff/_virtual_sites.py b/openff/interchange/smirnoff/_virtual_sites.py index 0629a88af..7f1955a98 100644 --- a/openff/interchange/smirnoff/_virtual_sites.py +++ b/openff/interchange/smirnoff/_virtual_sites.py @@ -1,3 +1,4 @@ +import logging import math from typing import Literal @@ -24,6 +25,8 @@ SMIRNOFFvdWCollection, ) +logger = logging.getLogger(__name__) + _DEGREES_TO_RADIANS = numpy.pi / 180.0 @@ -179,6 +182,12 @@ def store_potentials( # type: ignore[override] ), }, ) + + logger.info( + "Charge section VirtualSites applied to virtual site with orientation atoms at topology indices " + f"{virtual_site_key.orientation_atom_indices}", + ) + electrostatics_collection.key_map[virtual_site_key] = electrostatics_key electrostatics_collection.potentials[electrostatics_key] = electrostatics_potential diff --git a/openff/interchange/warnings.py b/openff/interchange/warnings.py index 67ca825b2..3a3dc4aab 100644 --- a/openff/interchange/warnings.py +++ b/openff/interchange/warnings.py @@ -21,5 +21,9 @@ class ForceFieldModificationWarning(UserWarning): """Warning for when a ForceField is modified.""" +class PresetChargesAndVirtualSitesWarning(UserWarning): + """Warning when possibly using preset charges and virtual sites together.""" + + class InterchangeCombinationWarning(UserWarning): """Warning for when combining Interchange objects.""" diff --git a/pyproject.toml b/pyproject.toml index fb79f2dca..6163e361a 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,13 +1,19 @@ [build-system] -requires = ["setuptools", "versioneer[toml]==0.29"] +requires = [ + "setuptools", + "versioningit", +] build-backend = "setuptools.build_meta" -[tool.versioneer] -VCS = "git" -style = "pep440" -versionfile_source = "openff/interchange/_version.py" -versionfile_build = "openff/interchange/_version.py" -tag_prefix = "v" +[project] +name="openff-interchange" +description = "A project (and object) for storing, manipulating, and converting molecular mechanics data." +readme = "README.md" +authors = [{name = "Open Force Field Initiative", email = "info@openforcefield.org"}] +license = {text = "MIT"} +dynamic = ["version"] + +[tool.versioningit] [tool.interrogate] ignore-init-method = true diff --git a/setup.cfg b/setup.cfg index 4aa3f4299..ed47a2469 100644 --- a/setup.cfg +++ b/setup.cfg @@ -1,6 +1,5 @@ [coverage:run] omit = - openff/interchange/_version.py */*/_tests/* [coverage:report] diff --git a/setup.py b/setup.py index 175e40feb..1b36eb65a 100644 --- a/setup.py +++ b/setup.py @@ -1,23 +1,6 @@ -import versioneer from setuptools import find_namespace_packages, setup -short_description = "A project (and object) for storing, manipulating, and converting molecular mechanics data." - -with open("README.md") as handle: - long_description = handle.read() - - setup( - name="openff-interchange", - author="Open Force Field Initiative", - author_email="info@openforcefield.org", - description=short_description[0], - long_description=long_description, - long_description_content_type="text/markdown", - version=versioneer.get_version(), - cmdclass=versioneer.get_cmdclass(), - license="MIT", packages=find_namespace_packages(), include_package_data=True, - setup_requires=[], ) diff --git a/stubs/foyer/__init__.py b/stubs/foyer/__init__.py new file mode 100644 index 000000000..826925080 --- /dev/null +++ b/stubs/foyer/__init__.py @@ -0,0 +1,3 @@ +from foyer.forcefield import Forcefield + +__all__ = ("Forcefield",)