use crate::error::ParseError;
#[derive(Clone, Copy, Debug, PartialEq)]
pub struct Dimension {
pub exponents: [f64; 5],
}
impl Dimension {
pub const ZERO: Self = Self {
exponents: [0.0; 5],
};
pub const LENGTH: Self = Self {
exponents: [1.0, 0.0, 0.0, 0.0, 0.0],
};
pub const TIME: Self = Self {
exponents: [0.0, 1.0, 0.0, 0.0, 0.0],
};
pub const MASS: Self = Self {
exponents: [0.0, 0.0, 1.0, 0.0, 0.0],
};
pub const CHARGE: Self = Self {
exponents: [0.0, 0.0, 0.0, 1.0, 0.0],
};
fn mul(self, other: Self) -> Self {
let mut e = [0.0; 5];
for (slot, (left, right)) in e
.iter_mut()
.zip(self.exponents.iter().zip(other.exponents.iter()))
{
*slot = left + right;
}
Self { exponents: e }
}
fn div(self, other: Self) -> Self {
let mut e = [0.0; 5];
for (slot, (left, right)) in e
.iter_mut()
.zip(self.exponents.iter().zip(other.exponents.iter()))
{
*slot = left - right;
}
Self { exponents: e }
}
fn pow(self, p: f64) -> Self {
let mut e = [0.0; 5];
for (slot, exp) in e.iter_mut().zip(self.exponents.iter()) {
*slot = exp * p;
}
Self { exponents: e }
}
fn compatible(self, other: Self) -> bool {
self.exponents
.iter()
.zip(other.exponents.iter())
.all(|(a, b)| (a - b).abs() < 1e-9)
}
}
#[derive(Clone, Copy, Debug)]
struct UnitValue {
si_factor: f64,
dim: Dimension,
}
const ENERGY: Dimension = Dimension {
exponents: [2.0, -2.0, 1.0, 0.0, 0.0],
};
struct BaseUnit {
name: &'static str,
value: UnitValue,
}
const fn unit(si_factor: f64, dim: Dimension) -> UnitValue {
UnitValue { si_factor, dim }
}
const BASE_UNITS: &[BaseUnit] = &[
BaseUnit {
name: "a",
value: unit(1e-10, Dimension::LENGTH),
},
BaseUnit {
name: "a0",
value: unit(5.291_772_109_03e-11, Dimension::LENGTH),
},
BaseUnit {
name: "amu",
value: unit(1.660_539_066_60e-27, Dimension::MASS),
},
BaseUnit {
name: "angstrom",
value: unit(1e-10, Dimension::LENGTH),
},
BaseUnit {
name: "bohr",
value: unit(5.291_772_109_03e-11, Dimension::LENGTH),
},
BaseUnit {
name: "da",
value: unit(1.660_539_066_60e-27, Dimension::MASS),
},
BaseUnit {
name: "dalton",
value: unit(1.660_539_066_60e-27, Dimension::MASS),
},
BaseUnit {
name: "e",
value: unit(1.602_176_634e-19, Dimension::CHARGE),
},
BaseUnit {
name: "electron_charge",
value: unit(1.602_176_634e-19, Dimension::CHARGE),
},
BaseUnit {
name: "ev",
value: unit(1.602_176_634e-19, ENERGY),
},
BaseUnit {
name: "femtosecond",
value: unit(1e-15, Dimension::TIME),
},
BaseUnit {
name: "fs",
value: unit(1e-15, Dimension::TIME),
},
BaseUnit {
name: "ha",
value: unit(4.359_744_722_207_1e-18, ENERGY),
},
BaseUnit {
name: "hartree",
value: unit(4.359_744_722_207_1e-18, ENERGY),
},
BaseUnit {
name: "j",
value: unit(1.0, ENERGY),
},
BaseUnit {
name: "joule",
value: unit(1.0, ENERGY),
},
BaseUnit {
name: "kcal",
value: unit(4184.0, ENERGY),
},
BaseUnit {
name: "kg",
value: unit(1.0, Dimension::MASS),
},
BaseUnit {
name: "kilogram",
value: unit(1.0, Dimension::MASS),
},
BaseUnit {
name: "kj",
value: unit(1000.0, ENERGY),
},
BaseUnit {
name: "m",
value: unit(1.0, Dimension::LENGTH),
},
BaseUnit {
name: "meter",
value: unit(1.0, Dimension::LENGTH),
},
BaseUnit {
name: "metre",
value: unit(1.0, Dimension::LENGTH),
},
BaseUnit {
name: "mev",
value: unit(1.602_176_634e-22, ENERGY),
},
BaseUnit {
name: "mol",
value: unit(6.022_140_76e23, Dimension::ZERO),
},
BaseUnit {
name: "nanometer",
value: unit(1e-9, Dimension::LENGTH),
},
BaseUnit {
name: "nanosecond",
value: unit(1e-9, Dimension::TIME),
},
BaseUnit {
name: "nm",
value: unit(1e-9, Dimension::LENGTH),
},
BaseUnit {
name: "ns",
value: unit(1e-9, Dimension::TIME),
},
BaseUnit {
name: "picosecond",
value: unit(1e-12, Dimension::TIME),
},
BaseUnit {
name: "ps",
value: unit(1e-12, Dimension::TIME),
},
BaseUnit {
name: "s",
value: unit(1.0, Dimension::TIME),
},
BaseUnit {
name: "second",
value: unit(1.0, Dimension::TIME),
},
BaseUnit {
name: "u",
value: unit(1.660_539_066_60e-27, Dimension::MASS),
},
BaseUnit {
name: "å",
value: unit(1e-10, Dimension::LENGTH),
},
];
#[cfg(test)]
pub(crate) fn cached_unit_table_addr() -> usize {
BASE_UNITS.as_ptr() as usize
}
fn lookup_base(name: &str) -> Option<&'static UnitValue> {
BASE_UNITS
.binary_search_by(|entry| entry.name.cmp(name))
.ok()
.map(|index| &BASE_UNITS[index].value)
}
pub fn parse_unit_expression(expr: &str) -> Result<(f64, Dimension), ParseError> {
let s: String = expr.chars().filter(|c| !c.is_whitespace()).collect();
if s.is_empty() {
return Err(ParseError::ValidationError("empty unit expression".into()));
}
parse_expr(&s.to_lowercase())
}
fn parse_expr(s: &str) -> Result<(f64, Dimension), ParseError> {
let mut factor = 1.0_f64;
let mut dim = Dimension::ZERO;
let mut i = 0;
let bytes = s.as_bytes();
let mut pending_div = false;
while i < bytes.len() {
if bytes[i] == b'*' {
i += 1;
continue;
}
if bytes[i] == b'/' {
pending_div = true;
i += 1;
continue;
}
let (atom_f, atom_d, consumed) = parse_atom(&s[i..])?;
i += consumed;
if pending_div {
factor /= atom_f;
dim = dim.div(atom_d);
pending_div = false;
} else {
factor *= atom_f;
dim = dim.mul(atom_d);
}
}
Ok((factor, dim))
}
fn parse_atom(s: &str) -> Result<(f64, Dimension, usize), ParseError> {
let bytes = s.as_bytes();
if bytes.is_empty() {
return Err(ParseError::ValidationError(
"trailing operator in unit".into(),
));
}
if bytes[0] == b'(' {
let mut depth = 0;
for (j, &b) in bytes.iter().enumerate() {
if b == b'(' {
depth += 1;
} else if b == b')' {
depth -= 1;
if depth == 0 {
let inner = &s[1..j];
let (f, d) = parse_expr(inner)?;
let mut end = j + 1;
let (f2, d2, extra) = apply_power(f, d, &s[end..])?;
end += extra;
return Ok((f2, d2, end));
}
}
}
return Err(ParseError::ValidationError("unbalanced '(' in unit".into()));
}
let mut j = 0;
while j < bytes.len()
&& (bytes[j].is_ascii_alphanumeric() || bytes[j] == b'_' || bytes[j] > 127)
{
j += 1;
}
if j == 0 {
return Err(ParseError::ValidationError(format!(
"expected unit name in '{s}'"
)));
}
let name = &s[..j];
let base = lookup_base(name)
.ok_or_else(|| ParseError::ValidationError(format!("unknown unit '{name}'")))?;
let mut end = j;
let (f, d, extra) = apply_power(base.si_factor, base.dim, &s[end..])?;
end += extra;
Ok((f, d, end))
}
fn apply_power(f: f64, d: Dimension, rest: &str) -> Result<(f64, Dimension, usize), ParseError> {
let bytes = rest.as_bytes();
if bytes.first() == Some(&b'^') {
let mut k = 1;
let neg = bytes.get(1) == Some(&b'-');
if neg {
k = 2;
}
let start = k;
while k < bytes.len() && (bytes[k].is_ascii_digit() || bytes[k] == b'.') {
k += 1;
}
if k == start {
return Err(ParseError::ValidationError(
"expected exponent after ^".into(),
));
}
let mut p: f64 = rest[start..k]
.parse()
.map_err(|_| ParseError::ValidationError("invalid unit exponent".into()))?;
if neg {
p = -p;
}
return Ok((f.powf(p), d.pow(p), k));
}
Ok((f, d, 0))
}
pub fn unit_conversion_factor(from_unit: &str, to_unit: &str) -> Result<f64, ParseError> {
let (f_from, d_from) = parse_unit_expression(from_unit)?;
let (f_to, d_to) = parse_unit_expression(to_unit)?;
if !d_from.compatible(d_to) {
return Err(ParseError::ValidationError(format!(
"incompatible units '{from_unit}' and '{to_unit}'"
)));
}
Ok(f_from / f_to)
}
pub fn velocity_file_time_unit(
spec_version: u32,
metadata: &serde_json::Map<String, serde_json::Value>,
override_unit: Option<&str>,
) -> Result<String, ParseError> {
if let Some(unit) = override_unit {
validate_unit_for_quantity("time", unit)?;
return Ok(unit.to_string());
}
let minor = metadata
.get(crate::types::meta::CON_SPEC_MINOR)
.and_then(|v| v.as_u64())
.unwrap_or(0);
if spec_version < 3 || minor == 0 {
return Ok(NATURAL_TIME_UNIT.to_string());
}
let time = metadata
.get(crate::types::meta::UNITS)
.and_then(|u| u.get("time"))
.and_then(|t| t.as_str())
.filter(|s| !s.is_empty());
match time {
Some(unit) => {
validate_unit_for_quantity("time", unit)?;
Ok(unit.to_string())
}
None => Err(ParseError::InvalidUnits(
"v3.1 velocities require units.time".into(),
)),
}
}
pub fn velocity_scale_to_natural(file_time: &str) -> Result<f64, ParseError> {
validate_unit_for_quantity("time", file_time)?;
let factor = unit_conversion_factor(file_time, NATURAL_TIME_UNIT)?;
if factor == 0.0 || !factor.is_finite() {
return Err(ParseError::InvalidUnits(format!(
"time unit '{file_time}' has a zero scale"
)));
}
Ok(1.0 / factor)
}
pub fn validate_unit_for_quantity(quantity: &str, unit: &str) -> Result<(), ParseError> {
let (_, dim) = parse_unit_expression(unit)?;
let expect = match quantity {
"length" => Dimension::LENGTH,
"time" => Dimension::TIME,
"mass" => Dimension::MASS,
"energy" => Dimension::MASS
.mul(Dimension::LENGTH.pow(2.0))
.div(Dimension::TIME.pow(2.0)),
"velocity" => Dimension::LENGTH.div(Dimension::TIME),
"force" => Dimension::MASS
.mul(Dimension::LENGTH)
.div(Dimension::TIME.pow(2.0)),
_ => {
return Err(ParseError::ValidationError(format!(
"unknown quantity '{quantity}'"
)));
}
};
if !dim.compatible(expect) {
return Err(ParseError::ValidationError(format!(
"unit '{unit}' is not valid for quantity '{quantity}'"
)));
}
Ok(())
}
pub fn validate_v3_units_metadata(units: &serde_json::Value) -> Result<(), ParseError> {
let obj = units
.as_object()
.ok_or_else(|| ParseError::ValidationError("units must be a JSON object".into()))?;
for key in ["length", "energy"] {
let Some(v) = obj.get(key) else {
return Err(ParseError::ValidationError(format!(
"v3 units must include non-empty '{key}'"
)));
};
let Some(s) = v.as_str() else {
return Err(ParseError::ValidationError(format!(
"units.{key} must be a string"
)));
};
if s.trim().is_empty() {
return Err(ParseError::ValidationError(format!(
"units.{key} must be non-empty"
)));
}
validate_unit_for_quantity(key, s)?;
}
for (key, qty) in [
("mass", "mass"),
("time", "time"),
("velocity", "velocity"),
("force", "force"),
] {
if let Some(v) = obj.get(key) {
let s = v.as_str().ok_or_else(|| {
ParseError::ValidationError(format!("units.{key} must be a string"))
})?;
validate_unit_for_quantity(qty, s)?;
}
}
Ok(())
}
pub const NATURAL_TIME_UNIT: &str = "(amu*angstrom^2/eV)^0.5";
pub fn default_v3_units_json() -> serde_json::Value {
serde_json::json!({
"length": "angstrom",
"energy": "eV",
"mass": "amu"
})
}
fn preferred_atom(alias: &str) -> Option<&'static str> {
let lower = alias.trim().to_lowercase();
Some(match lower.as_str() {
"a" | "å" | "angstrom" | "ångstrom" => "angstrom",
"nm" | "nanometer" | "nanometre" => "nm",
"m" | "meter" | "metre" => "m",
"bohr" | "a0" => "bohr",
"fs" | "femtosecond" | "femtoseconds" => "fs",
"ps" | "picosecond" | "picoseconds" => "ps",
"ns" | "nanosecond" | "nanoseconds" => "ns",
"s" | "sec" | "second" | "seconds" => "s",
"ev" => "eV",
"mev" => "meV",
"hartree" | "ha" => "hartree",
"j" | "joule" => "J",
"kj" => "kJ",
"kcal" => "kcal",
"amu" | "u" | "dalton" | "da" => "amu",
"kg" | "kilogram" => "kg",
"mol" => "mol",
_ => return None,
})
}
pub fn canonicalize_unit_expression(expr: &str) -> Result<String, ParseError> {
let trimmed = expr.trim();
if trimmed.is_empty() {
return Err(ParseError::ValidationError("empty unit".into()));
}
let _ = unit_conversion_factor(trimmed, trimmed)?;
let compact: String = trimmed.chars().filter(|c| !c.is_whitespace()).collect();
if compact == NATURAL_TIME_UNIT {
return Ok(NATURAL_TIME_UNIT.to_string());
}
let mut out = String::new();
let bytes = compact.as_bytes();
let mut i = 0;
while i < bytes.len() {
let b = bytes[i];
if matches!(b, b'*' | b'/' | b'(' | b')' | b'^') {
if !out.is_empty() && !out.ends_with(' ') && b != b')' && b != b'^' {
out.push(' ');
}
out.push(b as char);
if b != b'(' && b != b'^' {
out.push(' ');
}
i += 1;
continue;
}
if b == b'-' || b.is_ascii_digit() || b == b'.' {
out.push(b as char);
i += 1;
continue;
}
let start = i;
while i < bytes.len() && (bytes[i].is_ascii_alphabetic() || bytes[i] > 127) {
i += 1;
}
if start == i {
return Err(ParseError::ValidationError(format!(
"bad unit token in '{expr}'"
)));
}
let atom = &compact[start..i];
out.push_str(preferred_atom(atom).unwrap_or(atom));
}
let rendered = out.split_whitespace().collect::<Vec<_>>().join(" ");
let rendered_compact: String = rendered.chars().filter(|c| !c.is_whitespace()).collect();
if rendered_compact == NATURAL_TIME_UNIT {
return Ok(NATURAL_TIME_UNIT.to_string());
}
Ok(rendered)
}
pub fn canonicalize_units_object(
units: &serde_json::Value,
) -> Result<serde_json::Value, ParseError> {
let obj = units
.as_object()
.ok_or_else(|| ParseError::ValidationError("units must be a JSON object".into()))?;
let mut out = serde_json::Map::new();
for (k, v) in obj {
let s = v
.as_str()
.ok_or_else(|| ParseError::ValidationError(format!("units.{k} must be a string")))?;
out.insert(
k.clone(),
serde_json::Value::String(canonicalize_unit_expression(s)?),
);
}
Ok(serde_json::Value::Object(out))
}
pub fn canonicalize_units_metadata(
units: &serde_json::Value,
) -> Result<serde_json::Value, ParseError> {
validate_v3_units_metadata(units)?;
canonicalize_units_object(units)
}
#[cfg(test)]
mod tests {
use super::*;
#[test]
fn base_table_is_shared_across_calls_and_threads() {
for pair in BASE_UNITS.windows(2) {
assert!(
pair[0].name < pair[1].name,
"{} is not before {}",
pair[0].name,
pair[1].name
);
}
let energy = Dimension::MASS
.mul(Dimension::LENGTH.pow(2.0))
.div(Dimension::TIME.pow(2.0));
assert!(
parse_unit_expression("eV").unwrap().1.compatible(energy),
"energy exponents drifted from M L^2 T^-2"
);
let first = cached_unit_table_addr();
let second = cached_unit_table_addr();
assert_eq!(first, second);
let factors = std::thread::scope(|scope| {
let mut joins = Vec::new();
for _ in 0..8 {
joins.push(scope.spawn(|| parse_unit_expression("Angstrom").unwrap().0));
}
joins
.into_iter()
.map(|join| join.join().unwrap())
.collect::<Vec<_>>()
});
assert!(factors.iter().all(|factor| (*factor - 1e-10).abs() < 1e-22));
}
#[test]
fn angstrom_to_nm() {
let f = unit_conversion_factor("angstrom", "nm").unwrap();
assert!((f - 0.1).abs() < 1e-12);
}
#[test]
fn ev_to_mev() {
let f = unit_conversion_factor("eV", "meV").unwrap();
assert!((f - 1000.0).abs() < 1e-6);
}
#[test]
fn dimensional_mismatch() {
assert!(unit_conversion_factor("angstrom", "eV").is_err());
}
#[test]
fn ev_to_kcal_per_mol_is_metatomic() {
let f = unit_conversion_factor("eV", "kcal / mol").unwrap();
assert!((f - 23.060_548).abs() < 1e-3, "got {f}");
let force = unit_conversion_factor("eV / angstrom", "kJ / mol / angstrom").unwrap();
assert!((force - 96.485_332).abs() < 1e-3, "got {force}");
let ns = unit_conversion_factor("ns", "ps").unwrap();
assert!((ns - 1000.0).abs() < 1e-9, "got {ns}");
}
#[test]
fn validate_length_energy() {
validate_unit_for_quantity("length", "nm").unwrap();
validate_unit_for_quantity("energy", "hartree").unwrap();
assert!(validate_unit_for_quantity("energy", "angstrom").is_err());
}
#[test]
fn v3_units_require_length_energy() {
assert!(validate_v3_units_metadata(&serde_json::json!({"length": "A"})).is_err());
validate_v3_units_metadata(&default_v3_units_json()).unwrap();
}
#[test]
fn default_v3_units_do_not_invent_time_fs() {
let units = default_v3_units_json();
assert!(units.get("time").is_none(), "{units}");
assert_eq!(units["length"], "angstrom");
assert_eq!(units["energy"], "eV");
assert_eq!(units["mass"], "amu");
}
#[test]
fn canonicalize_aliases() {
assert_eq!(canonicalize_unit_expression("A").unwrap(), "angstrom");
assert_eq!(canonicalize_unit_expression("ev").unwrap(), "eV");
assert_eq!(canonicalize_unit_expression("femtosecond").unwrap(), "fs");
assert_eq!(
canonicalize_unit_expression("eV/angstrom").unwrap(),
"eV / angstrom"
);
let u = canonicalize_units_metadata(&serde_json::json!({
"length": "A",
"energy": "ev",
"time": "femtosecond"
}))
.unwrap();
assert_eq!(u["length"], "angstrom");
assert_eq!(u["energy"], "eV");
assert_eq!(u["time"], "fs");
}
#[test]
fn natural_time_spelling_is_pinned() {
let spaced = "( amu * angstrom ^ 2 / eV ) ^ 0.5";
assert_eq!(
canonicalize_unit_expression(spaced).unwrap(),
NATURAL_TIME_UNIT
);
assert!(!NATURAL_TIME_UNIT.contains(' '));
let factor = unit_conversion_factor(NATURAL_TIME_UNIT, "fs").unwrap();
assert!((factor - 10.180_505_710_759_414).abs() < 1e-9, "{factor}");
let mixed = "(AMU * angstrom^2 / ev)^0.5";
assert_eq!(
canonicalize_unit_expression(mixed).unwrap(),
NATURAL_TIME_UNIT
);
}
}