Generating Canonical AM1-BCC Charges

Problem

You want to generate general, canonical AM1-BCC charges from a single input structure. The atomic partial charges need to be appropriate for a wide spectrum of energy-based methods including minimization, posing in an active site, SZMAP, EON, or molecular dynamics.

Ingredients

Difficulty Level

🌶️ 🌶️

Download

Download code

can_am1_bcc.py

See also the Usage subsection.

Source Code

can_am1_bcc
#!/usr/bin/env python3
# (C) 2026 Cadence Design Systems, Inc. (Cadence)
# All rights reserved.
# TERMS FOR USE OF SAMPLE CODE The software below ("Sample Code") is
# provided to current licensees or subscribers of Cadence products or
# SaaS offerings (each a "Customer").
# Customer is hereby permitted to use, copy, and modify the Sample Code,
# subject to these terms. Cadence claims no rights to Customer's
# modifications. Modification of Sample Code is at Customer's sole and
# exclusive risk. Sample Code may require Customer to have a then
# current license or subscription to the applicable Cadence offering.
# THE SAMPLE CODE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
# EXPRESS OR IMPLIED.  OPENEYE DISCLAIMS ALL WARRANTIES, INCLUDING, BUT
# NOT LIMITED TO, WARRANTIES OF MERCHANTABILITY, FITNESS FOR A
# PARTICULAR PURPOSE AND NONINFRINGEMENT. In no event shall Cadence be
# liable for any damages or liability in connection with the Sample Code
# or its use.

"""Generate canonical AM1-BCC charges for molecules."""

import argparse
import os
import pathlib
import sys

from openeye import oechem, oeomega, oequacpac
from rich.console import Console
from rich_argparse import HelpPreviewAction, RichHelpFormatter

__SCRIPT_NAME__ = pathlib.Path(__file__).absolute().stem
__SCRIPT_DESC__ = "Generate canonical AM1-BCC charges for molecules."
__SCRIPT_TOOLKITS__ = ["oechem", "oeomega", "oequacpac"]


def parse_options() -> argparse.Namespace:
    """Set up command line options."""
    parser = argparse.ArgumentParser(
        add_help=True,
        formatter_class=RichHelpFormatter,
        description="[yellow]" + __SCRIPT_DESC__ + "[/yellow]",
    )

    io_group = parser.add_argument_group("Input/output options")
    io_group.add_argument(
        "--in",
        dest="input",
        metavar="INPUT-FILE",
        type=str,
        required=True,
        help="input molecule file with 3D coordinates",
    )
    io_group.add_argument(
        "--out",
        dest="output",
        metavar="OUTPUT-FILE",
        type=str,
        required=True,
        help="output molecule file (mol2 or oeb)",
    )

    omega_group = parser.add_argument_group("Omega options")
    omega_group.add_argument(
        "--energy-window",
        type=float,
        default=15.0,
        help="energy window for conformer generation (default: %(default)s)",
    )
    omega_group.add_argument(
        "--max-confs",
        type=int,
        default=800,
        help="maximum number of conformers (default: %(default)s)",
    )
    omega_group.add_argument(
        "--rms-threshold",
        type=float,
        default=1.0,
        help="RMS threshold for conformer deduplication (default: %(default)s)",
    )

    parser.add_argument("--help-image", action=HelpPreviewAction)
    return parser.parse_args()


def _setup_omega(args: argparse.Namespace) -> oeomega.OEOmega:
    """Configure and return an Omega conformer generator."""
    opts = oeomega.OEOmegaOptions()
    opts.SetIncludeInput(True)
    opts.SetCanonOrder(False)
    opts.SetSampleHydrogens(True)
    opts.SetEnergyWindow(args.energy_window)
    opts.SetMaxConfs(args.max_confs)
    opts.SetRMSThreshold(args.rms_threshold)
    return oeomega.OEOmega(opts)


def assign_canonical_charges(
    mol: oechem.OEMCMolBase, omega: oeomega.OEOmega, console: Console
) -> bool:
    """Generate conformers and assign AM1-BCC ELF10 charges to a molecule."""
    if not omega(mol):
        console.print(f"[red]Failed to generate conformers for {mol.GetTitle()}[/red]")
        return False

    oequacpac.OEAssignCharges(mol, oequacpac.OEAM1BCCELF10Charges())

    sum_formal = 0
    abs_formal = 0
    sum_partial = 0.0
    for atom in mol.GetAtoms():
        sum_formal += atom.GetFormalCharge()
        abs_formal += abs(atom.GetFormalCharge())
        sum_partial += atom.GetPartialCharge()

    console.print(
        f"[green]{mol.GetTitle()}[/green]: {abs_formal} formal charges give "
        f"total charge {sum_formal}; sum of partial charges {sum_partial:5.4f}"
    )
    return True


def main() -> int:
    """Generate canonical AM1-BCC charges."""
    args = parse_options()

    console = Console()

    ifs = oechem.oemolistream()
    if not ifs.open(args.input):
        console.print(f"[red]Unable to open {args.input} for reading.[/red]")
        return os.EX_USAGE

    if not oechem.OEIs3DFormat(ifs.GetFormat()):
        console.print("[red]Invalid input format: need 3D coordinates.[/red]")
        return os.EX_USAGE

    ofs = oechem.oemolostream()
    if not ofs.open(args.output):
        console.print(f"[red]Unable to open {args.output} for writing.[/red]")
        return os.EX_USAGE

    if ofs.GetFormat() not in [oechem.OEFormat_MOL2, oechem.OEFormat_OEB]:
        console.print("[red]MOL2 or OEB output file is required.[/red]")
        return os.EX_USAGE

    omega = _setup_omega(args)

    for mol in ifs.GetOEMols():
        if assign_canonical_charges(mol, omega, console):
            conf = mol.GetConf(oechem.OEHasConfIdx(0))
            oechem.OEWriteMolecule(ofs, conf)

    console.print(f"Output written to [green]{args.output}[/green]")
    return os.EX_OK


setattr(main, "__SCRIPT_NAME__", __SCRIPT_NAME__)
setattr(main, "__SCRIPT_DESC__", __SCRIPT_DESC__)
setattr(main, "__SCRIPT_TOOLKITS__", __SCRIPT_TOOLKITS__)

if __name__ == "__main__":
    sys.exit(main())

Solution

This script avoids the two main problems that can occur when computing AM1-BCC charges from a single structure: asymmetric charges and distorted charges caused by short-range polar interactions. Conformers are generated with Omega and the ELF10 method automatically selects the best conformer ensemble used to compute well-behaved, symmetric AM1-BCC charges. The input structure is then written out with these atomic partial charges, ready for use.

The assign_canonical_charges function generates conformers for a molecule using Omega and then assigns AM1-BCC ELF10 charges via the OEAssignCharges function with the OEAM1BCCELF10Charges method.

def assign_canonical_charges(
    mol: oechem.OEMCMolBase, omega: oeomega.OEOmega, console: Console
) -> bool:
    """Generate conformers and assign AM1-BCC ELF10 charges to a molecule."""
    if not omega(mol):
        console.print(f"[red]Failed to generate conformers for {mol.GetTitle()}[/red]")
        return False

    oequacpac.OEAssignCharges(mol, oequacpac.OEAM1BCCELF10Charges())

    sum_formal = 0
    abs_formal = 0
    sum_partial = 0.0
    for atom in mol.GetAtoms():
        sum_formal += atom.GetFormalCharge()
        abs_formal += abs(atom.GetFormalCharge())
        sum_partial += atom.GetPartialCharge()

    console.print(
        f"[green]{mol.GetTitle()}[/green]: {abs_formal} formal charges give "
        f"total charge {sum_formal}; sum of partial charges {sum_partial:5.4f}"
    )
    return True

Usage

See the Download section to download the script.

> can_am1_bcc --help
../_images/can_am1_bcc-help.svg

The command below show how to run the script with an input file input.mol2.

> can_am1_bcc --in input.mol2 --out charged.oeb

See also in Omega TK manual

API

See also in Quacpac TK manual

API