Multi-Molecule Constrained Parameterization¶

Transferability across multiple molecules that share common atom types is convenient for classical force field partial charge parameterization as well as machine learning interatomic potential (MLIP) feature generation. To this end, this tutorial highlights an automated atom typing scheme using Morgan fingerprints to constrain shared chemical environments across molecules to carry identical charges. The full runnable script is at examples/tutorials/multi_constrained_fitting.py.

1a. Molecule, Conformer, and Record Generation¶

First, as before, we select four ionic liquid cations for this example: MMIM, EMIM, BMIM, and C6MIM. As in the previous tutorials, we define the SMILES, generate conformers, and create the GDMA records. Here, it is useful to save the .sqlite record file for later use.

import pathlib

from openff.recharge.utilities.molecule import extract_conformers
from openff.toolkit import Molecule

from pympfit import GDMASettings, MoleculeGDMARecord, Psi4GDMAGenerator
from pympfit.gdma.storage import MoleculeGDMAStore

smiles_map = {
    "MMIM": "CN1C=C[N+](=C1)C",
    "EMIM": "CCN1C=C[N+](=C1)C",
    "BMIM": "CCCCN1C=C[N+](=C1)C",
    "C6MIM": "CCCCCCN1C=C[N+](=C1)C",
}

molecules = [Molecule.from_smiles(smi) for smi in smiles_map.values()]
names = list(smiles_map.keys())

settings = GDMASettings()
records = []
for mol in molecules:
    mol.generate_conformers(n_conformers=1)
    [conf] = extract_conformers(mol)
    coords, multipoles = Psi4GDMAGenerator.generate(
        mol, conf, settings, minimize=True,
    )
    records.append(
        MoleculeGDMARecord.from_molecule(mol, coords, multipoles, settings)
    )

1b. [Optional] Import GDMA SQLite Records¶

Instead of calculating expensive Psi4/GDMA multipoles each time, one can also import the .sqlite data from a previous run or as provided in this repo at pympfit/pympfit/tests/data/gdma/ionic_liquids.sqlite:

SQLITE_FILE = pathlib.Path("pympfit/tests/data/gdma/ionic_liquids.sqlite")

store = MoleculeGDMAStore(str(SQLITE_FILE))
records = []
for smi in smiles_map.values():
    recs = store.retrieve(smiles=smi)
    records.append(recs[0])

2. Generate Atom Type Labels¶

Morgan fingerprints are molecular fingerprints — each fingerprint bit (provided as a unique hash) corresponds to a fragment of the molecule based on a radius parameter, measured in number of bonds beyond the selected atom.

generate_global_atom_type_labels combines these fingerprints with within-molecule symmetry detection to assign labels. Atoms with the same label across different molecules will be constrained to carry the same charge.

from pympfit import generate_global_atom_type_labels

labels = generate_global_atom_type_labels(
    molecules,
    radius=2,
    equivalize_between_methyl_carbons=True,
    equivalize_between_methyl_hydrogens=True,
    equivalize_between_other_heavy_atoms=True,
    equivalize_between_other_hydrogen_atoms=True,
)

The shared atom types can be visualized using the tutorial utility:

from util import _generate_rdkit_imgs

_generate_rdkit_imgs(molecules, labels, names, "figures/", prefix="multi_molecule_")

The generated atom labels highlight the transferable atom types across the molecules. Atoms with the same color share the same label and will be constrained to carry identical charges. Unhighlighted atoms are molecule-specific.

3. Obtain Charges¶

Pass all records and molecules together. The solver fits a single reduced parameter vector across all four molecules simultaneously, enforcing that shared atom types carry identical total charges. Each molecule’s total charge is constrained to +1.0 (cation) via the conchg penalty.

from pympfit import ConstrainedSciPySolver, generate_constrained_mpfit_charge_parameter

solver = ConstrainedSciPySolver(conchg=10.0)
parameters = generate_constrained_mpfit_charge_parameter(
    records, molecules, solver=solver,
)

Each molecule’s total charge is +1.0, enforced by the conservation penalty. Shared atom types carry identical charges across molecules — for example, all N2 atoms have charge +0.0750, all H5 atoms have +0.0605, and all C3 atoms have +0.0687, regardless of which molecule they appear in.