Force Field Assignment

Force Field Assignment—the assignment of partial charges and force field terms (bonds, angles, dihedrals)—is traditionally the bottleneck in high-throughput molecular dynamics.

ChemFAST overcomes this by treating assignment as a hierarchical data problem. It employs a Hybrid Strategy that combines the rigorous accuracy of experimental databases with the infinite chemical coverage of Graph Neural Networks (GNNs).

1. The Hybrid Strategy

ChemFAST does not rely on a single method. Instead, for every molecule in your system, the engine provides an API allowing you to make a dynamic decision based on the chemical familiarity of the substructures.

The default workflow follows a strict priority logic:

  1. Tier 1: Database Exact Match (High Fidelity) The engine first queries the internal OPLS-AA database. If a chemical fragment (defined by its graph hash) is found, the validated parameters are retrieved directly. This ensures that well-known molecules (e.g., common solvents, amino acids, standard monomers) use parameters identical to established literature.

  2. Tier 2: Machine Learning Prediction (High Generalizability) If a fragment is “unknown” to the database (e.g., a novel conjugated polymer backbone or a transition state mimic), ChemFAST falls back to its pre-trained Graph Attention Network (GAT) to predict the parameters.

2. Tier 1: Database Matching via WL-Hash

ChemFAST structures its force field library as a relational SQL database. To rigorously identify chemical environments, the engine employs a multi-modal query system utilizing SMARTS patterns, atom-type names, and Weisfeiler-Lehman (WL) Graph Hashing. Unlike traditional ‘Atom Typing’ which relies on ambiguous heuristic labels (e.g., CA, CT, O_3), WL-Hashing mathematically encodes the precise topological neighborhood of each atom, ensuring high-fidelity parameter retrieval.

  • The Algorithm:

    For every atom, the algorithm iteratively aggregates information from its neighbors (up to 2~3 hops). This generates a unique hash string that represents the atom’s “chemical environment DNA”.

  • The Query:

    ChemFAST compares this hash against its internal SQL database (opls.db), which contains millions of pre-computed hashes from the standard OPLS-AA and BOSS/LigParGen datasets.

  • Key Advantages:

    1. Robustness: By relying on exact topological hashing rather than arbitrary nomenclature, DoMD significantly reduces “missing atom type” errors caused by naming mismatches in input files.

    2. High Automation: The process eliminates the laborious human cost of manual parameter lookup and “atom-typing,” thereby preventing user-induced transcription errors.

    3. Extensibility: The SQL architecture is designed to be open and modular. Users can freely expand the parameter space by injecting custom databases—such as those derived from ab initio calculations (e.g., using BOSS) or future iterations of the OPLS force field—without modifying the core code.

3. Tier 2: GAT-Based ML Prediction

When the database search returns no results for a specific chemical environment, ChemFAST activates the Graph Attention Network (GAT) module.

Model Architecture: Hierarchical Graph Representation

A molecule transcends a mere collection of atoms; it is a hierarchy of geometric entities ranging from atoms to bonds, angles, and torsions. To accurately predict parameters for these distinct multi-body terms, DoMD employs a Hierarchical Graph Attention Network (HGAT).

This architecture operates on three topologically coupled graph representations, successively transforming the molecular graph \(\mathcal{G}\) into higher-order line graphs to explicitly represent complex geometric dependencies.

Level 1: The Atom-Centric Graph (\(\mathcal{G}\))

At the foundational level, the model constructs the primary molecular graph where atoms are nodes and chemical bonds are edges.

  • Nodes: Atoms (\(u\)). Input features include elemental identity, hybridization state, aromaticity, and ring membership.

  • Edges: Chemical Bonds (\(e_{uv}\)).

  • Predictions:

    • Partial Charges (\(q_i\)): A physics-informed regression head predicts the partial charge for each atom-node.

    • LJ Types (\(\sigma, \epsilon\)): A multi-class classification head predicts the atom type index, mapping directly to Lennard-Jones parameters.

    • Improper Dihedrals: In OPLS, impropers maintain planarity and are defined by a central “hub” atom. The model predicts improper types via node classification on the central atom.

    • Bond Types (\(K_b, r_0\)): The model performs edge classification on this graph to assign parameters to each chemical bond.

Level 2: The Bond-Centric Graph (\(\mathcal{L}(\mathcal{G})\))

This level represents the Line Graph of \(\mathcal{G}\), effectively shifting the focus from atoms to interactions.

  • Transformation: Each bond \(e_{uv}\) in the atom graph becomes a Node in this representation.

  • Edges: Two bond-nodes are connected if they share a common atom (i.e., they form an angle geometry).

  • Predictions:

    • Angle Types (\(K_\theta, \theta_0\)): An edge classification head predicts parameters for each connection between bonds, which corresponds physically to a bond angle.

Level 3: The Angle-Centric Graph (\(\mathcal{L}(\mathcal{L}(\mathcal{G}))\))

This level represents the Line Graph of the Bond Graph (the 2nd-order Line Graph of \(\mathcal{G}\)), capturing long-range dependencies.

  • Transformation: Each angle (a triad of atoms \(i-j-k\), or an edge from Level 2) becomes a Node in this graph.

  • Edges: Two angle-nodes are connected if they share a common central bond (i.e., they form a dihedral geometry \(i-j-k-l\)).

  • Predictions:

    • Proper Dihedrals (\(C_0...V_5\)): A classification head operates on the Edges of this graph. Since an edge here represents the coupling between two angles (the torsion), this explicitly captures the 1-4 interaction environment required for accurate dihedral assignment.

Message Passing Flow

Information flows bi-directionally across these hierarchies. High-level geometric context (e.g., “this angle is constrained within a rigid ring”) propagates down to the atom level, while local electronic features (e.g., “this atom is highly electronegative”) propagate up to influence bond, angle, and dihedral predictions.

Physics-Informed Charge equilibration

A common issue with ML charge prediction is non-neutrality (e.g., total charge = +0.02e). ChemFAST enforces physical validity through a post-processing layer:

  • Total Charge Constraint: The predicted charges are globally corrected to ensure the sum equals the target integer charge of the molecule (usually 0, or \(\pm n\) for ions).

  • Symmetry Constraint: Atoms that are topologically symmetric (detected via graph automorphism) are forced to have identical charges.

4. User Control: The strategies Argument

While ChemFAST defaults to the ‘hybrid’ approach, users can explicitly control the assignment strategy for each molecule type via the Python API.

This is particularly useful for systems containing both simple ions (which must be exact) and complex polymers (which need ML).

# Example: Defining strategies for a mixed system
# mol_graphs contains [Polymer, Li+, TFSI-]

strategies = []
for mol in mol_graphs:
    if mol.number_of_nodes() > 20:
        # Complex Polymer -> Use ML (faster, handles novel connectivity)
        strategies.append('ml')
    elif mol.number_of_nodes() == 1:
        # Single Ion (Li+) -> Use Template (Must match standard OPLS ion params)
        strategies.append('tpl')
    else:
        # Small Molecule (TFSI) -> Try DB first, fallback to ML
        strategies.append('hybrid')

# Apply the strategies
ffs = assign_ff_parameters(mol_graphs, strategies=strategies)

Available Strategies:

  • 'hybrid' (Default): Try Database first; if failed, use ML. Recommended for most small molecules.

  • 'ml': Force usage of the GAT model. Recommended for large, irregular macromolecules where DB fragmentation might be slow or incomplete.

  • 'tpl' (Template): Strictly force Database lookup. Will raise an error if parameters are missing. Recommended for ions (Li+, Na+, Cl-) or standard solvents (Water, Methanol).

5. Scope and Force Field Compatibility

Target Force Field: OPLS-AA ChemFAST is currently optimized for the OPLS-AA (Optimized Potentials for Liquid Simulations - All Atom) force field family.

  • The internal database is built from OPLS/BOSS datasets.

  • The ML models are trained to reproduce OPLS/BOSS parameters.

Applicability Domains:

  • Well-Supported: Organic synthesis, polymers (thermoplastics, thermosets, conjugated systems), ionic liquids, organic solvents.

  • Requires Caution:

    • Biomolecules: While OPLS-AA/M works for proteins, specialized force fields (CHARMM/AMBER) are often preferred in bio-physics. ChemFAST can generate the topology, but users might want to validate parameters for highly specific protein residues.

    • Inorganics: ChemFAST assigns Lennard-Jones parameters to metals (e.g., Gold) based on interface compatibility, but it does not generate EAM potentials for bulk metal physics.