chematic

Pure-Rust cheminformatics library for Python — SMILES parsing, descriptors, fingerprints, substructure search, formats, reactions, and bounded 3D tools.

The current published line is 1.1.3. It includes molecule pickling, indexed SDF access, V3000 bond-stereo collections, and reusable Python MMFF94/UFF force fields. Compatibility and performance claims remain operation-, artifact-, and corpus-scoped; see the validation and migration guides before replacing a production RDKit workflow.

Installation

pip install chematic

Quick Start

import chematic

mol = chematic.from_smiles("CC(=O)Oc1ccccc1C(=O)O")  # aspirin

print(mol.mw)      # 180.16
print(mol.logp)    # 1.31
print(mol.tpsa)    # 63.6
print(mol.qed)     # 0.55

# New descriptors
print(mol.vabc)              # van der Waals volume (no 3D needed)
print(mol.schultz_mti)       # Schultz MTI
print(mol.gutman_mti)        # Gutman MTI*
print(mol.gravitational_index)  # gravitational index

# pKa prediction
print(mol.pka())   # {"most_acidic": 3.49, "most_basic": None}

# ADMET profile
print(mol.admet())
# {"bbb": False, "bbb_score": ..., "caco2": ..., "herg_risk": ..., "cyp3a4_risk": ...}

# Fingerprints (bytes, 2048-bit ECFP4)
fp = mol.ecfp4()

# Tanimoto similarity
mol2 = chematic.from_smiles("c1ccccc1")
sim = chematic.tanimoto(mol.ecfp4(), mol2.ecfp4())

# Mol objects can be cached or sent to worker processes with pickle. Graph
# chemistry, stereo, R-groups, source metadata, and SD properties survive.
import pickle
restored = pickle.loads(pickle.dumps(mol))

# Natural-language property summary (for LLM / MCP agents)
print(mol.describe())

# Structural diff between two molecules
ibuprofen = chematic.from_smiles("CC(C)Cc1ccc(CC(C)C(=O)O)cc1")
d = mol.diff(ibuprofen)  # {"summary": "...", "delta_mw": 66.1, "delta_logp": 2.75, ...}

# SVG / PDF / EPS depiction
svg = mol.to_svg()
pdf_bytes = mol.to_pdf()   # bytes; requires pdf feature
eps_str   = mol.to_eps()   # PostScript string

# ChemicalJSON (Avogadro 2 / MolSSI)
cjson_str = mol.to_cjson(coords=[])   # coords: list of (x,y,z) tuples, optional
mol2, coords = chematic.from_cjson(cjson_str)

# Template-based retrosynthesis (60 retro-SMIRKS templates)
mol3 = chematic.from_smiles("CC(=O)Nc1ccccc1")  # acetanilide
results = mol3.retro_disconnect(max_results=5)
for r in results:
    print(r["template"], "→", r["precursors"])
# amide_secondary → ['CC(=O)O', 'Nc1ccccc1']

# Filter by reaction class
amides = mol3.retro_disconnect(reaction_class="AmideBond")

# Bulk substructure match against a pre-parsed Mol list (returns indices)
mols = [chematic.from_smiles(s) for s in ["CCO", "c1ccccc1O", "CC(=O)O"]]
hits = chematic.bulk.substructure_match("[OH]", mols)  # → [0, 1, 2]

# All descriptors as a dict (for Pandas)
import pandas as pd
smiles = ["CCO", "c1ccccc1", "CC(=O)O"]
df = pd.DataFrame([chematic.from_smiles(s).descriptors() for s in smiles])

# Opt-in v2 embedding pipeline: torsion-knowledge-aware distance geometry +
# stereo verification/repair + policy-gated force field, with full per-stage
# evidence (never just final coordinates)
config = chematic.PipelineV2Config.safe(
    force_field="mmff94_with_uff_fallback",
    stereo_policy="repair_and_verify",
    ring_torsion_policy="fail_closed",
)
try:
    result = mol.embed_pipeline_v2(config)
    coords = result["coords"]                       # same atom order as mol
    print(result["force_field"]["actual_force_field_used"])  # fallback if MMFF94 lacked params
    print(result["final_validation"]["sound"])
except chematic.PipelineV2Error as e:
    print(e.diagnostics["stage"], e.diagnostics["cause"])   # structured, not just a message

# Reuse prepared MMFF94 terms across iterative energy/gradient calls
field = mol.mmff_force_field(coords)
energy, gradient = field.energy_gradient()
field.set_position_(0, [coords[0][0] + 0.01, *coords[0][1:]])
status, optimized_energy = field.minimize_()

Molecular GNN input

Mol.ogb_graph_features() returns the OGB categorical graph schema as NumPy int64 arrays. It works with molecules from SMILES, SDF records, and MOL2. The keys are node_feat (N, 9), edge_index (2, 2E), edge_feat (2E, 3), and num_nodes. The edges and their features have the same order, with one edge in each direction per bond. PyTorch and PyG are optional:

import chematic

for source_index, record, error in chematic.iter_sdf_checked("ligands.sdf"):
    if error is not None:
        print(f"Rejected SDF record {source_index}: {error}")
        continue
    assert record is not None
    graph = record.ogb_graph_features()
    # Attach record.name and record.properties() as labels or IDs as needed.
    # With PyTorch Geometric installed:
    # import torch
    # from torch_geometric.data import Data
    # data = Data(
    #     x=torch.from_numpy(graph["node_feat"]),
    #     edge_index=torch.from_numpy(graph["edge_index"]),
    #     edge_attr=torch.from_numpy(graph["edge_feat"]),
    #     num_nodes=graph["num_nodes"],
    # )

The default removes RDKit-removable hydrogen graph atoms, following OGB's default SMILES conversion. Pass remove_hs=False to retain the atoms already present in the molecule; this option does not add hydrogens. Unsupported chemistry raises ValueError. This export contains categorical graph features, not coordinates or SD properties. iter_sdf_checked reports a rejected record at its source index. Its default reader limits are 1 GiB and 100,000 records; pass positive max_input_bytes and max_records overrides for larger files. Use SdfRecord.ogb_graph_features() for V2000 or V3000 SDF input: it retains explicitly unknown double-bond stereo (V2000 code 3 or V3000 CFG=2). The generic Mol.ogb_graph_features() has no record-level stereo diagnostic to read. For a 3D SDF, chirality inferred by RDKit from coordinates can differ from the chirality retained in record.mol from the CTAB wedge representation.

Features

• Zero C/C++ dependencies — pure Rust, no RDKit or OpenBabel required
• SMILES / MOL / SDF / ChemicalJSON parsing and writing
• Indexed SDF datasets — SdfDataset provides cheap record counts, metadata, random access, and bounded range reads without materializing every molecule
• 190+ descriptor values (70+ functions; several return vectors): MW, LogP, TPSA, QED, Fsp3, SA Score, HBD, vabc, schultz_mti, gutman_mti, and gravitational_index. RDKit agreement is metric-specific; see the validation report.
• 14 fingerprint algorithms: ECFP2/4/6, FCFP4/6, MACCS, AtomPair, Torsion, …
• pKa prediction (15 SMARTS rules — unique to chematic)
• ADMET profile: BBB, Caco-2, hERG, CYP3A4
• Template-based retrosynthesis: mol.retro_disconnect() — 60 retro-SMIRKS templates, SA Score ranked
• SMARTS substructure search — chematic.smarts_match() and bulk.substructure_match(smarts, mols) (pre-parsed Mol list, returns indices)
• SVG / PDF / EPS depiction: mol.to_svg(), mol.to_pdf(), mol.to_eps()
• ChemicalJSON: mol.to_cjson(coords=[]), chematic.from_cjson(s) — Avogadro 2 / MolSSI compatible
• Opt-in v2 embedding pipeline: mol.embed_pipeline_v2(config) — torsion-knowledge-aware distance geometry, stereo verify/repair, and policy-gated force field (PipelineV2Config), returning full per-stage evidence (embed stats, torsion knowledge/optimization reports, stereo before/after, force-field actual policy and fallback, final geometry validation, stage timings) instead of just coordinates; raises chematic.PipelineV2Error with structured .diagnostics on failure
• Reusable MMFF94/UFF objects: mol.mmff_force_field(coords) and mol.uff_force_field(coords) retain prepared terms across energy, gradient, coordinate-update, and in-place minimization calls
• Charge-conserving EEM: mol.eem_charges(coords) uses the published NEEMP model and preserves the formal-charge sum for ions

RDKit compatibility

The operation/profile boundaries and measured oracle lanes are listed in the Compatibility Contract dashboard. For persisted fingerprints, saved indices, canonical identity, stereo, and browser/Worker migration decisions, read the RDKit migration guide before changing a production workflow. Native ECFP and RDKit-compatible Morgan are separate profiles; rebuild stored fingerprints and indices under the chosen profile.

chematic.rdkit_compat provides a lightweight RDKit-compatible subset for environments where RDKit is unavailable (WASM, serverless, conda-free CI):

from chematic import rdkit_compat as Chem
from chematic.rdkit_compat import Descriptors, rdMolDescriptors, DataStructs

mol = Chem.MolFromSmiles("CC(=O)Oc1ccccc1C(=O)O")

# Descriptors
Descriptors.MolWt(mol)          # 180.16
rdMolDescriptors.CalcTPSA(mol)  # 63.6

# Fingerprint (ExplicitBitVect) with bitInfo
bitInfo = {}
fp = rdMolDescriptors.GetMorganFingerprintAsBitVect(mol, 2, nBits=2048, bitInfo=bitInfo)
fp.GetNumBits()                         # 2048
bitInfo                                 # {bit: ((atom_idx, radius), ...)}
DataStructs.TanimotoSimilarity(fp, fp)  # 1.0
DataStructs.BulkTanimotoSimilarity(fp, [fp])  # [1.0]

import numpy as np
arr = DataStructs.ConvertToNumpyArray(fp)  # (2048,) int8 for sklearn / PyTorch

# Explain RDKit atom-pair / topological-torsion codes
code = mol.rdkit_atom_code(0)
chematic.explain_rdkit_atom_code(code)  # (symbol, branches, pi_electrons)
path_code = mol.rdkit_topological_torsion_path_code([0, 1, 2, 3])
chematic.explain_rdkit_topological_torsion_code(path_code)

# Atom / Bond traversal
for atom in mol.GetAtoms():
    atom.GetSymbol(), atom.GetAtomicNum(), atom.IsInRing()
for bond in mol.GetBonds():
    bond.GetBondType(), bond.GetBondTypeAsDouble(), bond.IsInRing()

# Ring information
ri = mol.GetRingInfo()
ri.NumRings()       # 1
ri.AtomRings()      # tuple of tuples of atom indices
ri.NumAtomRings(0)  # rings containing atom 0

# SDF I/O with SD properties
with Chem.SDWriter("out.sdf") as w:
    mol.SetProp("ID", "aspirin")
    w.write(mol)
for m in Chem.SDMolSupplier("out.sdf"):
    print(m.GetProp("ID"))

# Large seekable SDF: index once, parse records only when requested
from chematic import SdfDataset
dataset = SdfDataset.open("out.sdf")
print(len(dataset), dataset.metadata(0).title)
print(dataset[-1].smiles)

Unsupported options raise NotImplementedError or TypeError — they are never silently ignored.

Compatibility matrix

A live differential suite (tests/test_rdkit_diff.py, auto-skipped when RDKit is absent) compares chematic against RDKit across descriptors, ring counts, SMARTS match counts, SDF round-trips, and Morgan self-similarity, writing an explainable diff to validation/results/rdkit_diff.jsonl.

chematic.rdkit_compat is not a full RDKit clone — it is a lightweight RDKit-compatible subset for common 2D cheminformatics workflows. See the full RDKit migration guide (compatibility matrix, differential-validation results, known divergences, and runnable examples).

License

MIT OR Apache-2.0