Skip to content

Units and Conversions

units

Canonical unit conversions for force field parameters.

ForceField stores all parameters in a canonical unit system:

============  ========================
Parameter     Canonical unit
============  ========================
bond_k        kcal/(mol·Å²)
bond_eq       Å
angle_k       kcal/(mol·rad²)
angle_eq      degrees
torsion_k     kcal/mol
vdw_radius    Å
vdw_epsilon   kcal/mol
============  ========================

Energy convention: E = k·(x − x₀)² (no ½ factor), matching AMBER.

Loaders convert from format-native units to canonical on load. Savers convert from canonical back to format-native on save. Engines convert from canonical to engine-native at the boundary.

mm3_bond_k_to_canonical

mm3_bond_k_to_canonical(k: float) -> float

Convert MM3 bond force constant (mdyn/Å) to canonical (kcal/mol/Ų).

Source code in q2mm/models/units.py
def mm3_bond_k_to_canonical(k: float) -> float:
    """Convert MM3 bond force constant (mdyn/Å) to canonical (kcal/mol/Ų)."""
    return k * MDYNA_TO_KCALMOLA2

canonical_to_mm3_bond_k

canonical_to_mm3_bond_k(k: float) -> float

Convert canonical bond force constant (kcal/mol/Ų) to MM3 (mdyn/Å).

Source code in q2mm/models/units.py
def canonical_to_mm3_bond_k(k: float) -> float:
    """Convert canonical bond force constant (kcal/mol/Ų) to MM3 (mdyn/Å)."""
    return k * KCALMOLA2_TO_MDYNA

mm3_angle_k_to_canonical

mm3_angle_k_to_canonical(k: float) -> float

Convert MM3 angle force constant (mdyn·Å/rad²) to canonical (kcal/mol/rad²).

Source code in q2mm/models/units.py
def mm3_angle_k_to_canonical(k: float) -> float:
    """Convert MM3 angle force constant (mdyn·Å/rad²) to canonical (kcal/mol/rad²)."""
    return k * MDYNA_RAD2_TO_KCALMOLRAD2

canonical_to_mm3_angle_k

canonical_to_mm3_angle_k(k: float) -> float

Convert canonical angle force constant (kcal/mol/rad²) to MM3 (mdyn·Å/rad²).

Source code in q2mm/models/units.py
def canonical_to_mm3_angle_k(k: float) -> float:
    """Convert canonical angle force constant (kcal/mol/rad²) to MM3 (mdyn·Å/rad²)."""
    return k * KCALMOLRAD2_TO_MDYNA_RAD2