use crate::types::{
ConFrame, SECTION_CHARGES, SECTION_DISPLACEMENTS, SECTION_ENERGIES, SECTION_FORCES,
SECTION_MAGMOMS, SECTION_SPINS, SECTION_SPREADS, SECTION_VELOCITIES, encode_fixed_bitmask,
meta,
};
use serde_json::{Value, json};
use std::collections::BTreeMap;
use std::fs::{File, OpenOptions};
use std::io::{self, BufWriter, Read, Seek, SeekFrom, Write};
use std::path::Path;
#[derive(Debug)]
pub enum WriteError {
Io(io::Error),
Format(crate::error::ParseError),
}
impl std::fmt::Display for WriteError {
fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
match self {
WriteError::Io(e) => write!(f, "{e}"),
WriteError::Format(e) => write!(f, "{e}"),
}
}
}
impl std::error::Error for WriteError {
fn source(&self) -> Option<&(dyn std::error::Error + 'static)> {
match self {
WriteError::Io(e) => Some(e),
WriteError::Format(e) => Some(e),
}
}
}
impl From<io::Error> for WriteError {
fn from(e: io::Error) -> Self {
WriteError::Io(e)
}
}
impl From<WriteError> for io::Error {
fn from(e: WriteError) -> Self {
match e {
WriteError::Io(e) => e,
WriteError::Format(err) => io::Error::new(io::ErrorKind::InvalidData, err.to_string()),
}
}
}
#[derive(Clone, Copy, Debug, PartialEq, Eq, Default)]
pub enum FloatFormat {
DecimalPlaces(usize),
#[default]
RoundTrip,
}
pub struct ConFrameWriter<W: Write> {
writer: BufWriter<W>,
float_format: FloatFormat,
canonical: bool,
metadata_cache: Option<MetadataCacheEntry>,
scratch: Vec<u8>,
ryu: ryu::Buffer,
fmt_slot: [FmtSlot; 3],
atom_workers: usize,
bound: Option<BoundHeader>,
stream_file: Option<File>,
force_map: bool,
keep_scratch: bool,
exact_scratch: bool,
keep_bound: bool,
bound_vel: Option<bool>,
bound_flags: Option<[bool; 8]>,
bound_line: Option<String>,
bound_no_lattice: Option<bool>,
}
struct BoundHeader {
meta: BTreeMap<String, Value>,
spec: u32,
}
#[derive(Clone, Copy)]
struct FmtSlot {
bits: u64,
text: [u8; 24],
len: u8,
state: u8,
}
impl FmtSlot {
const fn new() -> Self {
Self {
bits: 0,
text: [0; 24],
len: 0,
state: 0,
}
}
}
#[derive(Debug)]
struct MetadataCacheEntry {
spec_version: u32,
has_velocities: bool,
has_forces: bool,
has_energies: bool,
has_charges: bool,
has_spins: bool,
has_magmoms: bool,
has_displacements: bool,
has_spreads: bool,
metadata: std::collections::BTreeMap<String, serde_json::Value>,
serialized: String,
}
impl MetadataCacheEntry {
fn matches(
&self,
spec_version: u32,
sections: [bool; 8],
metadata: &std::collections::BTreeMap<String, serde_json::Value>,
) -> bool {
let [
has_velocities,
has_forces,
has_energies,
has_charges,
has_spins,
has_magmoms,
has_displacements,
has_spreads,
] = sections;
self.spec_version == spec_version
&& self.has_velocities == has_velocities
&& self.has_forces == has_forces
&& self.has_energies == has_energies
&& self.has_charges == has_charges
&& self.has_spins == has_spins
&& self.has_magmoms == has_magmoms
&& self.has_displacements == has_displacements
&& self.has_spreads == has_spreads
&& self.metadata.len() == metadata.len()
&& self
.metadata
.iter()
.zip(metadata)
.all(|((ka, a), (kb, b))| ka == kb && metadata_value_eq(a, b))
}
}
fn metadata_value_eq(a: &serde_json::Value, b: &serde_json::Value) -> bool {
use serde_json::Value;
match (a, b) {
(Value::Number(a), Value::Number(b)) if a.is_f64() && b.is_f64() => {
a.as_f64().map(f64::to_bits) == b.as_f64().map(f64::to_bits)
}
(Value::Array(a), Value::Array(b)) => {
a.len() == b.len() && a.iter().zip(b).all(|(a, b)| metadata_value_eq(a, b))
}
(Value::Object(a), Value::Object(b)) => {
a.len() == b.len()
&& a.iter()
.all(|(key, a)| b.get(key).is_some_and(|b| metadata_value_eq(a, b)))
}
_ => a == b,
}
}
fn require_minor2_section_arity(
minor: u64,
sections: &[serde_json::Value],
metadata: &serde_json::Map<String, serde_json::Value>,
) -> Result<(), WriteError> {
if minor < 2 {
return Ok(());
}
for name in sections.iter().filter_map(|value| value.as_str()) {
if crate::types::reserved_section_value_arity(name).is_some() {
continue;
}
let listed = metadata
.get(meta::SECTION_ARITY)
.and_then(|value| value.as_object())
.and_then(|obj| obj.get(name))
.and_then(|value| value.as_u64())
.is_some_and(|n| n >= 1);
if !listed {
return Err(WriteError::Format(
crate::error::ParseError::ValidationError(format!(
"minor >= 2 must list {name} in section_arity"
)),
));
}
}
Ok(())
}
struct StreamOut {
file: File,
start: u64,
armed: bool,
}
impl Drop for StreamOut {
fn drop(&mut self) {
if !self.armed {
return;
}
let _ = self.file.set_len(self.start);
let _ = self.file.seek(SeekFrom::Start(self.start));
}
}
impl<W: Write> ConFrameWriter<W> {
pub fn new(writer: W) -> Self {
Self::with_float_format(writer, FloatFormat::default())
}
pub fn with_precision(writer: W, precision: usize) -> Self {
Self::with_float_format(writer, FloatFormat::DecimalPlaces(precision))
}
pub fn with_float_format(writer: W, float_format: FloatFormat) -> Self {
Self {
writer: BufWriter::new(writer),
float_format,
canonical: false,
metadata_cache: None,
scratch: Vec::with_capacity(16 * 1024),
ryu: ryu::Buffer::new(),
fmt_slot: [FmtSlot::new(); 3],
atom_workers: 1,
bound: None,
stream_file: None,
force_map: false,
keep_scratch: false,
exact_scratch: false,
keep_bound: false,
bound_vel: None,
bound_flags: None,
bound_line: None,
bound_no_lattice: None,
}
}
fn retain_scratch(mut self, on: bool) -> Self {
self.keep_scratch = on;
self
}
pub fn parallel_atoms(mut self, n: usize) -> Self {
self.atom_workers = n.max(1);
self
}
pub fn map_rows(mut self, on: bool) -> Self {
self.force_map = on;
self
}
pub fn mapped_file(mut self, file: File) -> Self {
self.stream_file = Some(file);
self
}
pub fn canonical(mut self, on: bool) -> Self {
self.set_canonical(on);
self
}
pub fn set_canonical(&mut self, on: bool) {
self.canonical = on;
if on {
self.metadata_cache = None;
}
}
pub fn is_canonical(&self) -> bool {
self.canonical
}
fn refresh_metadata_cache(
&mut self,
spec_version: u32,
metadata: &BTreeMap<String, Value>,
has_vel: bool,
has_frc: bool,
has_eng: bool,
has_chg: bool,
has_spn: bool,
has_mm: bool,
has_dsp: bool,
has_spr: bool,
) -> Result<(), WriteError> {
let cache_hit = !self.canonical
&& self.metadata_cache.as_ref().is_some_and(|c| {
c.matches(
spec_version,
[
has_vel, has_frc, has_eng, has_chg, has_spn, has_mm, has_dsp, has_spr,
],
metadata,
)
});
if cache_hit {
return Ok(());
}
let mut meta_obj = serde_json::Map::new();
meta_obj.insert(meta::CON_SPEC_VERSION.into(), json!(spec_version));
let mut sections = Vec::new();
if has_vel {
sections.push(json!(SECTION_VELOCITIES));
}
if has_frc {
sections.push(json!(SECTION_FORCES));
}
if has_eng {
sections.push(json!(SECTION_ENERGIES));
}
if has_chg {
sections.push(json!(SECTION_CHARGES));
}
if has_spn {
sections.push(json!(SECTION_SPINS));
}
if has_mm {
sections.push(json!(SECTION_MAGMOMS));
}
if has_dsp {
sections.push(json!(SECTION_DISPLACEMENTS));
}
if has_spr {
sections.push(json!(SECTION_SPREADS));
}
let validate = metadata
.get(meta::VALIDATE)
.and_then(|value| value.as_bool())
.unwrap_or(false);
if !sections.is_empty() || validate {
meta_obj.insert(meta::SECTIONS.into(), json!(sections));
}
for (k, v) in metadata {
if k == meta::CON_SPEC_VERSION || k == meta::SECTIONS {
continue;
}
meta_obj.insert(k.clone(), v.clone());
}
if let Some(u) = meta_obj.get(meta::UNITS).cloned() {
let c = crate::units::canonicalize_units_object(&u).map_err(|e| {
WriteError::Format(crate::error::ParseError::InvalidUnits(e.to_string()))
})?;
meta_obj.insert(meta::UNITS.into(), c);
}
if spec_version >= 3 {
match meta_obj.get(meta::UNITS) {
None => {
return Err(WriteError::Format(crate::error::ParseError::InvalidUnits(
"con_spec_version >= 3 requires units.length and units.energy".into(),
)));
}
Some(u) => {
crate::units::validate_v3_units_metadata(u).map_err(|e| {
WriteError::Format(crate::error::ParseError::InvalidUnits(e.to_string()))
})?;
}
}
}
let minor = metadata
.get(meta::CON_SPEC_MINOR)
.and_then(|v| v.as_u64())
.unwrap_or(0);
if spec_version >= 3 && minor >= 1 && has_vel {
let time_ok = meta_obj
.get(meta::UNITS)
.and_then(|u| u.get("time"))
.and_then(|v| v.as_str())
.is_some_and(|s| !s.trim().is_empty());
if !time_ok {
return Err(WriteError::Format(crate::error::ParseError::InvalidUnits(
"minor >= 1 requires units.time when a velocities section is present".into(),
)));
}
}
require_minor2_section_arity(minor, §ions, &meta_obj)?;
let serialized = serde_json::Value::Object(meta_obj).to_string();
self.metadata_cache = Some(MetadataCacheEntry {
spec_version,
has_velocities: has_vel,
has_forces: has_frc,
has_energies: has_eng,
has_charges: has_chg,
has_spins: has_spn,
has_magmoms: has_mm,
has_displacements: has_dsp,
has_spreads: has_spr,
metadata: metadata.clone(),
serialized,
});
Ok(())
}
pub fn write_frame_binding(&mut self, frame: &mut ConFrame) -> Result<(), WriteError> {
self.write_frame_bound(frame)
}
pub fn write_frame_bound(&mut self, frame: &ConFrame) -> Result<(), WriteError> {
let spec = frame.header.spec_version.max(crate::CON_SPEC_VERSION);
let flags = [
frame.has_velocities(),
frame.has_forces(),
frame.has_energies(),
frame.has_charges(),
frame.has_spins(),
frame.has_magmoms(),
frame.has_displacements(),
frame.has_spreads(),
];
let line = frame.header.prebox_header.metadata_line();
let reuse = self.exact_scratch
&& self.bound.is_some()
&& self.bound.as_ref().is_some_and(|item| item.spec == spec)
&& self.bound_flags == Some(flags)
&& self.bound_line.as_deref() == Some(line);
if !reuse {
self.bound = Some(BoundHeader {
meta: binding_metadata(frame),
spec,
});
if self.exact_scratch {
self.bound_flags = Some(flags);
self.bound_vel = Some(flags[0]);
self.bound_line = Some(line.to_string());
self.bound_no_lattice = None;
}
}
self.keep_bound = self.exact_scratch;
self.write_frame(frame)
}
pub fn write_frame(&mut self, frame: &ConFrame) -> Result<(), WriteError> {
let zero_frame = frame.header.natm_types == 0
&& frame.atom_data.is_empty()
&& frame.header.natms_per_type.is_empty()
&& frame.header.masses_per_type.is_empty();
if zero_frame && !frame.header.sections.is_empty() {
return Err(WriteError::Format(crate::error::ParseError::InvalidCount(
"a zero-atom frame cannot declare sections".into(),
)));
}
if !zero_frame
&& (frame.header.natm_types == 0
|| frame.atom_data.is_empty()
|| frame.header.natms_per_type.contains(&0))
{
return Err(WriteError::Format(crate::error::ParseError::InvalidCount(
"each per-type count must be >= 1".into(),
)));
}
#[cfg(feature = "parallel")]
let parallel_rows = self.atom_workers >= 2 && frame.atom_data.len() >= WIDE_ATOM_GRAIN;
#[cfg(not(feature = "parallel"))]
let parallel_rows = false;
if !self.exact_scratch {
if parallel_rows {
if !ids_are_dense(&frame.atom_data, self.atom_workers) {
ensure_unique_ids(&frame.atom_data, self.atom_workers)?;
}
} else {
ensure_unique_ids(&frame.atom_data, 1)?;
}
}
let prec = self.float_format;
let stream_file = self.stream_file.take();
let keep_bound = self.keep_bound;
let bound = self.bound.take();
let spec_version = bound
.as_ref()
.map(|item| item.spec)
.unwrap_or(frame.header.spec_version);
let skip_refresh = keep_bound && self.metadata_cache.is_some();
if !skip_refresh {
let owned_meta;
let metadata = if let Some(item) = bound.as_ref() {
&item.meta
} else {
owned_meta = frame.header.metadata.to_map();
&owned_meta
};
if spec_version >= 2 {
self.refresh_metadata_cache(
spec_version,
metadata,
frame.has_velocities(),
frame.has_forces(),
frame.has_energies(),
frame.has_charges(),
frame.has_spins(),
frame.has_magmoms(),
frame.has_displacements(),
frame.has_spreads(),
)?;
}
}
if keep_bound {
self.bound = bound;
}
let (boxl, angles) = if self.keep_bound && self.bound_no_lattice == Some(true) {
(frame.header.boxl, frame.header.angles)
} else {
match frame.header.lattice_vectors() {
Some(vecs) => {
self.bound_no_lattice = Some(false);
crate::types::lengths_and_angles(vecs)
}
None => {
self.bound_no_lattice = Some(true);
(frame.header.boxl, frame.header.angles)
}
}
};
let atoms = frame.atom_data.len();
let [
has_vel,
has_frc,
has_eng,
has_chg,
has_spn,
has_mm,
has_dsp,
has_spr,
] = if self.keep_bound {
self.bound_flags.unwrap_or([
frame.has_velocities(),
frame.has_forces(),
frame.has_energies(),
frame.has_charges(),
frame.has_spins(),
frame.has_magmoms(),
frame.has_displacements(),
frame.has_spreads(),
])
} else {
[
frame.has_velocities(),
frame.has_forces(),
frame.has_energies(),
frame.has_charges(),
frame.has_spins(),
frame.has_magmoms(),
frame.has_displacements(),
frame.has_spreads(),
]
};
let mut blocks = 1usize;
if has_vel {
blocks += 1;
}
if has_frc {
blocks += 1;
}
if has_eng {
blocks += 1;
}
if has_chg {
blocks += 1;
}
if has_spn {
blocks += 1;
}
if has_mm {
blocks += 1;
}
if has_dsp {
blocks += 1;
}
if has_spr {
blocks += 1;
}
let mut stream_out = match stream_file {
Some(mut file) => {
let start = file.stream_position().map_err(WriteError::Io)?;
Some(StreamOut {
file,
start,
armed: true,
})
}
None => None,
};
if !self.exact_scratch {
self.scratch
.reserve(256 + atoms.saturating_mul(64).saturating_mul(blocks));
}
if !self.keep_scratch {
self.scratch.clear();
}
push_str_line(&mut self.scratch, &frame.header.prebox_header.user);
if spec_version < 2 {
push_str_line(
&mut self.scratch,
frame.header.prebox_header.metadata_line(),
);
} else {
let line = self
.metadata_cache
.as_ref()
.expect("metadata_cache populated above")
.serialized
.as_bytes();
self.scratch.extend_from_slice(line);
self.scratch.push(b'\n');
}
{
let buf = &mut self.scratch;
let ryu = &mut self.ryu;
let slots = &mut self.fmt_slot;
push_f64(buf, ryu, boxl[0], prec);
buf.push(b' ');
push_f64(buf, ryu, boxl[1], prec);
buf.push(b' ');
push_f64(buf, ryu, boxl[2], prec);
buf.push(b'\n');
push_f64(buf, ryu, angles[0], prec);
buf.push(b' ');
push_f64(buf, ryu, angles[1], prec);
buf.push(b' ');
push_f64(buf, ryu, angles[2], prec);
buf.push(b'\n');
push_str_line(buf, &frame.header.postbox_header[0]);
push_str_line(buf, &frame.header.postbox_header[1]);
push_u64(buf, frame.header.natm_types as u64);
buf.push(b'\n');
for (i, n) in frame.header.natms_per_type.iter().enumerate() {
if i > 0 {
buf.push(b' ');
}
push_u64(buf, *n as u64);
}
buf.push(b'\n');
for (i, m) in frame.header.masses_per_type.iter().enumerate() {
if i > 0 {
buf.push(b' ');
}
push_f64(buf, ryu, *m, prec);
}
buf.push(b'\n');
let save_velocity = has_vel;
let direct_velocity = save_velocity
&& self.atom_workers < 2
&& !self.force_map
&& frame
.header
.natms_per_type
.iter()
.all(|count| *count < WIDE_ATOM_GRAIN);
std::thread_local! {
static VEL_ROWS: std::cell::RefCell<Vec<Vec<u8>>> =
const { std::cell::RefCell::new(Vec::new()) };
}
struct RestoreRows(Vec<Vec<u8>>);
impl Drop for RestoreRows {
fn drop(&mut self) {
let rows = std::mem::take(&mut self.0);
VEL_ROWS.with(|cell| *cell.borrow_mut() = rows);
}
}
let mut velocity_rows = if direct_velocity {
None
} else {
Some(RestoreRows(VEL_ROWS.with(|cell| cell.take())))
};
if let Some(velocity_rows) = velocity_rows.as_mut() {
let ntypes = frame.header.natms_per_type.len();
if velocity_rows.0.len() < ntypes {
velocity_rows.0.resize_with(ntypes, Vec::new);
}
for row in &mut velocity_rows.0 {
row.clear();
}
}
let mut atom_idx_offset = 0;
for (type_idx, &num_atoms_in_type) in frame.header.natms_per_type.iter().enumerate() {
let symbol = &frame.atom_data[atom_idx_offset].symbol;
push_str_line(buf, symbol);
push_component_line(buf, "Coordinates of Component", type_idx);
let mapped = self.force_map
|| (self.atom_workers >= 2 && num_atoms_in_type >= WIDE_ATOM_GRAIN);
append_coord_rows(
buf,
ryu,
slots,
&frame.atom_data[atom_idx_offset..atom_idx_offset + num_atoms_in_type],
prec,
spec_version,
self.atom_workers,
mapped,
self.force_map,
stream_out.as_mut().map(|out| &mut out.file),
if save_velocity && !direct_velocity {
Some(&mut velocity_rows.as_mut().unwrap().0[type_idx])
} else {
None
},
)?;
atom_idx_offset += num_atoms_in_type;
}
if has_vel {
buf.push(b'\n');
let mut vel_idx_offset = 0;
for (type_idx, &num_atoms_in_type) in frame.header.natms_per_type.iter().enumerate()
{
let symbol = &frame.atom_data[vel_idx_offset].symbol;
push_str_line(buf, symbol);
push_component_line(buf, "Velocities of Component", type_idx);
let atoms =
&frame.atom_data[vel_idx_offset..vel_idx_offset + num_atoms_in_type];
if direct_velocity {
append_velocity_lines(buf, ryu, atoms, prec, spec_version);
} else {
buf.extend_from_slice(&velocity_rows.as_ref().unwrap().0[type_idx]);
}
vel_idx_offset += num_atoms_in_type;
}
}
if has_frc {
buf.push(b'\n');
let mut force_idx_offset = 0;
for (type_idx, &num_atoms_in_type) in frame.header.natms_per_type.iter().enumerate()
{
let symbol = &frame.atom_data[force_idx_offset].symbol;
push_str_line(buf, symbol);
push_component_line(buf, "Forces of Component", type_idx);
for i in 0..num_atoms_in_type {
let atom = &frame.atom_data[force_idx_offset + i];
let [fx, fy, fz] = atom.force().unwrap_or([0.0; 3]);
push_xyz_line(
buf,
ryu,
slots,
[fx, fy, fz],
prec,
written_fixed_flag(atom.fixed(), spec_version),
atom.atom_id,
);
}
force_idx_offset += num_atoms_in_type;
}
}
if has_eng {
buf.push(b'\n');
let mut energy_idx_offset = 0;
for (type_idx, &num_atoms_in_type) in frame.header.natms_per_type.iter().enumerate()
{
let symbol = &frame.atom_data[energy_idx_offset].symbol;
push_str_line(buf, symbol);
push_component_line(buf, "Energies of Component", type_idx);
for i in 0..num_atoms_in_type {
let atom = &frame.atom_data[energy_idx_offset + i];
let e = atom.energy().unwrap_or(0.0);
push_scalar_line(
buf,
ryu,
slots,
e,
prec,
written_fixed_flag(atom.fixed(), spec_version),
atom.atom_id,
);
}
energy_idx_offset += num_atoms_in_type;
}
}
if has_chg {
buf.push(b'\n');
let mut off = 0;
for (type_idx, &num_atoms_in_type) in frame.header.natms_per_type.iter().enumerate()
{
let symbol = &frame.atom_data[off].symbol;
push_str_line(buf, symbol);
push_component_line(buf, "Charges of Component", type_idx);
for i in 0..num_atoms_in_type {
let atom = &frame.atom_data[off + i];
let q = atom.charge().unwrap_or(0.0);
push_scalar_line(
buf,
ryu,
slots,
q,
prec,
written_fixed_flag(atom.fixed(), spec_version),
atom.atom_id,
);
}
off += num_atoms_in_type;
}
}
if has_spn {
buf.push(b'\n');
let mut off = 0;
for (type_idx, &num_atoms_in_type) in frame.header.natms_per_type.iter().enumerate()
{
let symbol = &frame.atom_data[off].symbol;
push_str_line(buf, symbol);
push_component_line(buf, "Spins of Component", type_idx);
for i in 0..num_atoms_in_type {
let atom = &frame.atom_data[off + i];
let s = atom.spin().unwrap_or(0.0);
push_scalar_line(
buf,
ryu,
slots,
s,
prec,
written_fixed_flag(atom.fixed(), spec_version),
atom.atom_id,
);
}
off += num_atoms_in_type;
}
}
if has_mm {
buf.push(b'\n');
let mut off = 0;
for (type_idx, &num_atoms_in_type) in frame.header.natms_per_type.iter().enumerate()
{
let symbol = &frame.atom_data[off].symbol;
push_str_line(buf, symbol);
push_component_line(buf, "Magmoms of Component", type_idx);
for i in 0..num_atoms_in_type {
let atom = &frame.atom_data[off + i];
let [mx, my, mz] = atom.magmom().unwrap_or([0.0; 3]);
push_xyz_line(
buf,
ryu,
slots,
[mx, my, mz],
prec,
written_fixed_flag(atom.fixed(), spec_version),
atom.atom_id,
);
}
off += num_atoms_in_type;
}
}
if has_dsp {
buf.push(b'\n');
let mut off = 0;
for (type_idx, &num_atoms_in_type) in frame.header.natms_per_type.iter().enumerate()
{
let symbol = &frame.atom_data[off].symbol;
push_str_line(buf, symbol);
push_component_line(buf, "Displacements of Component", type_idx);
for i in 0..num_atoms_in_type {
let atom = &frame.atom_data[off + i];
let [dx, dy, dz] = atom.displacement().unwrap_or([0.0; 3]);
push_xyz_line(
buf,
ryu,
slots,
[dx, dy, dz],
prec,
written_fixed_flag(atom.fixed(), spec_version),
atom.atom_id,
);
}
off += num_atoms_in_type;
}
}
if has_spr {
buf.push(b'\n');
let mut off = 0;
for (type_idx, &num_atoms_in_type) in frame.header.natms_per_type.iter().enumerate()
{
let symbol = &frame.atom_data[off].symbol;
push_str_line(buf, symbol);
push_component_line(buf, "Spreads of Component", type_idx);
for i in 0..num_atoms_in_type {
let atom = &frame.atom_data[off + i];
let [sx, sy, sz] = atom.spread().unwrap_or([0.0; 3]);
if !crate::types::spread_row_ok([sx, sy, sz]) {
let err = if sx.is_finite() && sy.is_finite() && sz.is_finite() {
crate::types::invalid_spread(off + i)
} else {
crate::error::ParseError::NonFiniteValue(format!(
"spreads atom {}",
off + i
))
};
return Err(WriteError::Format(err));
}
push_xyz_line(
buf,
ryu,
slots,
[sx, sy, sz],
prec,
written_fixed_flag(atom.fixed(), spec_version),
atom.atom_id,
);
}
off += num_atoms_in_type;
}
}
}
if self.keep_scratch {
if let Some(out) = stream_out.as_mut() {
out.armed = false;
}
return Ok(());
}
if let Some(out) = stream_out.as_mut() {
if !self.scratch.is_empty() {
out.file.write_all(&self.scratch)?;
self.scratch.clear();
}
out.file.flush()?;
out.armed = false;
} else if self.scratch.len() >= 64 * 1024 {
self.writer.flush()?;
self.writer.get_mut().write_all(&self.scratch)?;
} else {
self.writer.write_all(&self.scratch)?;
}
Ok(())
}
pub fn flush(&mut self) -> io::Result<()> {
self.writer.flush()
}
pub fn extend<'a>(
&mut self,
frames: impl Iterator<Item = &'a ConFrame>,
) -> Result<(), WriteError> {
for frame in frames {
self.write_frame(frame)?;
}
Ok(())
}
}
fn ids_are_dense(atoms: &[crate::types::AtomDatum], workers: usize) -> bool {
let n = atoms.len();
if workers < 2 || n < WIDE_ATOM_GRAIN {
return atoms
.iter()
.enumerate()
.all(|(index, atom)| atom.atom_id == index as u64);
}
#[cfg(feature = "parallel")]
{
let chunk = n.div_ceil(workers);
return crate::parser::on_pool(workers, || {
use rayon::prelude::*;
atoms.par_chunks(chunk).enumerate().all(|(index, group)| {
let base = index * chunk;
group
.iter()
.enumerate()
.all(|(offset, atom)| atom.atom_id == (base + offset) as u64)
})
});
}
#[cfg(not(feature = "parallel"))]
{
let _ = workers;
atoms
.iter()
.enumerate()
.all(|(index, atom)| atom.atom_id == index as u64)
}
}
fn ensure_unique_ids(atoms: &[crate::types::AtomDatum], workers: usize) -> Result<(), WriteError> {
let n = atoms.len();
if n == 0 {
return Ok(());
}
if ids_are_dense(atoms, workers) {
return Ok(());
}
let mut seen_ids = std::collections::HashSet::with_capacity(n);
for atom in atoms {
if !seen_ids.insert(atom.atom_id) {
return Err(WriteError::Format(
crate::error::ParseError::ValidationError(format!(
"atom id {} is repeated in the frame",
atom.atom_id
)),
));
}
}
Ok(())
}
fn written_fixed_flag(fixed: [bool; 3], spec_version: u32) -> u8 {
if spec_version < 2 {
if fixed.iter().any(|axis| *axis) { 1 } else { 0 }
} else {
encode_fixed_bitmask(fixed)
}
}
fn push_u64(buf: &mut Vec<u8>, n: u64) {
const DIGITS: &[u8; 200] = b"00010203040506070809101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899";
let mut tmp = [0u8; 20];
let mut i = 20;
let mut n = n;
while n >= 100 {
let q = n / 100;
let r = (n % 100) as usize;
i -= 2;
tmp[i] = DIGITS[r * 2];
tmp[i + 1] = DIGITS[r * 2 + 1];
n = q;
}
if n < 10 {
i -= 1;
tmp[i] = b'0' + n as u8;
} else {
let r = n as usize;
i -= 2;
tmp[i] = DIGITS[r * 2];
tmp[i + 1] = DIGITS[r * 2 + 1];
}
buf.extend_from_slice(&tmp[i..]);
}
fn push_u128(buf: &mut Vec<u8>, mut n: u128) {
if n <= u128::from(u64::MAX) {
push_u64(buf, n as u64);
return;
}
let mut tmp = [0u8; 40];
let mut i = 40;
while n > 0 {
i -= 1;
tmp[i] = b'0' + (n % 10) as u8;
n /= 10;
}
buf.extend_from_slice(&tmp[i..]);
}
fn push_str_line(buf: &mut Vec<u8>, text: &str) {
buf.extend_from_slice(text.as_bytes());
buf.push(b'\n');
}
fn push_component_line(buf: &mut Vec<u8>, label: &str, index: usize) {
buf.extend_from_slice(label.as_bytes());
buf.push(b' ');
push_u64(buf, (index + 1) as u64);
buf.push(b'\n');
}
fn push_f64(buf: &mut Vec<u8>, ryu: &mut ryu::Buffer, value: f64, format: FloatFormat) {
match format {
FloatFormat::DecimalPlaces(precision) => push_f64_prec(buf, value, precision),
FloatFormat::RoundTrip => push_round_trip(buf, ryu, None, value),
}
}
fn push_round_trip(
buf: &mut Vec<u8>,
ryu: &mut ryu::Buffer,
slot: Option<&mut FmtSlot>,
value: f64,
) {
if !value.is_finite() {
let text = ryu.format(value);
buf.extend_from_slice(text.as_bytes());
return;
}
if let Some(slot) = slot {
let bits = value.to_bits();
if slot.bits == bits {
if slot.state == 2 {
buf.extend_from_slice(&slot.text[..slot.len as usize]);
return;
}
if slot.state == 1 {
let n = unsafe { ryu::raw::format64(value, slot.text.as_mut_ptr()) };
debug_assert!(n <= 24);
slot.len = n as u8;
slot.state = 2;
buf.extend_from_slice(&slot.text[..n]);
return;
}
}
write_finite(buf, value);
slot.bits = bits;
slot.state = 1;
return;
}
write_finite(buf, value);
}
fn write_finite(buf: &mut Vec<u8>, value: f64) {
if let Some((bytes, n)) = encode_fixed_decimal_into(value) {
buf.extend_from_slice(&bytes[..n]);
return;
}
if buf.capacity() - buf.len() < 24 {
let mut tmp = [0u8; 24];
let n = unsafe { ryu::raw::format64(value, tmp.as_mut_ptr()) };
debug_assert!(n <= 24);
buf.extend_from_slice(&tmp[..n]);
return;
}
let n = unsafe { ryu::raw::format64(value, buf.as_mut_ptr().add(buf.len())) };
debug_assert!(n <= 24);
unsafe { buf.set_len(buf.len() + n) };
}
fn round_trip_width(value: f64) -> usize {
if !value.is_finite() {
return ryu::Buffer::new().format(value).len();
}
let bits = value.to_bits();
let sign = bits >> 63 != 0;
let mag = f64::from_bits(bits & !(1u64 << 63));
if let Some(n) = fixed_decimal_width(mag, sign) {
return n;
}
if let Some((_, n)) = encode_fixed_decimal_into(value) {
return n;
}
let mut tmp = [0u8; 24];
unsafe { ryu::raw::format64(value, tmp.as_mut_ptr()) }
}
fn fixed_decimal_width(mag: f64, sign: bool) -> Option<usize> {
let sign_len = usize::from(sign);
if mag == 0.0 {
return Some(sign_len + 3);
}
if mag <= 1e16 && mag.fract() == 0.0 {
let n = mag as u64;
if n < 10u64.pow(16) && n as f64 == mag {
return Some(sign_len + count_digits(n) as usize + 2);
}
}
let scaled = mag * 1_000_000.0;
if !scaled.is_finite() || scaled > 1e16 {
return None;
}
let nearest = scaled.round();
let mut n = nearest as u64;
if nearest != n as f64 || (n as f64) / 1_000_000.0 != mag {
return None;
}
let mut places = 6u32;
while places > 0 && n.is_multiple_of(10) {
n /= 10;
places -= 1;
}
if places == 0 {
if n >= 10u64.pow(16) {
return None;
}
return Some(sign_len + count_digits(n) as usize + 2);
}
let digits = count_digits(n);
let kk = digits as i32 - places as i32;
if !(kk > -5 && kk <= 16) {
return None;
}
let body = if kk > 0 {
digits as usize + 1
} else {
2 + (-kk) as usize + digits as usize
};
Some(sign_len + body)
}
fn encode_fixed_decimal_into(value: f64) -> Option<([u8; 32], usize)> {
let mut out = [0u8; 32];
let n = encode_fixed_decimal(value, &mut out)?;
Some((out, n))
}
fn encode_fixed_decimal(value: f64, out: &mut [u8; 32]) -> Option<usize> {
let bits = value.to_bits();
let sign = bits >> 63 != 0;
let mag = f64::from_bits(bits & !(1u64 << 63));
if mag == 0.0 {
let text: &[u8] = if sign { b"-0.0" } else { b"0.0" };
out[..text.len()].copy_from_slice(text);
return Some(text.len());
}
if mag <= 1e16 && mag.fract() == 0.0 {
let n = mag as u64;
if n < 10u64.pow(16) && n as f64 == mag {
return Some(write_fixed_digits(out, sign, n, 0));
}
}
let mut scale = 10.0f64;
for d in 1u32..=6 {
let scaled = mag * scale;
if !scaled.is_finite() || scaled > 1e16 {
break;
}
let nearest = scaled.round();
let n = nearest as u64;
if nearest == n as f64 && (n as f64) / scale == mag {
let mut n = n;
let mut places = d;
while places > 0 && n.is_multiple_of(10) {
n /= 10;
places -= 1;
}
if places == 0 {
if n >= 10u64.pow(16) {
return None;
}
return Some(write_fixed_digits(out, sign, n, 0));
}
let digits = count_digits(n);
let kk = digits as i32 - places as i32;
if !(kk > -5 && kk <= 16) {
return None;
}
return Some(write_fixed_digits(out, sign, n, places));
}
scale *= 10.0;
}
None
}
fn count_digits(mut n: u64) -> u32 {
let mut digits = 1u32;
while n >= 10 {
n /= 10;
digits += 1;
}
digits
}
fn write_fixed_digits(out: &mut [u8; 32], sign: bool, mut n: u64, places: u32) -> usize {
let mut tmp = [0u8; 20];
let mut i = tmp.len();
if n == 0 {
tmp[i - 1] = b'0';
i -= 1;
} else {
while n > 0 {
i -= 1;
tmp[i] = b'0' + (n % 10) as u8;
n /= 10;
}
}
let digs = &tmp[i..];
let mut at = 0usize;
if sign {
out[at] = b'-';
at += 1;
}
if places == 0 {
out[at..at + digs.len()].copy_from_slice(digs);
at += digs.len();
out[at] = b'.';
out[at + 1] = b'0';
return at + 2;
}
let kk = digs.len() as i32 - places as i32;
if kk > 0 {
let cut = kk as usize;
out[at..at + cut].copy_from_slice(&digs[..cut]);
at += cut;
out[at] = b'.';
at += 1;
out[at..at + digs.len() - cut].copy_from_slice(&digs[cut..]);
at += digs.len() - cut;
return at;
}
out[at] = b'0';
out[at + 1] = b'.';
at += 2;
for _ in 0..(-kk) {
out[at] = b'0';
at += 1;
}
out[at..at + digs.len()].copy_from_slice(digs);
at + digs.len()
}
fn push_f64_prec(buf: &mut Vec<u8>, v: f64, prec: usize) {
if !v.is_finite() || prec != 6 {
let _ = write!(buf, "{v:.prec$}");
return;
}
if v.is_sign_negative() {
buf.push(b'-');
}
let ax = v.abs();
if prec == 0 {
push_u128(buf, ax.round() as u128);
return;
}
let scale = 10u128.pow(prec as u32);
let n = (ax * scale as f64).round() as u128;
let int_part = n / scale;
let frac = n % scale;
push_u128(buf, int_part);
buf.push(b'.');
let mut tmp = [b'0'; 20];
let mut x = frac;
let mut i = prec;
while i > 0 {
i -= 1;
tmp[i] = b'0' + (x % 10) as u8;
x /= 10;
}
buf.extend_from_slice(&tmp[..prec]);
}
const WIDE_ATOM_GRAIN: usize = 4096;
fn take_format_buf(cap: usize) -> Vec<u8> {
let mut pool = FORMAT_BUFS.lock().unwrap_or_else(|err| err.into_inner());
if let Some(mut buf) = pool.pop() {
if buf.capacity() < cap {
buf.reserve(cap);
}
buf.clear();
return buf;
}
drop(pool);
Vec::with_capacity(cap)
}
fn give_format_buf(buf: Vec<u8>) {
let mut pool = FORMAT_BUFS.lock().unwrap_or_else(|err| err.into_inner());
if pool.len() < 128 {
pool.push(buf);
}
}
static FORMAT_BUFS: std::sync::Mutex<Vec<Vec<u8>>> = std::sync::Mutex::new(Vec::new());
fn coord_slab(atoms: &[crate::types::AtomDatum]) -> Option<*const f64> {
let n = atoms.len();
if n < 2 {
return None;
}
let base = atoms[0].coord_ptr();
if base.is_null() {
return None;
}
let aligned = |index: usize| unsafe { atoms[index].coord_ptr() == base.add(index * 3) };
if aligned(1) && aligned(n / 2) && aligned(n - 1) {
Some(base)
} else {
None
}
}
fn load_xyz(atom: &crate::types::AtomDatum, slab: Option<*const f64>, index: usize) -> [f64; 3] {
if let Some(base) = slab {
unsafe {
let ptr = base.add(index * 3);
[*ptr, *ptr.add(1), *ptr.add(2)]
}
} else {
[atom.x(), atom.y(), atom.z()]
}
}
fn push_atom_motion(
coord: &mut Vec<u8>,
vel: Option<&mut Vec<u8>>,
ryu: &mut ryu::Buffer,
slots: &mut [FmtSlot; 3],
vel_slots: &mut [FmtSlot; 3],
xyz: [f64; 3],
atom: &crate::types::AtomDatum,
prec: FloatFormat,
spec_version: u32,
) {
let fixed = written_fixed_flag(atom.fixed(), spec_version);
push_xyz_line(coord, ryu, slots, xyz, prec, fixed, atom.atom_id);
if let Some(vel) = vel {
let [vx, vy, vz] = atom.velocity().unwrap_or([0.0; 3]);
push_xyz_line(vel, ryu, vel_slots, [vx, vy, vz], prec, fixed, atom.atom_id);
}
}
pub(crate) fn binding_header_bytes(frame: &ConFrame) -> Result<Vec<u8>, WriteError> {
let spec = frame.header.spec_version.max(crate::CON_SPEC_VERSION);
let meta = binding_metadata(frame);
let mut writer = ConFrameWriter::new(Vec::<u8>::new());
writer.refresh_metadata_cache(
spec,
&meta,
frame.has_velocities(),
frame.has_forces(),
frame.has_energies(),
frame.has_charges(),
frame.has_spins(),
frame.has_magmoms(),
frame.has_displacements(),
frame.has_spreads(),
)?;
let (boxl, angles) = match frame.header.lattice_vectors() {
Some(vecs) => crate::types::lengths_and_angles(vecs),
None => (frame.header.boxl, frame.header.angles),
};
let prec = FloatFormat::RoundTrip;
let mut buf = Vec::with_capacity(256);
push_str_line(&mut buf, &frame.header.prebox_header.user);
if spec < 2 {
push_str_line(&mut buf, frame.header.prebox_header.metadata_line());
} else {
let line = writer
.metadata_cache
.as_ref()
.expect("metadata")
.serialized
.as_bytes();
buf.extend_from_slice(line);
buf.push(b'\n');
}
let ryu = &mut writer.ryu;
push_f64(&mut buf, ryu, boxl[0], prec);
buf.push(b' ');
push_f64(&mut buf, ryu, boxl[1], prec);
buf.push(b' ');
push_f64(&mut buf, ryu, boxl[2], prec);
buf.push(b'\n');
push_f64(&mut buf, ryu, angles[0], prec);
buf.push(b' ');
push_f64(&mut buf, ryu, angles[1], prec);
buf.push(b' ');
push_f64(&mut buf, ryu, angles[2], prec);
buf.push(b'\n');
push_str_line(&mut buf, &frame.header.postbox_header[0]);
push_str_line(&mut buf, &frame.header.postbox_header[1]);
push_u64(&mut buf, frame.header.natm_types as u64);
buf.push(b'\n');
for (i, n) in frame.header.natms_per_type.iter().enumerate() {
if i > 0 {
buf.push(b' ');
}
push_u64(&mut buf, *n as u64);
}
buf.push(b'\n');
for (i, m) in frame.header.masses_per_type.iter().enumerate() {
if i > 0 {
buf.push(b' ');
}
push_f64(&mut buf, ryu, *m, prec);
}
buf.push(b'\n');
Ok(buf)
}
pub(crate) fn format_coord_range(frame: &ConFrame, begin: usize, end: usize) -> Vec<u8> {
let spec = frame.header.spec_version.max(crate::CON_SPEC_VERSION);
let prec = FloatFormat::RoundTrip;
let mut buf = Vec::with_capacity((end - begin).saturating_mul(32) + 64);
let mut ryu = ryu::Buffer::new();
let mut slots = [FmtSlot::new(); 3];
let mut vel_slots = [FmtSlot::new(); 3];
let mut offset = 0usize;
for (type_idx, &count) in frame.header.natms_per_type.iter().enumerate() {
let type_end = offset + count;
let from = begin.max(offset);
let to = end.min(type_end);
if from < to {
if from == offset {
let symbol = frame.atom_data[offset].symbol.as_ref();
push_str_line(&mut buf, symbol);
push_component_line(&mut buf, "Coordinates of Component", type_idx);
}
let atoms = &frame.atom_data[from..to];
let slab = coord_slab(atoms);
for (index, atom) in atoms.iter().enumerate() {
push_atom_motion(
&mut buf,
None,
&mut ryu,
&mut slots,
&mut vel_slots,
load_xyz(atom, slab, index),
atom,
prec,
spec,
);
}
}
offset = type_end;
}
buf
}
pub fn format_coordinate_lines(atoms: &[crate::types::AtomDatum]) -> Result<Vec<u8>, WriteError> {
let mut buf = take_format_buf(atoms.len().saturating_mul(48));
let mut ryu = ryu::Buffer::new();
let mut slots = [FmtSlot::new(); 3];
append_coord_rows(
&mut buf,
&mut ryu,
&mut slots,
atoms,
FloatFormat::RoundTrip,
2,
1,
false,
false,
None,
None,
)?;
Ok(buf)
}
pub fn recycle_format_buffer(buf: Vec<u8>) {
give_format_buf(buf);
}
fn append_velocity_lines(
buf: &mut Vec<u8>,
ryu: &mut ryu::Buffer,
atoms: &[crate::types::AtomDatum],
prec: FloatFormat,
spec_version: u32,
) {
let mut vel_slots = [FmtSlot::new(); 3];
for atom in atoms {
let [vx, vy, vz] = atom.velocity().unwrap_or([0.0; 3]);
push_xyz_line(
buf,
ryu,
&mut vel_slots,
[vx, vy, vz],
prec,
written_fixed_flag(atom.fixed(), spec_version),
atom.atom_id,
);
}
}
fn append_coord_rows(
buf: &mut Vec<u8>,
ryu: &mut ryu::Buffer,
slots: &mut [FmtSlot; 3],
atoms: &[crate::types::AtomDatum],
prec: FloatFormat,
spec_version: u32,
workers: usize,
mapped: bool,
force_map: bool,
stream: Option<&mut File>,
mut vel: Option<&mut Vec<u8>>,
) -> Result<(), WriteError> {
let want_vel = vel.is_some();
if let Some(file) = stream {
if mapped && !atoms.is_empty() {
file.write_all(buf).map_err(WriteError::Io)?;
buf.clear();
let vel_bytes = stream_coord_rows(
file,
atoms,
workers,
prec,
spec_version,
want_vel,
force_map,
)?;
if let Some(dest) = vel.as_mut() {
dest.extend_from_slice(&vel_bytes);
}
return Ok(());
}
}
if workers < 2 || atoms.len() < WIDE_ATOM_GRAIN {
let mut vel_slots = [FmtSlot::new(); 3];
let slab = coord_slab(atoms);
for (index, atom) in atoms.iter().enumerate() {
push_atom_motion(
buf,
vel.as_deref_mut(),
ryu,
slots,
&mut vel_slots,
load_xyz(atom, slab, index),
atom,
prec,
spec_version,
);
}
return Ok(());
}
#[cfg(feature = "parallel")]
{
let groups = workers.saturating_mul(4).max(workers);
let chunk = atoms.len().div_ceil(groups).max(1);
let nparts = atoms.len().div_ceil(chunk);
let bad_ids = std::sync::atomic::AtomicBool::new(false);
let next = std::sync::atomic::AtomicUsize::new(0);
let parts: Vec<std::sync::Mutex<Option<(Vec<u8>, Vec<u8>)>>> =
(0..nparts).map(|_| std::sync::Mutex::new(None)).collect();
crate::parser::on_pool(workers, || {
use rayon::prelude::*;
(0..workers).into_par_iter().for_each(|_| {
let mut ryu = ryu::Buffer::new();
let mut slots = [FmtSlot::new(); 3];
let mut vel_slots = [FmtSlot::new(); 3];
loop {
let index = next.fetch_add(1, std::sync::atomic::Ordering::Relaxed);
if index >= nparts {
break;
}
let start = index * chunk;
let end = (start + chunk).min(atoms.len());
let group = &atoms[start..end];
let mut local = Vec::with_capacity(group.len().saturating_mul(48));
let mut local_vel = if want_vel {
Vec::with_capacity(group.len().saturating_mul(48))
} else {
Vec::new()
};
let slab = coord_slab(group);
for (offset, atom) in group.iter().enumerate() {
if atom.atom_id != (start + offset) as u64 {
bad_ids.store(true, std::sync::atomic::Ordering::Relaxed);
}
push_atom_motion(
&mut local,
if want_vel { Some(&mut local_vel) } else { None },
&mut ryu,
&mut slots,
&mut vel_slots,
load_xyz(atom, slab, offset),
atom,
prec,
spec_version,
);
}
*parts[index].lock().expect("format part") = Some((local, local_vel));
}
});
});
if bad_ids.load(std::sync::atomic::Ordering::Relaxed) {
ensure_unique_ids(atoms, 1)?;
}
let mut coord_len = 0usize;
let mut vel_len = 0usize;
for part in &parts {
let guard = part.lock().expect("format part");
if let Some((coord, rows)) = guard.as_ref() {
coord_len += coord.len();
vel_len += rows.len();
}
}
let coord_at = buf.len();
buf.reserve(coord_len);
if let Some(dest) = vel.as_deref_mut() {
dest.reserve(vel_len);
}
if coord_len > 0 {
let coord_ptr = buf.as_mut_ptr() as usize;
let mut jobs = Vec::with_capacity(nparts);
let mut at = coord_at;
for part in &parts {
let guard = part.lock().expect("format part");
let Some((coord, _)) = guard.as_ref() else {
continue;
};
let n = coord.len();
if n == 0 {
continue;
}
jobs.push((coord.as_ptr() as usize, coord_ptr + at, n));
at += n;
}
let next_copy = std::sync::atomic::AtomicUsize::new(0);
crate::parser::on_pool(workers, || {
use rayon::prelude::*;
(0..workers).into_par_iter().for_each(|_| {
loop {
let index = next_copy.fetch_add(1, std::sync::atomic::Ordering::Relaxed);
if index >= jobs.len() {
break;
}
let (src, dest, n) = jobs[index];
unsafe {
std::ptr::copy_nonoverlapping(src as *const u8, dest as *mut u8, n);
}
}
});
});
unsafe { buf.set_len(coord_at + coord_len) };
}
if let Some(dest) = vel.as_deref_mut() {
for part in &parts {
let guard = part.lock().expect("format part");
if let Some((_, rows)) = guard.as_ref() {
dest.extend_from_slice(rows);
}
}
}
return Ok(());
}
#[cfg(not(feature = "parallel"))]
{
let _ = workers;
let mut vel_slots = [FmtSlot::new(); 3];
let slab = coord_slab(atoms);
for (index, atom) in atoms.iter().enumerate() {
push_atom_motion(
buf,
vel.as_deref_mut(),
ryu,
slots,
&mut vel_slots,
load_xyz(atom, slab, index),
atom,
prec,
spec_version,
);
}
return Ok(());
}
}
fn physical_cores() -> usize {
let logical = std::thread::available_parallelism()
.map(|n| n.get())
.unwrap_or(1);
let siblings =
std::fs::read_to_string("/sys/devices/system/cpu/cpu0/topology/thread_siblings_list")
.ok()
.map(|text| {
text.split([',', '-'])
.filter(|piece| !piece.trim().is_empty())
.count()
})
.filter(|n| *n > 0)
.unwrap_or(1);
(logical / siblings).max(1)
}
fn side_copy_workers(workers: usize) -> bool {
#[cfg(miri)]
{
let _ = workers;
return false;
}
#[cfg(not(miri))]
{
workers >= 2 && workers.saturating_mul(2) <= physical_cores()
}
}
pub(crate) fn advise_shared_store(file: &File, ptr: *mut u8, len: usize) {
#[cfg(all(target_os = "linux", not(miri)))]
if len >= 2 * 1024 * 1024 && mount_requests_huge(file) {
unsafe {
let _ = libc::madvise(ptr.cast::<libc::c_void>(), len, libc::MADV_HUGEPAGE);
}
}
#[cfg(any(not(target_os = "linux"), miri))]
{
let _ = (file, ptr, len);
}
}
#[cfg(all(target_os = "linux", not(miri)))]
fn mount_requests_huge(file: &File) -> bool {
use std::os::fd::AsRawFd;
let Ok(linked) = std::fs::read_link(format!("/proc/self/fd/{}", file.as_raw_fd())) else {
return false;
};
let Ok(table) = std::fs::read_to_string("/proc/mounts") else {
return false;
};
let path = linked.to_string_lossy();
let mut best_len = 0usize;
let mut huge = false;
for line in table.lines() {
let mut fields = line.split_whitespace();
let _device = fields.next();
let Some(raw_mount) = fields.next() else {
continue;
};
let mount = raw_mount.replace("\\040", " ");
if mount.len() >= best_len && path.starts_with(&mount) {
let _fstype = fields.next();
let options = fields.next().unwrap_or("");
best_len = mount.len();
huge = options.split(',').any(|opt| {
matches!(
opt,
"huge=within_size" | "huge=always" | "huge=advise" | "huge=force"
)
});
}
}
huge
}
pub(crate) fn file_is_tmpfs(file: &File) -> bool {
#[cfg(all(not(miri), target_os = "linux"))]
{
use std::os::fd::AsRawFd;
let mut stat = unsafe { std::mem::zeroed::<libc::statfs>() };
let rc = unsafe { libc::fstatfs(file.as_raw_fd(), &mut stat) };
return rc == 0 && stat.f_type == 0x01021994;
}
#[cfg(any(miri, not(target_os = "linux")))]
{
let _ = file;
false
}
}
fn stream_coord_pipeline(
file: &mut File,
atoms: &[crate::types::AtomDatum],
workers: usize,
prec: FloatFormat,
spec_version: u32,
want_vel: bool,
) -> Result<Vec<u8>, WriteError> {
let workers = workers.max(2).min(atoms.len().max(1));
let formatters = workers - 1;
let groups = if file_is_tmpfs(file) {
workers.saturating_mul(32)
} else {
workers.saturating_mul(8)
}
.max(workers);
let chunk = atoms.len().div_ceil(groups).max(1);
let nparts = atoms.len().div_ceil(chunk);
let base = file.stream_position().map_err(WriteError::Io)?;
let out = file.try_clone().map_err(WriteError::Io)?;
let next = std::sync::atomic::AtomicUsize::new(0);
let bad_ids = std::sync::atomic::AtomicBool::new(false);
let slots: Vec<std::sync::Mutex<Option<(Vec<u8>, Vec<u8>)>>> =
(0..nparts).map(|_| std::sync::Mutex::new(None)).collect();
let flag = &bad_ids;
let (end, mut vel_bytes) = std::thread::scope(|scope| -> Result<_, WriteError> {
let mut handles = Vec::with_capacity(formatters);
for worker in 0..formatters {
let next = &next;
let slots = &slots;
handles.push(scope.spawn(move || {
crate::parser::pin_worker(worker + 1);
let mut ryu = ryu::Buffer::new();
let mut slots_fmt = [FmtSlot::new(); 3];
let mut vel_slots = [FmtSlot::new(); 3];
loop {
let index = next.fetch_add(1, std::sync::atomic::Ordering::Relaxed);
if index >= nparts {
break;
}
let start = index * chunk;
let end = (start + chunk).min(atoms.len());
let group = &atoms[start..end];
let mut coord = Vec::with_capacity(group.len().saturating_mul(48));
let mut vel = if want_vel {
Vec::with_capacity(group.len().saturating_mul(48))
} else {
Vec::new()
};
let slab = coord_slab(group);
for (offset, atom) in group.iter().enumerate() {
if atom.atom_id != (start + offset) as u64 {
flag.store(true, std::sync::atomic::Ordering::Relaxed);
}
push_atom_motion(
&mut coord,
if want_vel { Some(&mut vel) } else { None },
&mut ryu,
&mut slots_fmt,
&mut vel_slots,
load_xyz(atom, slab, offset),
atom,
prec,
spec_version,
);
}
*slots[index].lock().expect("format slot") = Some((coord, vel));
}
}));
}
let mut offset = base;
let mut vel_bytes = Vec::new();
#[cfg(all(not(miri), target_os = "linux"))]
let mut old_affinity = unsafe { std::mem::zeroed::<libc::cpu_set_t>() };
#[cfg(all(not(miri), target_os = "linux"))]
unsafe {
libc::sched_getaffinity(0, std::mem::size_of::<libc::cpu_set_t>(), &mut old_affinity);
}
crate::parser::pin_worker(0);
let written = (|| -> Result<_, WriteError> {
for index in 0..nparts {
let (coord, vel) = loop {
if let Some(part) = slots[index].lock().expect("format slot").take() {
break part;
}
std::hint::spin_loop();
};
pwrite_full(&out, &coord, offset).map_err(WriteError::Io)?;
offset = offset.checked_add(coord.len() as u64).ok_or_else(|| {
WriteError::Io(io::Error::new(
io::ErrorKind::InvalidData,
"coordinate block is too large",
))
})?;
if want_vel {
vel_bytes.extend_from_slice(&vel);
}
}
Ok(())
})();
#[cfg(all(not(miri), target_os = "linux"))]
unsafe {
libc::sched_setaffinity(0, std::mem::size_of::<libc::cpu_set_t>(), &old_affinity);
}
written?;
for handle in handles {
if handle.join().is_err() {
return Err(WriteError::Io(io::Error::other("formatter panicked")));
}
}
Ok((offset, vel_bytes))
})?;
file.set_len(end).map_err(WriteError::Io)?;
file.seek(SeekFrom::Start(end)).map_err(WriteError::Io)?;
if bad_ids.load(std::sync::atomic::Ordering::Relaxed) {
ensure_unique_ids(atoms, 1)?;
}
if !want_vel {
vel_bytes.clear();
}
Ok(vel_bytes)
}
fn mapped_line_budget(prec: FloatFormat) -> usize {
match prec {
FloatFormat::RoundTrip => 128,
FloatFormat::DecimalPlaces(places) => {
let float = 1usize
.saturating_add(309)
.saturating_add(1)
.saturating_add(places);
float.saturating_mul(3).saturating_add(4 + 1 + 1 + 20 + 1)
}
}
}
fn stream_coord_rows(
file: &mut File,
atoms: &[crate::types::AtomDatum],
workers: usize,
prec: FloatFormat,
spec_version: u32,
want_vel: bool,
force_map: bool,
) -> Result<Vec<u8>, WriteError> {
let workers = workers.max(1).min(atoms.len().max(1));
let side_copies = side_copy_workers(workers);
if !side_copies && !force_map {
return stream_coord_pipeline(file, atoms, workers, prec, spec_version, want_vel);
}
let groups = if side_copies {
workers.saturating_mul(8).max(workers)
} else {
workers.saturating_mul(4).max(workers)
};
let chunk = atoms.len().div_ceil(groups).max(1);
let nparts = atoms.len().div_ceil(chunk);
let base = file.stream_position().map_err(WriteError::Io)?;
let line_budget = mapped_line_budget(prec);
let upper = atoms.len().checked_mul(line_budget).ok_or_else(|| {
WriteError::Io(io::Error::new(
io::ErrorKind::InvalidData,
"coordinate block is too large",
))
})?;
let reserve_end = base.checked_add(upper as u64).ok_or_else(|| {
WriteError::Io(io::Error::new(
io::ErrorKind::InvalidData,
"coordinate block is too large",
))
})?;
let mut map = {
file.set_len(reserve_end).map_err(WriteError::Io)?;
Some(unsafe { memmap2::MmapMut::map_mut(&*file).map_err(WriteError::Io)? })
};
let dest_base = map
.as_mut()
.map(|mapping| mapping.as_mut_ptr() as usize + base as usize)
.unwrap_or(0);
let lengths: Vec<std::sync::atomic::AtomicU64> = (0..nparts)
.map(|_| std::sync::atomic::AtomicU64::new(0))
.collect();
let bad_ids = std::sync::atomic::AtomicBool::new(false);
let failed = std::sync::atomic::AtomicBool::new(false);
let next_job = std::sync::atomic::AtomicUsize::new(0);
let produced = std::sync::atomic::AtomicUsize::new(0);
let copies: std::sync::Mutex<Vec<(Vec<u8>, usize)>> = std::sync::Mutex::new(Vec::new());
let ready = std::sync::Condvar::new();
let vel_parts: Vec<std::sync::Mutex<Option<Vec<u8>>>> =
(0..nparts).map(|_| std::sync::Mutex::new(None)).collect();
let flag = &bad_ids;
std::thread::scope(|scope| -> Result<(), WriteError> {
let mut copiers = Vec::with_capacity(workers);
if side_copies {
for copier in 0..workers {
let copies = &copies;
let ready = &ready;
let produced = &produced;
copiers.push(scope.spawn(move || {
crate::parser::pin_worker(workers + copier);
loop {
let job = {
let mut queue = copies.lock().expect("copies");
loop {
if let Some(job) = queue.pop() {
break Some(job);
}
if produced.load(std::sync::atomic::Ordering::Acquire) == nparts {
break None;
}
queue = ready.wait(queue).expect("copies");
}
};
let Some((coord, offset)) = job else {
break;
};
let n = coord.len();
if n > 0 {
unsafe {
std::ptr::copy_nonoverlapping(
coord.as_ptr(),
(dest_base + offset) as *mut u8,
n,
);
}
}
give_format_buf(coord);
}
}));
}
}
let mut formatters = Vec::with_capacity(workers);
for worker in 0..workers {
let lengths = &lengths;
let next_job = &next_job;
let vel_parts = &vel_parts;
let copies = &copies;
let ready = &ready;
let produced = &produced;
let failed = &failed;
formatters.push(scope.spawn(move || -> Result<(), WriteError> {
crate::parser::pin_worker(worker);
let mut ryu = ryu::Buffer::new();
let mut slots = [FmtSlot::new(); 3];
let mut vel_slots = [FmtSlot::new(); 3];
loop {
let index = next_job.fetch_add(1, std::sync::atomic::Ordering::Relaxed);
if index >= nparts {
break;
}
let start = index * chunk;
let end = (start + chunk).min(atoms.len());
let group = &atoms[start..end];
let mut coord = if side_copies {
take_format_buf(group.len().saturating_mul(48))
} else {
Vec::with_capacity(group.len().saturating_mul(48))
};
let mut vel = if want_vel {
Vec::with_capacity(group.len().saturating_mul(48))
} else {
Vec::new()
};
let slab = coord_slab(group);
for (offset, atom) in group.iter().enumerate() {
if atom.atom_id != (start + offset) as u64 {
flag.store(true, std::sync::atomic::Ordering::Relaxed);
}
push_atom_motion(
&mut coord,
if want_vel { Some(&mut vel) } else { None },
&mut ryu,
&mut slots,
&mut vel_slots,
load_xyz(atom, slab, offset),
atom,
prec,
spec_version,
);
}
let n = coord.len();
if n > group.len().saturating_mul(line_budget) {
failed.store(true, std::sync::atomic::Ordering::Release);
lengths[index].store(1, std::sync::atomic::Ordering::Release);
return Err(WriteError::Io(io::Error::other(
"coordinate line exceeded the mapped reserve",
)));
}
lengths[index].store(n as u64 + 1, std::sync::atomic::Ordering::Release);
let mut offset = 0usize;
for prev in 0..index {
loop {
if failed.load(std::sync::atomic::Ordering::Acquire) {
return Err(WriteError::Io(io::Error::other(
"coordinate line exceeded the mapped reserve",
)));
}
let known = lengths[prev].load(std::sync::atomic::Ordering::Acquire);
if known != 0 {
if failed.load(std::sync::atomic::Ordering::Acquire) {
return Err(WriteError::Io(io::Error::other(
"coordinate line exceeded the mapped reserve",
)));
}
offset += (known - 1) as usize;
break;
}
std::hint::spin_loop();
}
}
if side_copies {
let mut queue = copies.lock().expect("copies");
queue.push((coord, offset));
} else {
let n = coord.len();
if n > 0 {
unsafe {
std::ptr::copy_nonoverlapping(
coord.as_ptr(),
(dest_base + offset) as *mut u8,
n,
);
}
}
}
if side_copies {
produced.fetch_add(1, std::sync::atomic::Ordering::Release);
ready.notify_one();
}
if want_vel {
*vel_parts[index].lock().expect("vel part") = Some(vel);
}
}
Ok(())
}));
}
let formatted = (|| -> Result<(), WriteError> {
for handle in formatters {
handle.join().unwrap_or_else(|_| {
Err(WriteError::Io(io::Error::other("formatter panicked")))
})?;
}
Ok(())
})();
produced.store(nparts, std::sync::atomic::Ordering::Release);
ready.notify_all();
for handle in copiers {
handle
.join()
.map_err(|_| WriteError::Io(io::Error::other("copy panicked")))?;
}
formatted?;
Ok(())
})?;
let mut cursor = 0u64;
for length in &lengths {
cursor += length
.load(std::sync::atomic::Ordering::Relaxed)
.saturating_sub(1);
}
drop(map);
let end = base + cursor;
file.set_len(end).map_err(WriteError::Io)?;
file.seek(SeekFrom::Start(end)).map_err(WriteError::Io)?;
if bad_ids.load(std::sync::atomic::Ordering::Relaxed) {
ensure_unique_ids(atoms, 1)?;
}
let mut vel_bytes = Vec::new();
if want_vel {
for part in &vel_parts {
if let Some(vel) = part.lock().expect("vel part").as_ref() {
vel_bytes.extend_from_slice(vel);
}
}
}
drop(vel_parts);
Ok(vel_bytes)
}
fn push_xyz_line(
buf: &mut Vec<u8>,
ryu: &mut ryu::Buffer,
slots: &mut [FmtSlot; 3],
xyz: [f64; 3],
prec: FloatFormat,
fixed: u8,
atom_id: u64,
) {
let [x, y, z] = xyz;
match prec {
FloatFormat::RoundTrip => {
push_round_trip(buf, ryu, Some(&mut slots[0]), x);
buf.push(b' ');
push_round_trip(buf, ryu, Some(&mut slots[1]), y);
buf.push(b' ');
push_round_trip(buf, ryu, Some(&mut slots[2]), z);
}
FloatFormat::DecimalPlaces(_) => {
push_f64(buf, ryu, x, prec);
buf.push(b' ');
push_f64(buf, ryu, y, prec);
buf.push(b' ');
push_f64(buf, ryu, z, prec);
}
}
buf.push(b' ');
buf.push(b'0' + fixed);
buf.push(b' ');
push_u64(buf, atom_id);
buf.push(b'\n');
}
fn push_scalar_line(
buf: &mut Vec<u8>,
ryu: &mut ryu::Buffer,
slots: &mut [FmtSlot; 3],
v: f64,
prec: FloatFormat,
fixed: u8,
atom_id: u64,
) {
match prec {
FloatFormat::RoundTrip => push_round_trip(buf, ryu, Some(&mut slots[0]), v),
FloatFormat::DecimalPlaces(_) => push_f64(buf, ryu, v, prec),
}
buf.push(b' ');
push_u64(buf, u64::from(fixed));
buf.push(b' ');
push_u64(buf, atom_id);
buf.push(b'\n');
}
pub(crate) fn stamp_legacy_velocity_time(
meta: &mut BTreeMap<String, Value>,
spec_version: u32,
has_velocity: bool,
) {
if !has_velocity {
return;
}
let minor = meta
.get(meta::CON_SPEC_MINOR)
.and_then(Value::as_u64)
.unwrap_or(0);
if spec_version >= 3 && minor >= 1 {
return;
}
let mut units = match meta.get(meta::UNITS).cloned() {
Some(Value::Object(map)) => Value::Object(map),
Some(_) => return,
None => crate::units::default_v3_units_json(),
};
if let Some(object) = units.as_object_mut() {
object.insert(
"time".to_string(),
Value::String(crate::units::NATURAL_TIME_UNIT.to_string()),
);
meta.insert(meta::UNITS.to_string(), units);
if meta.contains_key(meta::CON_SPEC_MINOR) {
meta.insert(meta::CON_SPEC_MINOR.to_string(), Value::from(1u64));
}
}
}
pub fn binding_metadata(frame: &ConFrame) -> BTreeMap<String, Value> {
let mut meta = frame.header.metadata.to_map();
meta.remove(meta::CON_SPEC_VERSION);
meta.remove(meta::SECTIONS);
let has_velocity = frame.has_velocities();
stamp_legacy_velocity_time(&mut meta, frame.header.spec_version, has_velocity);
if !meta.contains_key(meta::UNITS) {
meta.insert(
meta::UNITS.to_string(),
crate::units::default_v3_units_json(),
);
}
if !meta.contains_key(meta::CON_SPEC_MINOR) {
meta.insert(meta::CON_SPEC_MINOR.to_string(), Value::from(1u64));
}
meta
}
pub fn binding_rewrite_frame(frame: &ConFrame) -> Result<ConFrame, crate::error::ParseError> {
let mut out = frame.clone();
out.header.metadata = crate::types::FrameMeta::from_map(binding_metadata(frame));
if out.header.spec_version < crate::CON_SPEC_VERSION {
out.header.spec_version = crate::CON_SPEC_VERSION;
}
Ok(out)
}
impl ConFrameWriter<File> {
pub fn from_path<P: AsRef<Path>>(path: P) -> io::Result<Self> {
let file = File::create(path)?;
Ok(Self::new(file))
}
pub fn from_path_with_precision<P: AsRef<Path>>(path: P, precision: usize) -> io::Result<Self> {
let file = File::create(path)?;
Ok(Self::with_precision(file, precision))
}
pub fn from_path_append<P: AsRef<Path>>(path: P) -> io::Result<Self> {
Ok(Self::new(open_uncompressed_append(path.as_ref())?))
}
pub fn from_path_append_with_precision<P: AsRef<Path>>(
path: P,
precision: usize,
) -> io::Result<Self> {
Ok(Self::with_precision(
open_uncompressed_append(path.as_ref())?,
precision,
))
}
}
pub(crate) fn open_uncompressed_append(path: &Path) -> io::Result<File> {
if let Some(ext) = path.extension().and_then(|ext| ext.to_str())
&& (ext == "gz" || ext == "zst")
{
return Err(io::Error::new(
io::ErrorKind::InvalidInput,
"append is only supported for an uncompressed file",
));
}
let mut file = OpenOptions::new()
.create(true)
.read(true)
.append(true)
.open(path)?;
let len = file.seek(SeekFrom::End(0))?;
if len > 0 {
file.seek(SeekFrom::End(-1))?;
let mut last = [0u8; 1];
file.read_exact(&mut last)?;
if last[0] != b'\n' {
file.write_all(b"\n")?;
}
}
Ok(file)
}
pub fn write_frames_to_vec(frames: &[ConFrame]) -> Result<Vec<u8>, WriteError> {
write_frames_to_vec_with_threads(frames, 1, None, false)
}
pub fn write_frames_to_vec_with_threads(
frames: &[ConFrame],
workers: usize,
precision: Option<usize>,
canonical: bool,
) -> Result<Vec<u8>, WriteError> {
let refs: Vec<&ConFrame> = frames.iter().collect();
write_frame_refs_with_threads(&refs, workers, precision, canonical)
}
fn temporary_sibling(path: &Path) -> io::Result<(std::path::PathBuf, File)> {
static NEXT: std::sync::atomic::AtomicU64 = std::sync::atomic::AtomicU64::new(0);
let dir = match path.parent() {
Some(dir) if !dir.as_os_str().is_empty() => dir.to_path_buf(),
_ => std::path::PathBuf::from("."),
};
let tmp = dir.join(format!(
".readcon-write-{}-{}.tmp",
std::process::id(),
NEXT.fetch_add(1, std::sync::atomic::Ordering::Relaxed)
));
let file = OpenOptions::new()
.read(true)
.write(true)
.create(true)
.truncate(true)
.open(&tmp)?;
Ok((tmp, file))
}
fn commit_temporary(
tmp: &Path,
dest: &Path,
result: Result<(), WriteError>,
) -> Result<(), WriteError> {
match result {
Ok(()) => match std::fs::rename(tmp, dest) {
Ok(()) => Ok(()),
Err(err) => {
let _ = std::fs::remove_file(tmp);
Err(WriteError::Io(err))
}
},
Err(err) => {
let _ = std::fs::remove_file(tmp);
Err(err)
}
}
}
pub fn write_frame_refs_to_path(
path: &std::path::Path,
frames: &[&ConFrame],
workers: usize,
precision: Option<usize>,
canonical: bool,
bind: bool,
) -> Result<(), WriteError> {
if workers < 2 {
return write_serial_frames(path, frames, precision, canonical, bind);
}
if workers >= 2 && frames.len() >= 2 {
return pipeline_frame_groups(path, frames, workers, precision, canonical, bind);
}
if workers >= 2 && frames.len() == 1 && frames[0].atom_data.len() >= WIDE_ATOM_GRAIN {
let mut writer = match precision {
Some(precision) => ConFrameWriter::with_precision(Vec::<u8>::new(), precision),
None => ConFrameWriter::new(Vec::<u8>::new()),
}
.canonical(canonical)
.parallel_atoms(workers);
let (tmp, file) = temporary_sibling(path).map_err(WriteError::Io)?;
writer.stream_file = Some(file);
let result = if bind {
writer.write_frame_bound(frames[0])
} else {
writer.write_frame(frames[0])
};
drop(writer);
return commit_temporary(&tmp, path, result);
}
let parts = formatted_parts(frames, workers, precision, canonical, bind)?;
if parts.len() == 1 {
std::fs::write(path, &parts[0]).map_err(WriteError::Io)
} else {
pwrite_parts(path, &parts).map_err(WriteError::Io)
}
}
pub fn write_frame_refs_with_threads(
frames: &[&ConFrame],
workers: usize,
precision: Option<usize>,
canonical: bool,
) -> Result<Vec<u8>, WriteError> {
let mut parts = formatted_parts(frames, workers, precision, canonical, false)?;
if parts.len() == 1 {
return Ok(parts.pop().unwrap_or_default());
}
let mut out = Vec::new();
for part in parts {
out.extend_from_slice(&part);
}
Ok(out)
}
pub(crate) fn frame_byte_weight(frame: &ConFrame) -> usize {
let atoms = frame.atom_data.len().max(1);
if !frame.header.sections.is_empty() {
return atoms.saturating_mul(1 + frame.header.sections.len());
}
let mut blocks = 1usize;
if frame.velocities.is_present() {
blocks += 1;
}
if frame.forces.is_present() {
blocks += 1;
}
if frame.atom_energies.is_present() {
blocks += 1;
}
if frame.charges.is_present() {
blocks += 1;
}
if frame.spins.is_present() {
blocks += 1;
}
if frame.magmoms.is_present() {
blocks += 1;
}
if frame.displacements.is_present() {
blocks += 1;
}
if frame.spreads.is_present() {
blocks += 1;
}
if blocks == 1
&& let Some(atom) = frame.atom_data.first()
{
if atom.has_velocity() {
blocks += 1;
}
if atom.has_forces() {
blocks += 1;
}
if atom.has_energy() {
blocks += 1;
}
if atom.has_charge() {
blocks += 1;
}
if atom.has_spin() {
blocks += 1;
}
if atom.has_magmom() {
blocks += 1;
}
if atom.has_displacement() {
blocks += 1;
}
if atom.has_spread() {
blocks += 1;
}
}
atoms.saturating_mul(blocks)
}
fn equal_frame_ranges(n: usize, groups: usize) -> Vec<(usize, usize)> {
let groups = groups.max(1).min(n.max(1));
let mut out = Vec::with_capacity(groups);
if n == 0 {
return out;
}
for group in 0..groups {
let start = group * n / groups;
let end = (group + 1) * n / groups;
if start < end {
out.push((start, end));
}
}
out
}
pub(crate) fn group_frame_ranges(frames: &[&ConFrame], groups: usize) -> Vec<(usize, usize)> {
let n = frames.len();
if n == 0 {
return Vec::new();
}
let groups = groups.max(1).min(n);
if n == 1 || !sampled_weights_vary(frames) {
return equal_frame_ranges(n, groups);
}
let weights: Vec<usize> = frames
.iter()
.map(|frame| frame_byte_weight(frame))
.collect();
pack_frames_by_weight(&weights, groups)
}
fn sampled_weights_vary(frames: &[&ConFrame]) -> bool {
let anchor = frame_byte_weight(frames[0]);
let last = frames.len() - 1;
let samples = 16.min(frames.len());
for i in 0..samples {
let index = if samples == 1 {
0
} else {
i * last / (samples - 1)
};
if frame_byte_weight(frames[index]) != anchor {
return true;
}
}
false
}
pub(crate) fn pack_frames_by_weight(weights: &[usize], workers: usize) -> Vec<(usize, usize)> {
let n = weights.len();
if n == 0 {
return Vec::new();
}
let workers = workers.max(1).min(n);
let mut prefix = vec![0u64; n + 1];
for (i, weight) in weights.iter().enumerate() {
prefix[i + 1] = prefix[i] + (*weight as u64).max(1);
}
let total = prefix[n];
let mut out = Vec::with_capacity(workers);
let mut start = 0usize;
for worker in 0..workers {
let end = if worker + 1 == workers {
n
} else {
let target = ((worker as u64 + 1) * total) / workers as u64;
let latest = n - (workers - worker - 1);
let mut end = start + 1;
while end < latest && prefix[end] < target {
end += 1;
}
end
};
out.push((start, end));
start = end;
}
out
}
pub(crate) fn pipeline_format_threads(workers: usize) -> usize {
workers.max(1)
}
pub(crate) fn direct_slice_write(workers: usize, total_weight: usize) -> bool {
workers >= 8 && total_weight >= 3_000_000
}
fn written_box(frame: &ConFrame) -> ([f64; 3], [f64; 3]) {
match frame.header.lattice_vectors() {
Some(vecs) => crate::types::lengths_and_angles(vecs),
None => (frame.header.boxl, frame.header.angles),
}
}
fn add_widths(sum: &mut usize, values: &[f64]) {
*sum += sum_spelling_widths(values);
}
fn milli_width_table() -> &'static [u8] {
const N: usize = 1_000_001;
static TABLE: std::sync::OnceLock<Box<[u8]>> = std::sync::OnceLock::new();
TABLE.get_or_init(|| {
let mut table = vec![0u8; N].into_boxed_slice();
for n in 0..N {
table[n] = spell_width((n as f64) / 10_000.0) as u8;
}
table
})
}
fn spell_width(value: f64) -> usize {
if !value.is_finite() {
return ryu::Buffer::new().format(value).len();
}
if let Some((_, n)) = encode_fixed_decimal_into(value) {
return n;
}
let mut tmp = [0u8; 24];
unsafe { ryu::raw::format64(value, tmp.as_mut_ptr()) }
}
fn sum_spelling_widths(values: &[f64]) -> usize {
let table = milli_width_table();
let mut sum = 0usize;
for &value in values {
let bits = value.to_bits();
let neg = bits >> 63 != 0;
let mag = f64::from_bits(bits & !(1u64 << 63));
let scaled = mag * 10_000.0;
let nearest = scaled.round();
if nearest <= 1_000_000.0 {
let n = nearest as usize;
if ((n as f64) / 10_000.0).to_bits() == mag.to_bits() {
sum += table[n] as usize + usize::from(neg);
continue;
}
}
sum += spell_width(value);
}
sum
}
fn packed_float_widths(frame: &ConFrame) -> usize {
let mut sum = 0usize;
let (boxl, angles) = written_box(frame);
add_widths(&mut sum, &boxl);
add_widths(&mut sum, &angles);
add_widths(&mut sum, &frame.header.masses_per_type);
if let Some(base) = coord_slab(&frame.atom_data) {
let n = frame.atom_data.len() * 3;
let coords = unsafe { std::slice::from_raw_parts(base, n) };
add_widths(&mut sum, coords);
} else {
for atom in &frame.atom_data {
sum += round_trip_width(atom.x());
sum += round_trip_width(atom.y());
sum += round_trip_width(atom.z());
}
}
sum += optional_rows_width(frame.velocities.is_present().then(|| &*frame.velocities));
sum += optional_rows_width(frame.forces.is_present().then(|| &*frame.forces));
sum += optional_rows_width(frame.magmoms.is_present().then(|| &*frame.magmoms));
sum += optional_rows_width(
frame
.displacements
.is_present()
.then(|| &*frame.displacements),
);
sum += optional_rows_width(frame.spreads.is_present().then(|| &*frame.spreads));
sum += optional_column_width(
frame
.atom_energies
.is_present()
.then(|| &*frame.atom_energies),
);
sum += optional_column_width(frame.charges.is_present().then(|| &*frame.charges));
sum += optional_column_width(frame.spins.is_present().then(|| &*frame.spins));
if !frame.velocities.is_present() && frame.has_velocities() {
for atom in &frame.atom_data {
let [x, y, z] = atom.velocity().unwrap_or([0.0; 3]);
sum += round_trip_width(x) + round_trip_width(y) + round_trip_width(z);
}
}
if !frame.forces.is_present() && frame.has_forces() {
for atom in &frame.atom_data {
let [x, y, z] = atom.force().unwrap_or([0.0; 3]);
sum += round_trip_width(x) + round_trip_width(y) + round_trip_width(z);
}
}
sum
}
fn optional_rows_width(storage: Option<&crate::storage_dtype::Array2Storage>) -> usize {
let Some(storage) = storage else {
return 0;
};
let Some(values) = storage.f64_slice() else {
return 0;
};
let mut sum = 0usize;
add_widths(&mut sum, values);
sum
}
fn optional_column_width(storage: Option<&crate::storage_dtype::Array1Storage>) -> usize {
let Some(storage) = storage else {
return 0;
};
let Some(values) = storage.f64_slice() else {
return 0;
};
let mut sum = 0usize;
add_widths(&mut sum, values);
sum
}
fn same_packed_shape(reference: &ConFrame, frame: &ConFrame) -> bool {
if frame.atom_data.len() != reference.atom_data.len()
|| frame.atom_data.is_empty()
|| frame.header.spec_version != reference.header.spec_version
|| frame.header.natm_types != reference.header.natm_types
|| frame.header.natms_per_type != reference.header.natms_per_type
|| frame.header.masses_per_type.len() != reference.header.masses_per_type.len()
|| frame.header.prebox_header.user != reference.header.prebox_header.user
|| frame.header.prebox_header.metadata_line()
!= reference.header.prebox_header.metadata_line()
|| frame.header.postbox_header != reference.header.postbox_header
|| frame.has_velocities() != reference.has_velocities()
|| frame.has_forces() != reference.has_forces()
|| frame.has_energies() != reference.has_energies()
|| frame.has_charges() != reference.has_charges()
|| frame.has_spins() != reference.has_spins()
|| frame.has_magmoms() != reference.has_magmoms()
|| frame.has_displacements() != reference.has_displacements()
|| frame.has_spreads() != reference.has_spreads()
|| !ids_are_dense(&frame.atom_data, 1)
{
return false;
}
let mut offset = 0usize;
for &count in &frame.header.natms_per_type {
if frame.atom_data[offset].symbol != reference.atom_data[offset].symbol {
return false;
}
offset += count;
}
true
}
fn packed_frame_lengths(
frames: &[&ConFrame],
workers: usize,
precision: Option<usize>,
canonical: bool,
bind: bool,
) -> Option<Vec<usize>> {
if precision.is_some() || frames.len() < 2 {
return None;
}
let reference = frames[0];
if reference.atom_data.is_empty() || !ids_are_dense(&reference.atom_data, 1) {
return None;
}
let workers = workers.max(1).min(frames.len());
let shape_ok = std::thread::scope(|scope| {
let mut guards = Vec::with_capacity(workers);
for thread in 0..workers {
let start = thread * frames.len() / workers;
let end = (thread + 1) * frames.len() / workers;
guards.push(scope.spawn(move || {
frames[start..end]
.iter()
.all(|frame| same_packed_shape(reference, frame))
}));
}
guards
.into_iter()
.all(|guard| guard.join().unwrap_or(false))
});
if !shape_ok {
return None;
}
let widths = std::thread::scope(|scope| {
let mut guards = Vec::with_capacity(workers);
for thread in 0..workers {
let start = thread * frames.len() / workers;
let end = (thread + 1) * frames.len() / workers;
guards.push(scope.spawn(move || {
#[cfg(feature = "parallel")]
crate::parser::pin_worker(thread);
frames[start..end]
.iter()
.map(|frame| packed_float_widths(frame))
.collect::<Vec<_>>()
}));
}
let mut widths = Vec::with_capacity(frames.len());
for guard in guards {
widths.extend(guard.join().unwrap_or_default());
}
widths
});
if widths.len() != frames.len() {
return None;
}
let sample = format_group_bytes(
std::slice::from_ref(&reference),
1,
precision,
canonical,
bind,
)
.ok()?;
let overhead = sample.len().checked_sub(widths[0])?;
let lengths: Vec<usize> = widths.iter().map(|width| overhead + width).collect();
for &index in [0, frames.len() / 2, frames.len() - 1].iter() {
let got = format_group_bytes(
std::slice::from_ref(&frames[index]),
1,
precision,
canonical,
bind,
)
.ok()?
.len();
if got != lengths[index] {
return None;
}
}
Some(lengths)
}
fn format_group_into(
group: &[&ConFrame],
dest: &mut [u8],
precision: Option<usize>,
canonical: bool,
bind: bool,
) -> Result<(), WriteError> {
if dest.is_empty() {
return Ok(());
}
struct Held {
writer: ConFrameWriter<Vec<u8>>,
}
impl Drop for Held {
fn drop(&mut self) {
std::mem::forget(std::mem::take(&mut self.writer.scratch));
}
}
let cap = dest.len();
let ptr = dest.as_mut_ptr();
let mut held = Held {
writer: match precision {
Some(precision) => ConFrameWriter::with_precision(Vec::new(), precision),
None => ConFrameWriter::new(Vec::new()),
}
.canonical(canonical)
.retain_scratch(true),
};
held.writer.exact_scratch = true;
let previous = std::mem::replace(&mut held.writer.scratch, unsafe {
Vec::from_raw_parts(ptr, 0, cap)
});
drop(previous);
let result = write_group(&mut held.writer, group, bind);
let written = held.writer.scratch.len();
let stable = held.writer.scratch.as_ptr() == ptr;
result?;
if written != cap || !stable {
return Err(WriteError::Io(io::Error::other(
"packed frame length drifted",
)));
}
Ok(())
}
fn place_packed(
file: &File,
frames: &[&ConFrame],
lengths: &[usize],
workers: usize,
precision: Option<usize>,
canonical: bool,
bind: bool,
) -> Result<(), WriteError> {
let total: u64 = lengths.iter().map(|len| *len as u64).sum();
file.set_len(total).map_err(WriteError::Io)?;
if total == 0 {
return Ok(());
}
let ranges = pack_frames_by_weight(lengths, workers);
let mut prefix = vec![0usize; lengths.len() + 1];
for (index, len) in lengths.iter().enumerate() {
prefix[index + 1] = prefix[index] + len;
}
if cfg!(miri) {
let mut parts = Vec::with_capacity(ranges.len());
for &(start, end) in &ranges {
parts.push(format_group_bytes(
&frames[start..end],
1,
precision,
canonical,
bind,
)?);
}
return store_disjoint_parts(file, &parts).map_err(WriteError::Io);
}
let result = std::thread::scope(|scope| {
let mut guards = Vec::with_capacity(ranges.len());
for (thread, &(start, end)) in ranges.iter().enumerate() {
let offset = prefix[start];
let span = prefix[end] - prefix[start];
let piece = &frames[start..end];
let cloned = file.try_clone().map_err(WriteError::Io)?;
guards.push(scope.spawn(move || -> Result<(), WriteError> {
#[cfg(feature = "parallel")]
crate::parser::pin_worker(thread);
if span == 0 {
return Ok(());
}
let page = 4096usize;
let map_off = offset & !(page - 1);
let delta = offset - map_off;
let map_len = delta + span;
let mut map = unsafe {
memmap2::MmapOptions::new()
.offset(map_off as u64)
.len(map_len)
.map_mut(&cloned)
.map_err(WriteError::Io)?
};
advise_shared_store(&cloned, map.as_mut_ptr(), map.len());
let dest = &mut map[delta..delta + span];
format_group_into(piece, dest, precision, canonical, bind)?;
drop(map);
Ok(())
}));
}
for guard in guards {
match guard.join() {
Ok(result) => result?,
Err(_) => {
return Err(WriteError::Io(io::Error::other("frame formatter panicked")));
}
}
}
Ok(())
});
result
}
fn write_serial_frames(
path: &std::path::Path,
frames: &[&ConFrame],
precision: Option<usize>,
canonical: bool,
bind: bool,
) -> Result<(), WriteError> {
let (tmp, file) = temporary_sibling(path).map_err(WriteError::Io)?;
let mut writer = match precision {
Some(precision) => ConFrameWriter::with_precision(file, precision),
None => ConFrameWriter::new(file),
}
.canonical(canonical);
let result =
write_group(&mut writer, frames, bind).and_then(|_| writer.flush().map_err(WriteError::Io));
drop(writer);
commit_temporary(&tmp, path, result)
}
fn pipeline_frame_groups(
path: &std::path::Path,
frames: &[&ConFrame],
workers: usize,
precision: Option<usize>,
canonical: bool,
bind: bool,
) -> Result<(), WriteError> {
if (2..=workers).contains(&frames.len()) {
let mut parts: Vec<Result<Vec<u8>, WriteError>> = Vec::with_capacity(frames.len());
std::thread::scope(|scope| {
let mut guards = Vec::with_capacity(frames.len());
for (index, frame) in frames.iter().enumerate() {
guards.push(scope.spawn(move || {
#[cfg(feature = "parallel")]
crate::parser::pin_worker(index);
format_group_bytes(std::slice::from_ref(frame), 1, precision, canonical, bind)
}));
}
for guard in guards {
parts.push(match guard.join() {
Ok(part) => part,
Err(_) => Err(WriteError::Io(io::Error::other("frame formatter panicked"))),
});
}
});
let mut formatted = Vec::with_capacity(parts.len());
for part in parts {
formatted.push(part?);
}
let (tmp, file) = temporary_sibling(path).map_err(WriteError::Io)?;
let result = store_disjoint_parts(&file, &formatted).map_err(WriteError::Io);
return commit_temporary(&tmp, path, result);
}
let heavy = direct_slice_write(
workers,
frame_byte_weight(frames[0]).saturating_mul(frames.len()),
);
if heavy
&& let Some(lengths) = cached_packed_lengths(frames, workers, precision, canonical, bind)
{
let (tmp, file) = temporary_sibling(path).map_err(WriteError::Io)?;
let result = place_packed(&file, frames, &lengths, workers, precision, canonical, bind);
return commit_temporary(&tmp, path, result);
}
let groups = workers.saturating_mul(4).max(workers);
let ranges = group_frame_ranges(frames, groups);
let pieces: Vec<&[&ConFrame]> = ranges
.iter()
.map(|&(start, end)| &frames[start..end])
.collect();
let format_threads = pipeline_format_threads(workers).min(pieces.len().max(1));
let (tx, rx) = std::sync::mpsc::sync_channel::<(usize, Result<Vec<u8>, WriteError>)>(2);
let next_job = std::sync::atomic::AtomicUsize::new(0);
let (tmp, file) = temporary_sibling(path).map_err(WriteError::Io)?;
let result = std::thread::scope(|scope| {
for _ in 0..format_threads {
let tx = tx.clone();
let pieces = &pieces;
let next_job = &next_job;
scope.spawn(move || {
loop {
let index = next_job.fetch_add(1, std::sync::atomic::Ordering::Relaxed);
if index >= pieces.len() {
break;
}
let message = format_group_bytes(pieces[index], 1, precision, canonical, bind);
if tx.send((index, message)).is_err() {
break;
}
}
});
}
drop(tx);
let mut file = file;
let mut next = 0usize;
let mut pending = BTreeMap::<usize, Vec<u8>>::new();
while next < pieces.len() {
let (index, part) = rx
.recv()
.map_err(|err| WriteError::Io(io::Error::other(err.to_string())))?;
pending.insert(index, part?);
while let Some(part) = pending.remove(&next) {
file.write_all(&part).map_err(WriteError::Io)?;
next += 1;
}
}
file.flush().map_err(WriteError::Io)?;
Ok(())
});
commit_temporary(&tmp, path, result)
}
struct PackedCache {
first: usize,
last: usize,
count: usize,
precision: Option<usize>,
canonical: bool,
bind: bool,
generation: u64,
lengths: Vec<usize>,
}
static PACKED_GEN: std::sync::atomic::AtomicU64 = std::sync::atomic::AtomicU64::new(1);
static PACKED_FIRST: std::sync::atomic::AtomicUsize = std::sync::atomic::AtomicUsize::new(0);
static PACKED_LAST: std::sync::atomic::AtomicUsize = std::sync::atomic::AtomicUsize::new(0);
pub(crate) fn invalidate_packed_lengths() {
PACKED_GEN.fetch_add(1, std::sync::atomic::Ordering::Relaxed);
}
pub(crate) fn invalidate_packed_lengths_if(frame: usize) {
let first = PACKED_FIRST.load(std::sync::atomic::Ordering::Relaxed);
if first != 0
&& (frame == first || frame == PACKED_LAST.load(std::sync::atomic::Ordering::Relaxed))
{
invalidate_packed_lengths();
PACKED_FIRST.store(0, std::sync::atomic::Ordering::Relaxed);
}
}
fn cached_packed_lengths(
frames: &[&ConFrame],
workers: usize,
precision: Option<usize>,
canonical: bool,
bind: bool,
) -> Option<Vec<usize>> {
let first = frames[0] as *const ConFrame as usize;
let last = *frames.last()? as *const ConFrame as usize;
let count = frames.len();
static CACHE: std::sync::Mutex<Option<PackedCache>> = std::sync::Mutex::new(None);
let generation = PACKED_GEN.load(std::sync::atomic::Ordering::Relaxed);
if let Some(hit) = CACHE.lock().unwrap_or_else(|err| err.into_inner()).as_ref()
&& hit.generation == generation
&& hit.first == first
&& hit.last == last
&& hit.count == count
&& hit.precision == precision
&& hit.canonical == canonical
&& hit.bind == bind
{
return Some(hit.lengths.clone());
}
let generation_before = PACKED_GEN.load(std::sync::atomic::Ordering::Relaxed);
let lengths = packed_frame_lengths(frames, workers, precision, canonical, bind)?;
let generation_after = PACKED_GEN.load(std::sync::atomic::Ordering::Relaxed);
if generation_before == generation_after {
PACKED_FIRST.store(first, std::sync::atomic::Ordering::Relaxed);
PACKED_LAST.store(last, std::sync::atomic::Ordering::Relaxed);
*CACHE.lock().unwrap_or_else(|err| err.into_inner()) = Some(PackedCache {
first,
last,
count,
precision,
canonical,
bind,
generation: generation_after,
lengths: lengths.clone(),
});
}
Some(lengths)
}
fn store_disjoint_parts(file: &File, parts: &[Vec<u8>]) -> io::Result<()> {
let mut offsets = Vec::with_capacity(parts.len());
let mut total = 0u64;
for part in parts {
offsets.push(total);
total += part.len() as u64;
}
file.set_len(total)?;
if total == 0 {
return Ok(());
}
if file_is_tmpfs(file) {
store_mapped_parts(file, parts, &offsets)
} else {
std::thread::scope(|scope| {
let mut guards = Vec::new();
for (part, at) in parts.iter().zip(offsets) {
let cloned = file.try_clone()?;
guards.push(scope.spawn(move || pwrite_full(&cloned, part, at)));
}
for guard in guards {
guard
.join()
.unwrap_or_else(|_| Err(io::Error::other("pwrite worker")))?;
}
Ok(())
})
}
}
fn store_mapped_parts(file: &File, parts: &[Vec<u8>], offsets: &[u64]) -> io::Result<()> {
let total = offsets
.last()
.copied()
.unwrap_or(0)
.saturating_add(parts.last().map(|part| part.len() as u64).unwrap_or(0));
let mut map = unsafe {
memmap2::MmapOptions::new()
.len(total as usize)
.map_mut(file)?
};
let base = map.as_mut_ptr() as usize;
std::thread::scope(|scope| {
for (part, at) in parts.iter().zip(offsets) {
let src = part.as_ptr() as usize;
let len = part.len();
let dest = base + *at as usize;
scope.spawn(move || {
if len == 0 {
return;
}
unsafe {
std::ptr::copy_nonoverlapping(src as *const u8, dest as *mut u8, len);
}
});
}
});
Ok(())
}
fn write_group<W: Write>(
writer: &mut ConFrameWriter<W>,
group: &[&ConFrame],
bind: bool,
) -> Result<(), WriteError> {
if bind {
for frame in group {
writer.write_frame_bound(frame)?;
}
return Ok(());
}
writer.extend(group.iter().copied())
}
fn format_group_bytes(
group: &[&ConFrame],
workers: usize,
precision: Option<usize>,
canonical: bool,
bind: bool,
) -> Result<Vec<u8>, WriteError> {
let mut buffer = Vec::new();
let mut writer = match precision {
Some(precision) => ConFrameWriter::with_precision(&mut buffer, precision),
None => ConFrameWriter::new(&mut buffer),
}
.canonical(canonical)
.retain_scratch(true);
if workers >= 2 {
writer = writer.parallel_atoms(workers);
}
write_group(&mut writer, group, bind)?;
Ok(std::mem::take(&mut writer.scratch))
}
fn formatted_parts(
frames: &[&ConFrame],
workers: usize,
precision: Option<usize>,
canonical: bool,
bind: bool,
) -> Result<Vec<Vec<u8>>, WriteError> {
let wide = frames.len() == 1 && frames[0].atom_data.len() >= WIDE_ATOM_GRAIN;
if workers < 2 || (frames.len() < 2 && !wide) {
return Ok(vec![format_group_bytes(
frames, 1, precision, canonical, bind,
)?]);
}
if wide {
return Ok(vec![format_group_bytes(
frames, workers, precision, canonical, bind,
)?]);
}
#[cfg(feature = "parallel")]
{
let chunk = frames.len().div_ceil(workers).max(1);
return crate::parser::on_pool(workers, || {
use rayon::prelude::*;
frames
.par_chunks(chunk)
.map(|group| format_group_bytes(group, 1, precision, canonical, bind))
.collect()
});
}
#[cfg(not(feature = "parallel"))]
{
let _ = workers;
Ok(vec![format_group_bytes(
frames, 1, precision, canonical, bind,
)?])
}
}
fn pwrite_parts(path: &std::path::Path, parts: &[Vec<u8>]) -> io::Result<()> {
let slices: Vec<&[u8]> = parts.iter().map(|part| part.as_slice()).collect();
pwrite_borrowed(path, &slices)
}
fn pwrite_borrowed(path: &std::path::Path, parts: &[&[u8]]) -> io::Result<()> {
let mut offset = 0u64;
let mut offsets = Vec::with_capacity(parts.len());
for part in parts {
offsets.push(offset);
offset += part.len() as u64;
}
let file = OpenOptions::new()
.write(true)
.create(true)
.truncate(true)
.open(path)?;
file.set_len(offset)?;
std::thread::scope(|scope| {
let mut handles = Vec::new();
for (part, at) in parts.iter().zip(offsets) {
let file = file.try_clone()?;
handles.push(scope.spawn(move || pwrite_full(&file, part, at)));
}
for handle in handles {
handle
.join()
.unwrap_or_else(|_| Err(io::Error::other("pwrite worker")))?;
}
Ok(())
})
}
fn pwrite_full(file: &File, buf: &[u8], mut offset: u64) -> io::Result<()> {
if buf.is_empty() {
return Ok(());
}
#[cfg(miri)]
{
use std::io::{Seek, SeekFrom, Write};
let mut out = file.try_clone()?;
out.seek(SeekFrom::Start(offset))?;
out.write_all(buf)?;
return Ok(());
}
#[cfg(all(not(miri), unix))]
{
use std::os::unix::io::AsRawFd;
let fd = file.as_raw_fd();
let mut written = 0usize;
while written < buf.len() {
let n = unsafe {
libc::pwrite(
fd,
buf[written..].as_ptr() as *const libc::c_void,
buf.len() - written,
offset as libc::off_t,
)
};
if n < 0 {
let err = io::Error::last_os_error();
if err.kind() == io::ErrorKind::Interrupted {
continue;
}
return Err(err);
}
if n == 0 {
return Err(io::Error::new(io::ErrorKind::WriteZero, "pwrite"));
}
written += n as usize;
offset += n as u64;
}
Ok(())
}
#[cfg(all(not(miri), windows))]
{
use std::os::windows::fs::FileExt;
file.seek_write(buf, offset)?;
Ok(())
}
}
fn write_at_offset(file: &File, buf: &[u8], offset: u64) -> io::Result<usize> {
#[cfg(unix)]
{
use std::os::unix::fs::FileExt;
return file.write_at(buf, offset);
}
#[cfg(windows)]
{
use std::os::windows::fs::FileExt;
file.seek_write(buf, offset)
}
}
pub fn write_frames_to_string(frames: &[ConFrame]) -> Result<String, WriteError> {
let bytes = write_frames_to_vec(frames)?;
String::from_utf8(bytes).map_err(|err| {
WriteError::Io(io::Error::new(
io::ErrorKind::InvalidData,
err.utf8_error().to_string(),
))
})
}
impl ConFrameWriter<flate2::write::GzEncoder<File>> {
pub fn from_path_gzip<P: AsRef<Path>>(path: P) -> io::Result<Self> {
let encoder = crate::compression::gzip_writer(path.as_ref())?;
Ok(Self::new(encoder))
}
pub fn from_path_gzip_with_precision<P: AsRef<Path>>(
path: P,
precision: usize,
) -> io::Result<Self> {
let encoder = crate::compression::gzip_writer(path.as_ref())?;
Ok(Self::with_precision(encoder, precision))
}
}
#[cfg(feature = "zstd")]
impl ConFrameWriter<zstd::stream::write::AutoFinishEncoder<'static, File>> {
pub fn from_path_zstd<P: AsRef<Path>>(path: P) -> io::Result<Self> {
let encoder = crate::compression::zstd_writer(path.as_ref())?;
Ok(Self::new(encoder))
}
pub fn from_path_zstd_with_precision<P: AsRef<Path>>(
path: P,
precision: usize,
) -> io::Result<Self> {
let encoder = crate::compression::zstd_writer(path.as_ref())?;
Ok(Self::with_precision(encoder, precision))
}
}
#[cfg(test)]
mod float_format_tests {
use super::{encode_fixed_decimal_into, push_f64_prec, push_u64, write_finite};
fn formatted(v: f64, prec: usize) -> String {
let mut buf = Vec::new();
push_f64_prec(&mut buf, v, prec);
String::from_utf8(buf).expect("utf8")
}
#[test]
fn matches_std_fixed_precision() {
let vals = [
0.0,
-0.0,
1.0,
-1.0,
0.9045,
6.975_299_999_999_995,
63.546,
1.008,
15.3456,
90.0,
218.0,
1e-6,
-1.23456789,
10.0,
];
for prec in [0usize, 6] {
for v in vals {
let got = formatted(v, prec);
let exp = format!("{v:.prec$}");
assert_eq!(got, exp, "v={v:?} prec={prec}");
}
}
}
#[test]
fn fixed_decimal_matches_ryu_including_edges() {
let mut ryu_buf = ryu::Buffer::new();
let mut values = vec![
0.0,
-0.0,
1.0,
-1.0,
10.0,
97.0,
-97.0,
63.546,
0.5,
0.125,
1.5,
-0.1,
0.1 + 0.2,
1.0 / 3.0,
1e15,
-1e15,
1e16,
1e20,
-1e20,
1e-6,
1e-10,
f64::MIN_POSITIVE,
f64::MAX,
f64::MIN,
f64::INFINITY,
f64::NEG_INFINITY,
f64::NAN,
];
let mut state = 0x1234_5678_9abc_def0u64;
for _ in 0..256 {
state = state.wrapping_mul(0x5851_f42d_4c95_7f2d).wrapping_add(1);
values.push(f64::from_bits(state));
}
for value in values {
let mut got = Vec::new();
let mut scratch = ryu::Buffer::new();
super::push_round_trip(&mut got, &mut scratch, None, value);
let expect = ryu_buf.format(value);
if value.is_nan() {
assert_eq!(std::str::from_utf8(&got).unwrap(), expect);
continue;
}
assert_eq!(std::str::from_utf8(&got).unwrap(), expect, "{value:?}");
}
}
#[test]
fn round_trip_spelling_matches_ryu_on_repeated_and_unique_rows() {
use crate::iterators::ConFrameIterator;
use crate::writer::ConFrameWriter;
let mut text = String::from(
"bench\n{\"con_spec_version\":2}\n10.0 10.0 10.0\n90.0 90.0 90.0\n0 0\n0 0 1\n1\n6\n12.0\nCu\nCoordinates of Component 1\n",
);
for i in 0..6 {
let x = (i % 2) as f64 * 1.5 + 0.25;
let y = 3.25;
let z = if i == 3 {
-0.0
} else {
(i as f64) * 0.0010000001
};
text.push_str(&format!(" {x:.17} {y:.17} {z:.17} 0 {i}\n"));
}
text.push_str("\nCu\nVelocities of Component 1\n");
for i in 0..6 {
let vx = 0.125;
let vy = if i % 2 == 0 { 0.5 } else { -0.5 };
let vz = (i as f64) * 1.0e-4;
text.push_str(&format!(" {vx:.17} {vy:.17} {vz:.17} 0 {i}\n"));
}
let frames: Vec<_> = ConFrameIterator::new(&text)
.map(|frame| frame.expect("parse"))
.collect();
let mut values = Vec::new();
for atom in &frames[0].atom_data {
values.extend([atom.x(), atom.y(), atom.z()]);
}
for atom in &frames[0].atom_data {
if let Some(v) = atom.velocity() {
values.extend(v);
}
}
let mut buf = Vec::new();
{
let mut writer = ConFrameWriter::new(&mut buf);
writer.extend(frames.iter()).unwrap();
}
let written = String::from_utf8(buf).unwrap();
let mut ryu = ryu::Buffer::new();
let mut spelled = Vec::new();
for line in written.lines() {
let toks: Vec<_> = line.split_whitespace().collect();
if toks.len() >= 5
&& toks[0].contains('.')
&& toks[3].bytes().all(|b| b.is_ascii_digit())
{
spelled.extend(toks[..3].iter().map(|tok| (*tok).to_string()));
}
}
assert_eq!(spelled.len(), values.len());
for (got, value) in spelled.iter().zip(values) {
let expect = if value.is_finite() {
ryu.format_finite(value).to_string()
} else {
ryu.format(value).to_string()
};
assert_eq!(got, &expect, "value {value:?}");
}
}
#[test]
fn fixed_decimal_matches_ryu_or_falls_back_to_it() {
let mut values = vec![
0.0,
-0.0,
1.0,
-2.0,
10.0,
99.0,
100.0,
0.1,
0.01,
0.001,
0.0001,
0.0002,
-0.0001,
1.5,
0.25,
0.125,
63.546,
12.0,
90.0,
3.25,
0.5,
-0.5,
1.0e-4,
1.25,
218.0,
6.975_299_999_999_995,
1.008,
15.3456,
-1.23456789,
1e-6,
1e16,
1e20,
];
for i in 0..200 {
values.push(i as f64);
values.push((i as f64) * 0.001);
values.push((i as f64) * 0.0001);
values.push(-(i as f64) * 0.001);
}
let mut ryu = ryu::Buffer::new();
for value in values {
let mut buf = Vec::new();
write_finite(&mut buf, value);
let expect = ryu.format_finite(value);
assert_eq!(std::str::from_utf8(&buf).unwrap(), expect, "{value:?}");
if encode_fixed_decimal_into(value).is_none() {
continue;
}
}
}
#[test]
fn u64_matches_decimal_display() {
let mut values = vec![0, 1, 9, 10, 11, 99, 100, 101, 999, 1000, 50000, u64::MAX];
values.extend(0..200);
for n in values {
let mut buf = Vec::new();
push_u64(&mut buf, n);
assert_eq!(std::str::from_utf8(&buf).unwrap(), n.to_string());
}
}
}
#[cfg(test)]
mod section_arity_tests {
use super::require_minor2_section_arity;
use serde_json::json;
#[test]
fn minor2_lists_a_new_section_in_section_arity() {
let widgets = vec![json!("widgets")];
let err = require_minor2_section_arity(2, &widgets, &serde_json::Map::new()).unwrap_err();
let msg = err.to_string();
assert!(msg.contains("section_arity"), "{msg}");
let mut listed = serde_json::Map::new();
listed.insert("section_arity".into(), json!({"widgets": 3}));
require_minor2_section_arity(2, &widgets, &listed).unwrap();
let reserved = vec![json!("velocities"), json!("energies")];
require_minor2_section_arity(2, &reserved, &serde_json::Map::new()).unwrap();
require_minor2_section_arity(1, &widgets, &serde_json::Map::new()).unwrap();
}
}
#[cfg(test)]
mod weight_partition {
use super::pack_frames_by_weight;
use proptest::prelude::*;
fn is_ordered_partition(n: usize, workers: usize, ranges: &[(usize, usize)]) -> bool {
if n == 0 {
return ranges.is_empty();
}
let workers = workers.max(1).min(n);
if ranges.len() != workers {
return false;
}
let mut cursor = 0usize;
for &(start, end) in ranges {
if start != cursor || end <= start || end > n {
return false;
}
cursor = end;
}
cursor == n
}
proptest! {
#![proptest_config(ProptestConfig {
cases: 64,
failure_persistence: None,
..ProptestConfig::default()
})]
#[test]
fn byte_weights_do_not_drop_or_reorder_frames(
weights in prop::collection::vec(0usize..50_000, 0..48),
workers in 1usize..17,
) {
let ranges = pack_frames_by_weight(&weights, workers);
prop_assert!(
is_ordered_partition(weights.len(), workers, &ranges),
"weights={weights:?} workers={workers} ranges={ranges:?}"
);
}
}
}
#[cfg(all(test, feature = "parallel"))]
mod parallel_write_tests {
use super::{
direct_slice_write, frame_byte_weight, pack_frames_by_weight, pipeline_format_threads,
};
#[test]
fn two_workers_both_format_frame_groups() {
assert_eq!(pipeline_format_threads(1), 1);
assert_eq!(pipeline_format_threads(2), 2);
assert_eq!(pipeline_format_threads(8), 8);
}
#[test]
fn serial_and_narrow_writes_skip_the_direct_slice() {
let long = 3_000_000;
assert!(!direct_slice_write(1, long));
assert!(!direct_slice_write(2, long));
assert!(!direct_slice_write(4, long));
assert!(direct_slice_write(8, long));
assert!(!direct_slice_write(8, long - 1));
}
#[test]
fn a_heavy_frame_is_not_glued_to_a_light_one() {
assert_eq!(
pack_frames_by_weight(&[1, 1, 1, 100], 2),
vec![(0, 3), (3, 4)]
);
assert_eq!(
pack_frames_by_weight(&[100, 1, 1, 1], 2),
vec![(0, 1), (1, 4)]
);
assert_eq!(
pack_frames_by_weight(&[4, 4, 4, 4, 4, 4, 4, 4], 4),
vec![(0, 2), (2, 4), (4, 6), (6, 8)]
);
}
#[test]
fn velocity_lines_weigh_more_than_the_same_coordinates() {
let con = std::fs::read_to_string(concat!(
env!("CARGO_MANIFEST_DIR"),
"/resources/test/tiny_cuh2.con"
))
.unwrap();
let vel = std::fs::read_to_string(concat!(
env!("CARGO_MANIFEST_DIR"),
"/resources/test/tiny_cuh2.convel"
))
.unwrap();
let plain = crate::iterators::frames_from_text(&con, Some(1)).unwrap();
let moving = crate::iterators::frames_from_text(&vel, Some(1)).unwrap();
assert!(frame_byte_weight(&moving[0]) > frame_byte_weight(&plain[0]));
assert_eq!(
frame_byte_weight(&moving[0]),
plain[0].atom_data.len().max(1) * 2
);
}
#[test]
fn packed_groups_match_the_serial_bytes() {
let con = std::fs::read_to_string(concat!(
env!("CARGO_MANIFEST_DIR"),
"/resources/test/tiny_cuh2.con"
))
.unwrap();
let vel = std::fs::read_to_string(concat!(
env!("CARGO_MANIFEST_DIR"),
"/resources/test/tiny_cuh2.convel"
))
.unwrap();
let text = format!("{con}{con}{vel}");
let frames = crate::iterators::frames_from_text(&text, Some(1)).unwrap();
assert!(frames.len() > 2);
let refs: Vec<&crate::types::ConFrame> = frames.iter().collect();
let mut serial = Vec::new();
{
let mut writer = super::ConFrameWriter::new(&mut serial);
for frame in &frames {
writer.write_frame_bound(frame).unwrap();
}
writer.flush().unwrap();
}
let path = std::env::temp_dir().join(format!("rc-pack-{}.con", std::process::id()));
super::write_frame_refs_to_path(&path, &refs, 2, None, false, true).unwrap();
let got = std::fs::read(&path).unwrap();
assert_eq!(got, serial);
let _ = std::fs::remove_file(&path);
}
#[test]
fn mixed_con_and_convel_groups_stay_within_ten_percent() {
let con = std::fs::read_to_string(concat!(
env!("CARGO_MANIFEST_DIR"),
"/resources/test/tiny_cuh2.con"
))
.unwrap();
let vel = std::fs::read_to_string(concat!(
env!("CARGO_MANIFEST_DIR"),
"/resources/test/tiny_cuh2.convel"
))
.unwrap();
let text = format!("{}{}", con.repeat(24), vel.repeat(24));
let frames = crate::iterators::frames_from_text(&text, Some(1)).unwrap();
assert_eq!(frames.len(), 48);
let refs: Vec<&crate::types::ConFrame> = frames.iter().collect();
let weights: Vec<usize> = refs.iter().map(|frame| frame_byte_weight(frame)).collect();
let groups = 4usize;
let ranges = super::group_frame_ranges(&refs, groups);
let group_weight = |ranges: &[(usize, usize)]| -> Vec<u64> {
ranges
.iter()
.map(|&(start, end)| weights[start..end].iter().map(|w| *w as u64).sum())
.collect()
};
let weighted = group_weight(&ranges);
let total: u64 = weighted.iter().sum();
let mean = total as f64 / weighted.len() as f64;
let max = *weighted.iter().max().unwrap() as f64;
assert!(
max <= mean * 1.10,
"heaviest group {max} mean {mean} weights {weighted:?}"
);
let equal = group_weight(&super::equal_frame_ranges(frames.len(), groups));
let equal_max = *equal.iter().max().unwrap() as f64;
let equal_mean = total as f64 / equal.len() as f64;
assert!(
equal_max > equal_mean * 1.10,
"equal frame counts should miss the 10% band: {equal:?}"
);
let mut serial = Vec::new();
{
let mut writer = super::ConFrameWriter::new(&mut serial);
for frame in &frames {
writer.write_frame_bound(frame).unwrap();
}
writer.flush().unwrap();
}
let path = std::env::temp_dir().join(format!("rc-mix-{}.con", std::process::id()));
super::write_frame_refs_to_path(&path, &refs, 4, None, false, true).unwrap();
assert_eq!(std::fs::read(&path).unwrap(), serial);
let _ = std::fs::remove_file(&path);
}
#[test]
fn disjoint_parts_rebuild_the_concatenation() {
let path = std::env::temp_dir().join(format!("rc-parts-{}.con", std::process::id()));
let file = std::fs::OpenOptions::new()
.read(true)
.write(true)
.create(true)
.truncate(true)
.open(&path)
.unwrap();
let parts = vec![
b"aa".to_vec(),
b"bbb".to_vec(),
Vec::new(),
b"cdef".to_vec(),
];
super::store_disjoint_parts(&file, &parts).unwrap();
drop(file);
assert_eq!(std::fs::read(&path).unwrap(), b"aabbbcdef");
let _ = std::fs::remove_file(&path);
}
#[test]
fn round_trip_width_matches_the_spelling() {
let mut samples = vec![
0.0,
-0.0,
1.0,
-2.0,
20.0,
90.0,
63.546,
0.001,
0.000001,
1e-7,
1e-20,
1e22,
-1.25,
9.983,
f64::MIN_POSITIVE,
f64::MAX,
f64::NAN,
f64::INFINITY,
f64::NEG_INFINITY,
];
for i in 0..2000 {
samples.push((i as f64) * 0.001);
samples.push(-((i as f64) + 0.5) * 0.017);
}
let mut buf = Vec::new();
let mut ryu = ryu::Buffer::new();
for value in samples {
buf.clear();
super::push_round_trip(&mut buf, &mut ryu, None, value);
assert_eq!(
super::round_trip_width(value),
buf.len(),
"spelling {buf:?} of {value}"
);
}
}
#[test]
fn packed_lengths_match_formatted_frames() {
use crate::types::ConFrameBuilder;
let frame = |shift: f64, weird: f64| {
let mut builder = ConFrameBuilder::new([20.0, 20.0, 20.0], [90.0, 90.0, 90.0]);
builder
.set_units(serde_json::json!({
"length": "angstrom",
"energy": "eV",
"mass": "amu",
"time": crate::units::NATURAL_TIME_UNIT,
}))
.unwrap();
builder
.add_atom("Cu", shift, shift + 0.001, weird, [false; 3], 0, 63.546)
.with_velocity([shift * 0.1, -0.25, 0.0]);
builder
.add_atom("Cu", 1.5, -2.0, 0.0, [false; 3], 1, 63.546)
.with_velocity([0.0, 1e-20, -0.0]);
builder.build().unwrap()
};
let frames: Vec<_> = (0..6)
.map(|i| {
frame(
i as f64 * 0.017,
if i == 3 { 1.0e22 } else { 0.125 + i as f64 },
)
})
.collect();
let refs: Vec<&crate::types::ConFrame> = frames.iter().collect();
let lengths = super::packed_frame_lengths(&refs, 4, None, false, true)
.expect("uniform frames have exact lengths");
for (frame, len) in refs.iter().zip(&lengths) {
let formatted =
super::format_group_bytes(std::slice::from_ref(frame), 1, None, false, true)
.unwrap();
assert_eq!(formatted.len(), *len);
}
let mut dest = vec![0u8; lengths.iter().sum()];
let ptr = dest.as_mut_ptr();
super::format_group_into(&refs, &mut dest, None, false, true).unwrap();
assert_eq!(dest.as_ptr(), ptr);
let mut serial = Vec::new();
{
let mut writer = super::ConFrameWriter::new(&mut serial);
for frame in &frames {
writer.write_frame_bound(frame).unwrap();
}
writer.flush().unwrap();
}
assert_eq!(dest, serial);
let path = std::env::temp_dir().join(format!("rc-packed-{}.con", std::process::id()));
super::write_frame_refs_to_path(&path, &refs, 4, None, false, true).unwrap();
assert_eq!(std::fs::read(&path).unwrap(), serial);
let _ = std::fs::remove_file(&path);
let mut other = ConFrameBuilder::new([10.0, 10.0, 10.0], [90.0, 90.0, 90.0]);
other
.set_units(serde_json::json!({
"length": "angstrom",
"energy": "eV",
"mass": "amu",
"time": crate::units::NATURAL_TIME_UNIT,
}))
.unwrap();
other
.add_atom("H", 0.0, 0.0, 0.0, [false; 3], 0, 1.0)
.with_velocity([0.0, 0.0, 0.0]);
let odd = other.build().unwrap();
let mixed = [refs[0], &odd];
assert!(super::packed_frame_lengths(&mixed, 2, None, false, true).is_none());
super::write_frame_refs_to_path(&path, &mixed, 2, None, false, true).unwrap();
let mut serial = Vec::new();
{
let mut writer = super::ConFrameWriter::new(&mut serial);
writer.write_frame_bound(mixed[0]).unwrap();
writer.write_frame_bound(mixed[1]).unwrap();
writer.flush().unwrap();
}
assert_eq!(std::fs::read(&path).unwrap(), serial);
let _ = std::fs::remove_file(&path);
}
#[test]
fn a_coordinate_edit_reprices_the_packed_lengths() {
use crate::types::ConFrameBuilder;
let mut frames = Vec::new();
for i in 0..4 {
let mut builder = ConFrameBuilder::new([20.0, 20.0, 20.0], [90.0, 90.0, 90.0]);
builder
.set_units(serde_json::json!({
"length": "angstrom",
"energy": "eV",
"mass": "amu",
"time": crate::units::NATURAL_TIME_UNIT,
}))
.unwrap();
builder
.add_atom("Cu", i as f64, 0.25, 0.5, [false; 3], 0, 63.546)
.with_velocity([0.1, 0.0, 0.0]);
builder
.add_atom("Cu", 1.5, -2.0, 0.0, [false; 3], 1, 63.546)
.with_velocity([0.0, 0.2, 0.0]);
frames.push(builder.build().unwrap());
}
let refs: Vec<&crate::types::ConFrame> = frames.iter().collect();
let before = super::cached_packed_lengths(&refs, 2, None, false, true).unwrap();
let old = before[2];
frames[2].atom_data[0].set_x(1.0e22);
let refs: Vec<&crate::types::ConFrame> = frames.iter().collect();
let after = super::cached_packed_lengths(&refs, 2, None, false, true).unwrap();
let formatted =
super::format_group_bytes(std::slice::from_ref(&refs[2]), 1, None, false, true)
.unwrap();
assert_ne!(old, formatted.len());
assert_eq!(after[2], formatted.len());
}
}
#[cfg(all(test, feature = "parallel"))]
mod parallel_write_bytes {
use super::{write_frames_to_vec, write_frames_to_vec_with_threads};
use crate::iterators::frames_from_text;
fn wide_text(n: usize) -> String {
let mut text = String::new();
text.push_str(
"grid\n{\"con_spec_version\":2}\n10.0 10.0 10.0\n90.0 90.0 90.0\n0 0\n0 0\n1\n",
);
text.push_str(&format!("{n}\n63.546\nCu\nCoordinates of Component 1\n"));
for i in 0..n {
text.push_str(&format!("{i}.125 0.25 0.5 0 {i}\n"));
}
text
}
#[test]
fn parallel_bytes_match_the_serial_writer() {
let wide = frames_from_text(&wide_text(5000), Some(1)).unwrap();
let serial = write_frames_to_vec(&wide).unwrap();
let parallel = write_frames_to_vec_with_threads(&wide, 4, None, false).unwrap();
assert_eq!(parallel, serial);
let mut with_vel = wide_text(5000);
with_vel.push_str("\nCu\nVelocities of Component 1\n");
for i in 0..5000 {
with_vel.push_str(&format!("{}.0001 0.0002 -0.0001 0 {i}\n", i % 7));
}
let moving = frames_from_text(&with_vel, Some(1)).unwrap();
let serial = write_frames_to_vec(&moving).unwrap();
let parallel = write_frames_to_vec_with_threads(&moving, 4, None, false).unwrap();
assert_eq!(parallel, serial);
let one = frames_from_text(
&std::fs::read_to_string(concat!(
env!("CARGO_MANIFEST_DIR"),
"/resources/test/cuh2.con"
))
.unwrap()
.repeat(8),
Some(1),
)
.unwrap();
let serial = write_frames_to_vec(&one).unwrap();
let parallel = write_frames_to_vec_with_threads(&one, 4, None, false).unwrap();
assert_eq!(parallel, serial);
}
#[test]
fn pwrite_matches_the_concatenated_writer() {
let one = frames_from_text(
&std::fs::read_to_string(concat!(
env!("CARGO_MANIFEST_DIR"),
"/resources/test/cuh2.con"
))
.unwrap()
.repeat(8),
Some(1),
)
.unwrap();
let serial = write_frames_to_vec(&one).unwrap();
let refs: Vec<&crate::types::ConFrame> = one.iter().collect();
let path = std::env::temp_dir().join(format!(
"readcon-pwrite-{}-{}.con",
std::process::id(),
one.len()
));
super::write_frame_refs_to_path(&path, &refs, 4, None, false, false).unwrap();
let got = std::fs::read(&path).unwrap();
let _ = std::fs::remove_file(&path);
assert_eq!(got, serial);
}
#[test]
fn stamped_pwrite_matches_the_binding_writer() {
let one = frames_from_text(
&std::fs::read_to_string(concat!(
env!("CARGO_MANIFEST_DIR"),
"/resources/test/cuh2.con"
))
.unwrap()
.repeat(8),
Some(1),
)
.unwrap();
let mut serial = Vec::new();
{
let mut writer = super::ConFrameWriter::new(&mut serial);
for frame in &one {
let mut stamped = frame.clone();
writer.write_frame_binding(&mut stamped).unwrap();
}
}
let refs: Vec<&crate::types::ConFrame> = one.iter().collect();
let path = std::env::temp_dir().join(format!(
"readcon-pwrite-bind-{}-{}.con",
std::process::id(),
one.len()
));
super::write_frame_refs_to_path(&path, &refs, 4, None, false, true).unwrap();
let got = std::fs::read(&path).unwrap();
let _ = std::fs::remove_file(&path);
assert_eq!(got, serial);
}
#[test]
fn a_wide_file_write_matches_the_serial_bytes_when_ids_are_not_dense() {
let n = 4096usize;
let mut text = String::from(
"grid\n{\"con_spec_version\":2}\n10.0 10.0 10.0\n90.0 90.0 90.0\n0 0\n0 0\n1\n",
);
text.push_str(&format!("{n}\n63.546\nCu\nCoordinates of Component 1\n"));
for i in 0..n {
text.push_str(&format!("{i}.125 0.25 0.5 0 {}\n", i + 10));
}
text.push_str("\nCu\nVelocities of Component 1\n");
for i in 0..n {
text.push_str(&format!("{}.0001 0.0002 -0.0001 0 {}\n", i % 7, i + 10));
}
let frames = frames_from_text(&text, Some(1)).unwrap();
let serial = write_frames_to_vec(&frames).unwrap();
let parallel = write_frames_to_vec_with_threads(&frames, 8, None, false).unwrap();
assert_eq!(parallel, serial);
let refs: Vec<&crate::types::ConFrame> = frames.iter().collect();
let path = std::env::temp_dir().join(format!(
"readcon-wide-stream-{}-{n}.con",
std::process::id()
));
super::write_frame_refs_to_path(&path, &refs, 8, None, false, false).unwrap();
let got = std::fs::read(&path).unwrap();
let _ = std::fs::remove_file(&path);
assert_eq!(got, serial);
}
#[test]
fn a_mixed_wide_frame_matches_the_serial_parse() {
let n = 7000usize;
let mut text = String::from(
"grid\n{\"con_spec_version\":2}\n10.0 10.0 10.0\n90.0 90.0 90.0\n0 0\n0 0\n2\n",
);
text.push_str(&format!("{n} {n}\n63.546 1.008\n"));
text.push_str("Cu\nCoordinates of Component 1\n");
for i in 0..n {
text.push_str(&format!("{i}.125 0.25 0.5 0 {i}\n"));
}
text.push_str("H\nCoordinates of Component 2\n");
for i in 0..n {
let id = n + i;
text.push_str(&format!("1.5 2.25 {i}.5 0 {id}\n"));
}
assert!(text.len() >= (256 << 10), "frame is {} bytes", text.len());
let serial = frames_from_text(&text, Some(1)).unwrap();
let parallel = frames_from_text(&text, Some(8)).unwrap();
assert_eq!(parallel.len(), 1);
assert_eq!(parallel[0].atom_data.len(), n * 2);
assert_eq!(parallel[0].atom_data[0].symbol.as_ref(), "Cu");
assert_eq!(parallel[0].atom_data[n].symbol.as_ref(), "H");
assert_eq!(
parallel[0].atom_data[n].atom_id,
serial[0].atom_data[n].atom_id
);
assert_eq!(parallel, serial);
}
#[test]
fn coordinate_lines_match_the_serial_frame_body() {
let frames = frames_from_text(&wide_text(32), Some(1)).unwrap();
let lines = super::format_coordinate_lines(&frames[0].atom_data).unwrap();
let serial = write_frames_to_vec(&frames).unwrap();
let text = std::str::from_utf8(&serial).unwrap();
let marker = "Coordinates of Component 1\n";
let start = text.find(marker).unwrap() + marker.len();
assert_eq!(lines, text.as_bytes()[start..]);
}
#[test]
fn pwrite_offsets_follow_an_exclusive_scan() {
let path =
std::env::temp_dir().join(format!("readcon-pwrite-scan-{}.bin", std::process::id()));
let parts: [&[u8]; 4] = [b"abc", b"", b"def\n", b"g"];
super::pwrite_borrowed(&path, &parts).unwrap();
let got = std::fs::read(&path).unwrap();
let _ = std::fs::remove_file(&path);
assert_eq!(got, b"abcdef\ng");
}
#[test]
fn mapped_write_matches_the_serial_write() {
use crate::types::ConFrameBuilder;
for n in [1usize, 3, 7, 4097, 8191] {
let mut builder = ConFrameBuilder::new([10.0, 10.0, 10.0], [90.0, 90.0, 90.0]);
for i in 0..n {
let x = (i % 97) as f64 * 0.125;
let y = (i % 13) as f64;
let z = if i % 2 == 0 { i as f64 * 0.001 } else { -0.5 };
builder.add_atom("Cu", x, y, z, [false, i % 5 == 0, false], i as u64, 63.546);
}
let frame = builder.build().expect("frame");
let serial = super::write_frames_to_vec(std::slice::from_ref(&frame)).expect("serial");
for workers in [1usize, 2, 3, 8] {
let path = std::env::temp_dir().join(format!(
"readcon-map-{n}-{workers}-{}.con",
std::process::id()
));
let mut writer = super::ConFrameWriter::new(Vec::<u8>::new())
.parallel_atoms(workers)
.map_rows(true);
writer.stream_file = Some(
std::fs::OpenOptions::new()
.read(true)
.write(true)
.create(true)
.truncate(true)
.open(&path)
.expect("create"),
);
writer.write_frame(&frame).expect("mapped");
drop(writer);
let got = std::fs::read(&path).expect("read");
let _ = std::fs::remove_file(&path);
assert_eq!(got, serial, "n={n} workers={workers}");
}
}
}
fn copper_frame(
n: usize,
symbol: &str,
mass: f64,
id_at: impl Fn(usize) -> u64,
) -> crate::types::ConFrame {
use crate::types::ConFrameBuilder;
let mut builder = ConFrameBuilder::new([10.0, 10.0, 10.0], [90.0, 90.0, 90.0]);
for i in 0..n {
builder.add_atom(
symbol,
i as f64 * 0.125,
0.25,
0.5,
[false, false, false],
id_at(i),
mass,
);
}
builder.build().expect("frame")
}
#[test]
fn two_wide_frames_on_two_threads_keep_up_with_one_frame() {
let n = 200_000usize;
let one = copper_frame(n, "Cu", 63.546, |i| i as u64);
let frames = [one.clone(), one.clone()];
let refs: Vec<&crate::types::ConFrame> = frames.iter().collect();
let serial = super::write_frames_to_vec(&frames).expect("serial");
let path =
std::env::temp_dir().join(format!("readcon-weak-write-{}.con", std::process::id()));
super::write_frame_refs_to_path(&path, &refs, 2, None, false, false).expect("parallel");
let got = std::fs::read(&path).expect("read");
let _ = std::fs::remove_file(&path);
assert_eq!(got, serial);
let time = |frames: &[crate::types::ConFrame], workers: usize| {
let refs: Vec<&crate::types::ConFrame> = frames.iter().collect();
let path = std::env::temp_dir().join(format!(
"readcon-weak-time-{}-{workers}.con",
std::process::id()
));
let _ = super::write_frame_refs_to_path(&path, &refs, workers, None, false, false);
let mut samples = Vec::with_capacity(5);
for _ in 0..5 {
let _ = std::fs::remove_file(&path);
let started = std::time::Instant::now();
super::write_frame_refs_to_path(&path, &refs, workers, None, false, false)
.expect("write");
samples.push(started.elapsed().as_secs_f64());
}
let _ = std::fs::remove_file(&path);
samples.sort_by(|left, right| left.partial_cmp(right).unwrap());
samples[samples.len() / 2]
};
let one_s = time(std::slice::from_ref(&one), 1);
let two_s = time(&frames, 2);
#[cfg(not(coverage))]
assert!(
two_s < one_s * 1.7,
"two frames on two threads took {two_s:.4}s against one frame {one_s:.4}s"
);
#[cfg(coverage)]
let _ = (one_s, two_s);
}
#[test]
fn parallel_write_rejects_an_id_repeated_in_another_type() {
use crate::types::ConFrameBuilder;
let n = 4096usize;
let mut builder = ConFrameBuilder::new([10.0, 10.0, 10.0], [90.0, 90.0, 90.0]);
for i in 0..n {
builder.add_atom(
"Cu",
i as f64 * 0.125,
0.25,
0.5,
[false, false, false],
i as u64,
63.546,
);
}
for i in 0..n {
builder.add_atom(
"H",
i as f64 * 0.125,
1.5,
0.5,
[false, false, false],
i as u64,
1.008,
);
}
let frame = builder.build().expect("frame");
let serial = write_frames_to_vec(std::slice::from_ref(&frame));
let parallel =
write_frames_to_vec_with_threads(std::slice::from_ref(&frame), 4, None, false);
assert!(serial.is_err(), "serial writer accepted a repeated id");
assert!(parallel.is_err(), "parallel writer accepted a repeated id");
let path =
std::env::temp_dir().join(format!("readcon-repeat-id-{}.con", std::process::id()));
std::fs::write(&path, b"keep").unwrap();
let refs = [&frame];
let err = super::write_frame_refs_to_path(&path, &refs, 4, None, false, false);
assert!(err.is_err(), "path write accepted a repeated id");
assert_eq!(std::fs::read(&path).unwrap(), b"keep");
let _ = std::fs::remove_file(&path);
}
#[test]
fn a_failed_later_group_leaves_the_previous_file() {
let good = copper_frame(2, "Cu", 63.546, |i| i as u64);
let bad = copper_frame(2, "Cu", 63.546, |_| 0);
let path =
std::env::temp_dir().join(format!("readcon-group-fail-{}.con", std::process::id()));
std::fs::write(&path, b"keep").unwrap();
let refs = [&good, &bad];
let err = super::write_frame_refs_to_path(&path, &refs, 4, None, false, false);
assert!(err.is_err());
assert_eq!(std::fs::read(&path).unwrap(), b"keep");
let _ = std::fs::remove_file(&path);
}
#[test]
fn mapped_decimal_lines_match_the_serial_writer() {
let frame = copper_frame(64, "Cu", 63.546, |i| i as u64);
let frames = [frame];
let serial = write_frames_to_vec_with_threads(&frames, 1, Some(40), false).unwrap();
let path =
std::env::temp_dir().join(format!("readcon-long-line-{}.con", std::process::id()));
let mut writer = super::ConFrameWriter::with_precision(Vec::<u8>::new(), 40)
.parallel_atoms(4)
.map_rows(true);
writer.stream_file = Some(
std::fs::OpenOptions::new()
.read(true)
.write(true)
.create(true)
.truncate(true)
.open(&path)
.unwrap(),
);
writer.write_frame(&frames[0]).unwrap();
drop(writer);
let got = std::fs::read(&path).unwrap();
let _ = std::fs::remove_file(&path);
assert_eq!(got, serial);
}
#[test]
fn a_rejected_spread_does_not_keep_the_reserved_file() {
use crate::types::ConFrameBuilder;
let n = 4096usize;
let mut builder = ConFrameBuilder::new([10.0, 10.0, 10.0], [90.0, 90.0, 90.0]);
for i in 0..n {
builder.add_atom(
"Cu",
i as f64 * 0.125,
0.25,
0.5,
[false, false, false],
i as u64,
63.546,
);
}
builder.with_spread([0.0, 0.0, 0.0]);
let mut frame = builder.build().expect("frame");
frame.atom_data[0].set_spread(Some([-1.0, 0.0, 0.0]));
let path =
std::env::temp_dir().join(format!("readcon-spread-fail-{}.con", std::process::id()));
let file = std::fs::OpenOptions::new()
.read(true)
.write(true)
.create(true)
.truncate(true)
.open(&path)
.unwrap();
let mut writer = super::ConFrameWriter::new(Vec::<u8>::new())
.parallel_atoms(4)
.map_rows(true)
.mapped_file(file);
let err = writer.write_frame(&frame);
drop(writer);
assert!(err.is_err());
assert_eq!(std::fs::metadata(&path).unwrap().len(), 0);
let _ = std::fs::remove_file(&path);
}
#[test]
fn overlapped_pwrite_does_not_keep_bytes_past_the_frame() {
use std::io::{Seek, Write};
let frame = copper_frame(4096, "Cu", 63.546, |i| i as u64);
let serial = write_frames_to_vec(std::slice::from_ref(&frame)).unwrap();
let path = std::env::current_dir()
.unwrap()
.join(format!("readcon-pwrite-tail-{}.con", std::process::id()));
let mut file = std::fs::OpenOptions::new()
.read(true)
.write(true)
.create(true)
.truncate(true)
.open(&path)
.unwrap();
file.write_all(&vec![b'Z'; serial.len() + 64]).unwrap();
file.seek(std::io::SeekFrom::Start(0)).unwrap();
let workers = super::physical_cores().max(2);
let mut writer = super::ConFrameWriter::new(Vec::<u8>::new())
.parallel_atoms(workers)
.mapped_file(file);
writer.write_frame(&frame).unwrap();
drop(writer);
let got = std::fs::read(&path).unwrap();
let _ = std::fs::remove_file(&path);
assert_eq!(got, serial);
}
}
#[cfg(test)]
mod gzip_stream_memory {
use super::ConFrameWriter;
use crate::iterators::{frames_from_text, read_all_frames};
fn parsed_frame() -> crate::types::ConFrame {
let mut text = String::from(
"grid\n{\"con_spec_version\":2,\"units\":{\"length\":\"angstrom\",\"energy\":\"eV\"}}\n20 20 20\n90 90 90\n0 0\n0 0\n1\n48\n63.546\nCu\nCoordinates of Component 1\n",
);
for i in 0..48 {
text.push_str(&format!("0.{i} 1 2 0 {i}\n"));
}
frames_from_text(&text, Some(1)).unwrap().remove(0)
}
#[test]
fn parsed_frames_gzip_write_drops_the_metadata_maps() {
let frames = vec![parsed_frame(); 4000];
assert!(
frames
.iter()
.all(|frame| !frame.header.metadata.has_cached_map())
);
let path = std::env::temp_dir().join(format!(
"readcon-gzip-peak-{}-{}.con.gz",
std::process::id(),
std::time::SystemTime::now()
.duration_since(std::time::UNIX_EPOCH)
.unwrap()
.as_nanos()
));
{
let mut writer = ConFrameWriter::from_path_gzip(&path).unwrap();
for frame in &frames {
writer.write_frame(frame).unwrap();
}
writer.flush().unwrap();
}
let cached = frames
.iter()
.filter(|frame| frame.header.metadata.has_cached_map())
.count();
assert_eq!(
cached, 0,
"the write retained {cached} metadata maps; each is about a kilobyte"
);
let back = read_all_frames(&path).unwrap();
assert_eq!(back.len(), frames.len());
assert_eq!(back[0].atom_data.len(), 48);
assert_eq!(back[0].atom_data[3].atom_id, 3);
assert_eq!(back[17].atom_data[0].symbol.as_ref(), "Cu");
let _ = std::fs::remove_file(&path);
}
#[test]
fn lattice_line_overrides_box_lengths_and_can_be_reread_from_the_cache() {
let text = "\
grid
{\"con_spec_version\":2,\"lattice_vectors\":[[10,0,0],[0,20,0],[0,0,30]]}
5 5 5
90 90 90
0 0
0 0
1
1
1.0
Cu
Coordinates of Component 1
0 0 0 0 1
";
let frame = frames_from_text(text, Some(1)).unwrap().remove(0);
assert!(!frame.header.metadata.has_cached_map());
let vecs = frame.header.lattice_vectors().unwrap();
assert_eq!(vecs, [[10.0, 0.0, 0.0], [0.0, 20.0, 0.0], [0.0, 0.0, 30.0]]);
assert!(
!frame.header.metadata.has_cached_map(),
"reading lattice vectors cached the whole metadata map"
);
let path =
std::env::temp_dir().join(format!("readcon-lattice-line-{}.con", std::process::id()));
{
let mut writer = ConFrameWriter::from_path(&path).unwrap();
writer.write_frame(&frame).unwrap();
writer.flush().unwrap();
}
let back = read_all_frames(&path).unwrap();
let _ = std::fs::remove_file(&path);
assert_eq!(back[0].header.boxl, [10.0, 20.0, 30.0]);
assert!(frame.header.metadata.get("lattice_vectors").is_some());
assert!(frame.header.metadata.has_cached_map());
assert_eq!(frame.header.lattice_vectors().unwrap()[1][1], 20.0);
}
#[test]
fn repeated_lattice_lookups_share_a_bounded_cache() {
let mut frames = Vec::new();
for side in [10.0, 11.0, 12.0, 13.0, 14.0] {
let text = format!(
"grid\n{{\"con_spec_version\":2,\"lattice_vectors\":[[{side},0,0],[0,{side},0],[0,0,{side}]]}}\n5 5 5\n90 90 90\n0 0\n0 0\n1\n1\n1.0\nCu\nCoordinates of Component 1\n0 0 0 0 1\n"
);
frames.push(frames_from_text(&text, Some(1)).unwrap().remove(0));
}
for _ in 0..2 {
for (index, frame) in frames.iter().enumerate() {
let side = 10.0 + index as f64;
let vecs = frame.header.lattice_vectors().unwrap();
assert_eq!(vecs[0][0], side);
assert_eq!(vecs[2][2], side);
assert!(!frame.header.metadata.has_cached_map());
}
}
let older = frames[1].header.lattice_vectors().unwrap();
assert_eq!(older[0][0], 11.0);
let text = "\
grid
{\"con_spec_version\":2,\"lattice_vectors\":[[14,0,0],[0,14,0],[0,0,14]]}
5 5 5
90 90 90
0 0
0 0
1
1
1.0
Cu
Coordinates of Component 1
0 0 0 0 1
";
let twin = frames_from_text(text, Some(1)).unwrap().remove(0);
assert_eq!(twin.header.lattice_vectors().unwrap()[0][0], 14.0);
assert!(!twin.header.metadata.has_cached_map());
}
#[test]
fn lattice_lookups_from_several_threads_keep_their_own_lines() {
let frames: Vec<_> = (1..17)
.map(|side| {
let text = format!(
"grid\n{{\"con_spec_version\":2,\"lattice_vectors\":[[{side},0,0],[0,{side},0],[0,0,{side}]]}}\n5 5 5\n90 90 90\n0 0\n0 0\n1\n1\n1.0\nCu\nCoordinates of Component 1\n0 0 0 0 1\n"
);
frames_from_text(&text, Some(1)).unwrap().remove(0)
})
.collect();
std::thread::scope(|scope| {
for _ in 0..4 {
scope.spawn(|| {
for _ in 0..8 {
for (index, frame) in frames.iter().enumerate() {
let side = 1.0 + index as f64;
let vecs = frame.header.lattice_vectors().unwrap();
assert_eq!(vecs[0][0], side);
assert_eq!(vecs[1][1], side);
assert_eq!(vecs[2][2], side);
}
}
});
}
});
}
}