Structure analysis

The amorphgen.analysis module provides ensemble structural analysis for amorphous structure files: pair distribution functions, structure factors, total correlation functions, coordination numbers, bond angles, bond-orientational order, ring statistics, Voronoi metrics, void distributions, oxygen speciation, elastic moduli, harmonic vibrational DOS, energy ranking, and validation against literature reference ranges.

The CLI entry point is amorphgen --analyse <FILE_OR_DIR>; the Python API is the StructureAnalyser class.

StructureAnalyser

class amorphgen.analysis.StructureAnalyser(source, cutoff='auto-rdf')[source]

Bases: object

Analyse amorphous structures: density, coordination, distances, angles, RDF, S(q), rings, Voronoi, voids, oxygen speciation, crystal-like bond order, elastic moduli, vibrational DOS and energy ranking.

Results retain per-structure observations and uncertainty summaries of equal-weight structure means (SEM, Student-t intervals and percentile bootstrap bounds). Pooled site/bond/angle std fields describe spread, not uncertainty. Curve bands are pointwise. Independent structures are assumed; fewer than two contributing structures leave intervals undefined.

Parameters:
  • source (str or list of str) – Path to a structure file, a directory, or a list of file paths.

  • cutoff (float, dict, or str) –

    • float: single cutoff (A) for all pairs

    • dict: pair-specific, e.g. {“Si-O”: 2.2}

    • ”auto-rdf”: from RDF first minimum (default — physically correct first coordination shell, matches the standard convention in neutron-diffraction analysis of glasses/liquids)

    • ”auto”: from bonding radii (fast; can truncate the first peak for systems with broad bond distributions such as a-Si, a-HfO2, chalcogenides — prefer “auto-rdf” unless you have a specific reason)

__init__(source, cutoff='auto-rdf')[source]
Parameters:
  • source (str, list of str, or list of Atoms) – Path to a structure file, a directory of structure files, a list of file paths, or a list of ASE Atoms objects.

  • cutoff (float, dict, or str) –

    Neighbour cutoff for CN, bond distances, and angles.

    • float: single cutoff (A) for all pairs

    • str with overrides: "In-O=2.6,Zn-O=2.3" keeps auto-rdf for every other pair; "auto,In-O=2.6" or "2.4,In-O=2.6" set the base rule explicitly

    • dict listing only some pairs: the rest come from auto-rdf (or from its "default" entry: a mode name or a number)

    • dict: pair-specific, e.g. {"Si-O": 2.2, "O-O": 2.8}

    • "auto-rdf" (default): first minimum of the partial RDF — the standard physical definition of the first coordination shell. Robust across material classes; prefer this for published numbers.

    • "auto": derived from Shannon/Cordero/Goldschmidt radii (minsep). Fast but can under-count coordination for systems with broad first-shell distributions (a-Si, a-HfO2, chalcogenides). Kept for backward compatibility and for the placement-stage workflows where no RDF is available yet.

Notes

Changed in v1.0.0: default switched from "auto" to "auto-rdf" to match neutron-diffraction analysis convention and to avoid systematically under-counting coordination for materials with broad bond-length distributions.

screen(config=True)[source]

Label candidates and record exclusion decisions without changing them.

Each screen has an independent exclude policy (default false). Coordination uses this analyser’s resolved chemical cutoffs; default allowed sets come from random generation’s inferred targets and tolerance. See amorphgen.analysis.screen_structures() for settings. The returned audit starts with zero analysed structures.

screened(config=True)[source]

Return (retained_analyser_or_None, audit) without mutating inputs.

Selection keeps the full candidate ensemble’s resolved cutoffs fixed and preserves file alignment. Labels alone do not remove a structure. Call mark_screening_analysed(audit, audit['retained_indices']) after successfully analysing the subset to record its contribution.

total_coordination(centre=None, partners=None)[source]

First-shell coordination of an element counting several partner types together.

Parameters:
  • centre (str, optional) – Element at the centre ("O"). Default: every element.

  • partners (iterable of str, optional) – Partner elements to count, such as In and Ga. Default: every bonded partner type (bonded = cation-anion or hetero covalent, the rule of the “Bonding coordination numbers” table). When given, the named partners are counted within their pair cutoffs whether or not the pair is classed as a bond.

Returns:

Keyed by element, with mean, std, min, max and distribution (CN to percent). For IGZO the default gives the total O coordination (Ga + In + Zn around O) that the per-pair O-Ga / O-In / O-Zn entries only give in parts; centre='O' with partners set to In and Ga gives the count over the two larger cations only.

Return type:

dict

density()[source]

Compute density (g/cm3) for each structure.

Returns:

{“values”: list[float], “mean”: float, “std”: float}

Return type:

dict

coordination(pair=None)[source]

Compute coordination numbers across all structures.

Parameters:

pair (str, optional) – Pair to analyse, e.g. “Si-O”. If None, all pairs computed.

Returns:

Keyed by pair string. Each value is a dict with “mean”, “std”, “min”, “max”, “distribution”, “total_atoms”.

Return type:

dict

cutoff_robustness(window=0.1, points=5)[source]

Report near-cutoff pair shares and directional coordination curves.

Sweep the resolved cutoffs over +/- window Angstrom (default 0.10) at points evenly spaced values (odd and >= 3, default 5). Auto cutoffs are held fixed, with no RDF refit. Positive cutoffs are shifted together, lower radii clipped to zero, and zero cutoffs kept at zero. The analyser and its structures are not modified.

pairs maps unordered element pairs to cutoff, actual cutoffs, near_pairs, n_pairs and near_fraction. The fraction divides contacts added across the window by contacts included at its upper endpoint; it is None for an empty denominator. Reciprocal contacts count once, distinct periodic images separately.

Each pair’s coordination has one curve per direction: mean pools central sites, ensemble_mean weights contributing structures equally, delta / ensemble_delta give upper-minus-lower changes, and per_structure retains aligned curves (None for absent centres). Pair counts also retain aligned per_structure observations. Zero-neighbour sites are included. The central point matches coordination(), including its exact distance-boundary rules. These perturbations measure sensitivity, not sampling uncertainty.

dimer_report(threshold_frac=0.85)[source]

Detect unphysical same-element / anion-anion close contacts.

Flags “dimers” — pairs closer than threshold_frac times the radii-derived minimum separation (O-O peroxide below ~1.9 A, Cl2, metal-metal dimers). These form when cold MLIP relaxation collapses under-coordinated seeds and mark structures that usually rank higher in energy; the recommended workflow is to generate an ensemble and keep the dimer-free, lowest-energy members.

Parameters:

threshold_frac (float) – Fraction of the pair minsep below which a contact counts as a dimer (default 0.85).

Returns:

{"pairs": {pair: {"count", "min_distance", "threshold"}}, "per_structure": list[int], "total": int, "n_structures": int, "threshold_frac": float}

Return type:

dict

bond_order(qbar6_threshold=0.3, min_neighbors=4, cutoff=None)[source]

Steinhardt q6, Lechner–Dellago qbar6 and ordered clusters.

Uses this analyser’s resolved neighbour cutoffs unless cutoff is provided. Ordered atoms have qbar6 at least qbar6_threshold and at least min_neighbors neighbour images. Calibrate the threshold against the material’s crystal and liquid structures. See amorphgen.analysis.compute_bond_order() for output fields.

bond_distances(pair=None)[source]

Compute bond distance statistics.

Parameters:

pair (str, optional) – Pair to analyse, e.g. “Si-O”. If None, all pairs computed.

Returns:

Keyed by pair string. Each value is a dict with “mean”, “std”, “min”, “max”, “count”.

Return type:

dict

bond_angles(triplet=None)[source]

Compute bond angle statistics.

Parameters:

triplet (str, optional) – Triplet to analyse, e.g. “O-Si-O”. If None, all triplets computed.

Returns:

Keyed by triplet string. Each value is a dict with “mean”, “std”, “min”, “max”, “count” (angles in degrees).

Return type:

dict

angle_distribution(triplet=None, bins=90, *, normalise=True, confidence=0.95, n_bootstrap=1000, seed=0)[source]

Per-structure angle histograms and uncertainty of their mean.

Integer bins spans 0–180 degrees, including linear angles. Explicit edges must cover that range. Each structure’s histogram is normalized independently before averaging; a structure without the triplet is missing, not a zero-density sample. Bands are pointwise percentile intervals from resampling complete structures.

rdf(pair=None, rmax=None, nbins=200, sigma=0.05, *, confidence=0.95, n_bootstrap=1000, seed=0)[source]

Compute the radial distribution function g(r).

Parameters:
  • pair (str, optional) – Pair to analyse, e.g. “Si-O”. If None, total RDF.

  • rmax (float, optional) – Maximum radius in A. Auto-detected from cell if None.

  • nbins (int) – Number of histogram bins (default 200).

  • sigma (float) – Gaussian smearing width in A (default DEFAULT_SMEARING = 0.05). Pass 0.0 for the raw histogram. S(q) methods are unaffected (they always use the raw g(r)).

Returns:

{“r”: list[float], “g_r”: list[float]}

Return type:

dict

structure_factor(pair=None, qmax=15.0, nq=300, rmax=None, weighting='unweighted', *, confidence=0.95, n_bootstrap=1000, seed=0)[source]

Compute the structure factor S(q) from g(r) via Fourier transform.

Parameters:
  • pair (str, optional) – Pair to analyse, e.g. "Ga-O". If None, computes the total S(q). The weighting argument only affects the total case; specific partials are always returned as their own Faber-Ziman S_AB(q).

  • qmax (float) – Maximum q in inverse-Angstrom (default 15.0).

  • nq (int) – Number of q points (default 300).

  • rmax (float, optional) – Max radius for the underlying g(r). None = auto (half cell).

  • weighting ({“unweighted”, “xray”, “neutron”}, default "unweighted") – How partials are combined into the total. "xray" uses q-dependent Waasmaier-Kirfel form factors (Faber-Ziman); "neutron" uses tabulated coherent scattering lengths. Match the weighting and normalization of the reference data.

Returns:

{"q": list[float], "s_q": list[float]}.

Return type:

dict

structure_factor_direct(qmax=15.0, nq=300, weighting='xray', sigma_q=0.0, partials=False, *, confidence=0.95, n_bootstrap=1000, seed=0)[source]

Compute Faber-Ziman S(q) from reciprocal-lattice scattering sums.

Avoids the finite-rmax integral in structure_factor(). The cell size still limits q sampling (q_min ~ 2*pi/L for a cubic cell); inspect n_per_bin and check convergence.

Parameters:
  • qmax (float) – Maximum q in inverse-Angstrom.

  • nq (int) – Number of q-bins for spherical averaging.

  • weighting ({“xray”, “neutron”, “unweighted”}, default "xray") – Per-element scattering factors (q-dependent Waasmaier-Kirfel form factors for x-rays, Sears scattering lengths for neutrons). See compute_structure_factor() for details.

  • sigma_q (float, default 0.0) – Gaussian re-binning width in 1/A, weighted by the number of q-vectors per shell. Reduces noise but can broaden features; compare with the raw values. 0 returns raw shell averages.

  • partials (bool, default False) – Also return the Faber-Ziman partial structure factors S_ab(q) computed from the per-species amplitudes.

Returns:

{"q": list[float], "s_q": list[float], "n_per_bin": list[int]}, plus "s_q_raw" when sigma_q > 0 and "partials" ({"A-B": list[float], ...}) when partials=True.

Return type:

dict

total_correlation(weighting='xray', qmin=0.3, qmax=20.0, nq=400, rmax=10.0, nr=600, window='lorch', *, confidence=0.95, n_bootstrap=1000, seed=0, sigma_q=0.05)[source]

Total correlation function T(r) = 4 pi r rho g(r), the curve a diffraction paper plots beside S(Q).

The g(r) behind it is SCATTERING-WEIGHTED and obtained by Fourier transforming the weighted S(Q) over the measured Q range, so it is directly comparable with published data and is NOT the same as rdf() with pair=None, which weights every pair equally. Set qmin/qmax/window to the experiment’s own values. sigma_q is the direct S(q) re-binning width in inverse Angstrom (default 0.05); use 0 for the unsmoothed shells.

Returns a dict with r, g_r, T_r, G_r (the reduced PDF), the q/s_q used, and rho.

xrd_pattern(wavelength=1.5406, qmax=None, nq=300, method='direct', rmax=None, sigma_q=0.0, *, q_batch=4096, confidence=0.95, n_bootstrap=1000, seed=0)[source]

Coherent X-ray intensity per atom versus 2theta in degrees.

wavelength is in Angstrom (default Cu K-alpha). Each structure’s X-ray S(q) is converted using its own composition before averaging. Returns q, two_theta, intensity, per_structure and pointwise uncertainty. Missing reciprocal shells stay missing. This profile includes no instrument, background or sample corrections. See amorphgen.analysis.compute_xrd_pattern() for conventions.

compare_experiment(experiment, kind=None, *, method='direct', weighting='xray', calculation_options=None, load_options=None, x_range=None, confidence=0.95, n_bootstrap=1000, seed=0)[source]

Compare a measured S(q) or T(r) file with this ensemble.

experiment is a path or a mapping returned by load_experiment. kind is "sq" (the default for paths) or "tr". Set calculation_options to the chosen scattering method’s options, e.g. {"qmax": 12, "sigma_q": 0.05}; for T(r), match the measured q range and window. load_options accepts delimiter, columns and skiprows. Coordinates use inverse Angstrom for S(q), Angstrom for T(r).

Each structure is interpolated onto the measured grid, with no extrapolation or interpolation across missing bins. The result has residuals, goodness-of-fit metrics and pointwise mean confidence bands. No scale, offset or other parameter is fitted. Scattering conventions and instrument resolution must already match the measured curve.

averaged_rdf(pair=None, rmax=None, nbins=200, *, confidence=0.95, n_bootstrap=1000, seed=0)[source]

Compute RDF per structure with mean and standard deviation.

Parameters:
  • pair (str, optional) – Pair to analyse, e.g. “Si-O”. If None, total RDF.

  • rmax (float, optional) – Maximum radius in A. Auto-detected from cell if None.

  • nbins (int) – Number of histogram bins (default 200).

Returns:

{“r”: list, “g_r_mean”: list, “g_r_std”: list, “n_structures”: int}

Return type:

dict

ring_statistics(bond_pair=None, cutoff=None, max_ring=12)[source]

Compute ring size statistics for network-forming structures.

Parameters:
  • bond_pair (tuple of str, optional) – Bond pair to trace, e.g. (“Si”, “O”). Auto-detected if None.

  • cutoff (float or dict, optional) – Positive bond cutoff in Angstrom, or an ASE pair-cutoff mapping. Uses analyser pair cutoff if None.

  • max_ring (int) – Maximum ring size to search (default 12).

Returns:

Size distribution, mean/spread, resolved and unresolved network edge counts, and per-structure statistics with uncertainty. Counts (including legacy total_rings) are shortest-cycle observations per edge, not the number of unique rings. Unresolved edges have no closure within max_ring. See amorphgen.analysis.rings.compute_ring_statistics().

Return type:

dict

polyhedral_connectivity(cation=None, anion=None)[source]

Corner/edge/face sharing between cation-centred polyhedra.

Two cations sharing one anion neighbour are corner-sharing, two anions edge-sharing, three or more face-sharing. Returns the linkage distribution, the percentage of cations in at least one edge- or face-sharing pair (near 0 in a-SiO2 / a-GeO2; tens of percent in unrelaxed random placements), per-species means and per-structure values. See amorphgen.analysis.structure.compute_polyhedral_connectivity().

voronoi(element=None)[source]

Compute Voronoi tessellation indices (n3, n4, n5, n6).

Parameters:

element (str, optional) – Element to analyse. If None, all atoms included.

Returns:

{“indices”: list, “distribution”: dict, “top_10”: list,

”mean_faces”: float, “total_atoms”: int}

Return type:

dict

void_distribution(n_samples=10000, probe_radius=0.0, radii=None, nbins=50, seed=0, *, probe_radii=None)[source]

Sample periodic free-space clearance and accessible volume.

Distances/radii are in Angstrom. This is a volume-weighted point clearance distribution; it does not identify connected pores. probe_radii selects thresholds for an accessible-volume curve evaluated from the same samples (default: histogram bin edges). Clearance quantiles describe points accessible to probe_radius. See amorphgen.analysis.voids.compute_void_distribution().

oxygen_speciation(network_formers=None)[source]

Count free/non-bridging/bridging and multiply shared oxygen.

Uses this analyser’s pair cutoffs. Defaults to Al, B, Ge, P and Si present; explicitly select the network formers for other oxides. See amorphgen.analysis.oxygen.compute_oxygen_speciation().

elastic_moduli(calculator=None, strain=0.005, relax=False, fmax=0.01, steps=200)[source]

Finite-strain stress response and Voigt/Reuss/Hill moduli (GPa).

Requires an active ASE stress calculator, supplied explicitly or attached to each structure. relax=True relaxes internal atomic positions at fixed cell. See amorphgen.analysis.elasticity.compute_elastic_moduli().

vibrational_dos(calculator=None, displacement=0.01, npoints=400, sigma=0.1)[source]

Harmonic DOS from finite-difference forces, in THz.

Uses Gamma-point normal modes of each supplied cell, with negative frequencies representing imaginary modes. Requires 6N force calls per structure. See amorphgen.analysis.vibrations.compute_vibrational_dos().

energy_ranking()[source]

Rank structures by potential energy per atom.

Reads energy from atoms.info or attached calculator.

Returns:

{“energies_per_atom”: dict, “ranking”: list, “best”: int,

”worst”: int, “best_energy”: float, “worst_energy”: float, “spread”: float}. Returns “warning” key if no energy data.

Return type:

dict

averaged_cn(pair=None)[source]

Compute per-structure mean CN with overall mean and std.

Parameters:

pair (str, optional) – Pair to analyse, e.g. “Si-O”. If None, all pairs computed.

Returns:

Keyed by pair string. Each value is a dict with “mean_per_structure”, “overall_mean”, “overall_std”, “n_structures”.

Return type:

dict

convergence_report(tolerances=None, *, descriptors=None, confidence=0.95, sizes=None, max_structures=1000000)[source]

Report ensemble precision against declared absolute tolerances.

Core descriptor names are density (g/cm3), coordination.Si-O and total_coordination.Si (neighbours), bond_distance.O-Si (angstrom), and bond_angle.O-Si-O (degrees), with the species present in this ensemble replacing these examples. Each observation is one structure’s mean, regardless of its atom count.

tolerances maps descriptor names to positive absolute confidence half-widths. Undeclared descriptors are reported without a pass/fail decision. descriptors can add named per-structure values or existing uncertainty summaries, for example {"rdf.total": self.rdf()["uncertainty"]}. Additional observations must align with this ensemble and may not replace core descriptors.

The curve uses full-ensemble variance and Student-t multipliers; reordering the same observations cannot change it. Projections assume independent structures with unchanged variance and descriptor availability. Curve-valued descriptors require every point’s half-width to meet the tolerance; intervals remain pointwise, not simultaneous. See amorphgen.analysis.convergence_report() for report fields, missing-data handling and the distinction from generation-prefix tests.

summary(show_angles=True, cutoff_window=0.1)[source]

Print and return structural statistics and cutoff sensitivity.

cutoff_window is the positive half-window in Angstrom used by cutoff_robustness() (default 0.10).

per_structure_summary(cutoff_window=0.1)[source]

Analyse each structure individually and produce a comparison table.

Appends ensemble cutoff sensitivity over +/- cutoff_window A.

Return type:

str — formatted table (also printed).

save_report(filepath, text=None, show_angles=True, cutoff_window=0.1)[source]

Save the summary report to a text file.

Parameters:
  • filepath (str) – Output file path.

  • text (str, optional) – Pre-computed report text. If None, calls summary().

  • show_angles (bool) – Include bond angles in the report (default True).

  • cutoff_window (float) – Half-window in Angstrom for cutoff robustness (default 0.10). Used only when generating report text here.

plot(**kwargs)[source]

Generate and save analysis plots (RDF, CN distribution, bond angles).

Parameters:
  • output_dir (str) – Directory for output files (default “.”).

  • prefix (str) – Filename prefix (default “analysis”).

  • rdf_pairs (list of str, optional) – Pairs to plot, e.g. [“Si-O”, “O-O”]. If None, all pairs.

  • rmax (float, optional) – Maximum radius for RDF. Auto-detected if None.

  • save_csv (bool) – Also save raw data as CSV files (default True).

Screening

StructureAnalyser.screen(config) returns per-candidate screening decisions without changing the analyser. StructureAnalyser.screened(config) returns (retained_analyser, report); the analyser is None when all candidates are excluded. The retained analyser uses the same resolved neighbour cutoffs as the original full ensemble. crystal_like.cutoff can independently override the order shell and is resolved once using that full ensemble. Neither method runs the subsequent analysis or marks structures as analysed.

from amorphgen.analysis import (
    StructureAnalyser, format_screening_report, mark_screening_analysed,
    write_screening_outputs,
)

sa = StructureAnalyser("ensemble/", cutoff="auto-rdf")
config = {
    "coordination": {"allowed": {"Si": [4], "O": [2]}},
    "close_contacts": {"threshold_frac": 0.7, "exclude": True},
    "unconverged": {"exclude": False},
}
retained, report = sa.screened(config)
if retained is not None:
    print(retained.summary())
    mark_screening_analysed(report, report["retained_indices"])
print(format_screening_report(report))
write_screening_outputs(report, "analysis/screening")

Each screen records its label and exclusion decision separately. exclude defaults to False; an unavailable metric yields an explicit unavailable label and follows that same exclusion policy. See Analysis for the six screens, bounds and thresholds, automatic coordination sets, and the generated / passed / labelled / analysed counts. The configuration is the screening mapping itself, not its enclosing analysis.screening YAML keys. True enables the default label-only screens.

The standalone screen_structures(atoms_list, config, *, cutoff="auto-rdf", source_names=None) accepts ASE structures directly. Reports initially have analysed=0; call mark_screening_analysed(report, indices) only after those original candidate indices have completed analysis. write_screening_outputs writes <prefix>.json, <prefix>_structures.csv and <prefix>_summary.csv. The CLI manages this bookkeeping and exports automatically.

amorphgen.analysis.validate_screening_config(config)[source]

Validate and normalize a screening mapping without mutating it.

True enables coordination, crystal-like order, close contacts and relaxation diagnostics. Density and energy need explicit min and/or max bounds. False, None and {} enable no screens. Individual screens accept a mapping or a boolean; omitted exclude is always false. Unknown keys and invalid thresholds raise ValueError before evaluation.

amorphgen.analysis.screen_structures(atoms_list, config, *, cutoff='auto-rdf', source_names=None)[source]

Evaluate configured screens and return JSON-safe records in input order.

Coordination uses chemically bonding neighbour-image counts and accepts auto_target_cn targets plus/minus its tolerance by default. Crystal-like order uses the existing neighbour-averaged qbar6 diagnostic. Close contacts include all pair types and periodic self-images; thresholds are strict lower bounds. Density (g/cm3) and energy (eV/atom) bounds are inclusive. Relaxation status comes only from explicit relaxation_converged or relax_converged booleans. Calculators are never evaluated.

Every nonpassing screen adds a label independently of its exclude flag. Missing measurements add <screen>_unavailable and follow that same flag. passed means no labels, retained_indices means not excluded, and analysed starts false until mark_screening_analysed() is called. Resolved neighbour cutoffs are shared over the input ensemble. An optional crystal_like.cutoff specifies its own neighbour definition, also resolved over the complete input ensemble before selection.

amorphgen.analysis.mark_screening_analysed(result, indices)[source]

Mark retained original input indices analysed after successful analysis.

amorphgen.analysis.format_screening_report(result)[source]

Render population accounting; columns overlap rather than form a funnel.

amorphgen.analysis.write_screening_outputs(result, prefix)[source]

Write .json, _structures.csv and _summary.csv; return paths.

Ring sizes and void clearance

from amorphgen.analysis import StructureAnalyser
from amorphgen.analysis.descriptors import save_descriptor

sa = StructureAnalyser("structures/", cutoff={"Si-O": 2.0})
rings = sa.ring_statistics(bond_pair=("Si", "O"), max_ring=16)
voids = sa.void_distribution(n_samples=20000, probe_radius=0.5, seed=42,
                             probe_radii=[0, 0.25, 0.5, 0.75, 1.0])
save_descriptor("rings", rings, "analysis/")
save_descriptor("voids", voids, "analysis/")

Ring counts and legacy total_rings count shortest-cycle observations per network edge, not unique cycles. mean_ring_size, std_ring_size (population spread), min_ring_size and max_ring_size summarize the resolved edges. n_network_edges, n_ring_edges, n_unresolved_edges and ring_edge_fraction report search coverage. An unresolved edge may close beyond max_ring; undefined size statistics and coverage are None. The resolved cutoff, counting convention and per-structure observations accompany the result.

Void probe_curve contains sorted unique radii and the aligned accessible_fraction, accessible_volume and corresponding *_stderr arrays. All thresholds use the same samples; these errors quantify Monte Carlo noise. Omitting probe_radii uses the histogram bin edges. The curve can include radii below the base probe_radius, independently of the histogram. clearance_quantiles gives empirical p10/p50/p90 clearances conditional on the base probe, using cell-volume weights across structures. No accessible samples gives None quantiles. Clearance describes local free space, not connected pores or maximal cavities.

Both results include per_structure observations and separate uncertainty summaries of equal-weight structure means. See Analysis for normalization, interpretation and exported files.

Cutoff robustness

StructureAnalyser.cutoff_robustness(window=0.1, points=5) measures contact and coordination sensitivity around the analyser’s resolved pair cutoffs. window is a finite positive half-width in Å; points is an odd integer of at least three. Automatic cutoffs are resolved once and frozen while the same offset is applied to each positive pair cutoff. Values are clipped at zero and zero cutoffs stay zero throughout the sweep.

The report includes pooled undirected periodic contact counts and the share whose inclusion changes from the lower to the upper endpoint, using the upper endpoint contact count as denominator. No contacts gives an undefined share. Coordination is directional, with pooled central-site means, per-structure means and equal-weight structure means retained. The ordinary distance <= pair cutoff and distance < largest cutoff boundary rules apply at every point. See Analysis for interpretation.

summary(show_angles=True, cutoff_window=0.1) and per_structure_summary(cutoff_window=0.1) include the five-point report by default. plot(..., cutoff_window=0.1) exports it when save_csv=True.

from amorphgen.analysis import (
    StructureAnalyser, format_cutoff_robustness, save_cutoff_robustness,
)

sa = StructureAnalyser("structures/", cutoff="auto-rdf")
report = sa.cutoff_robustness(window=0.15, points=7)
print(format_cutoff_robustness(report))
paths = save_cutoff_robustness(report, output_dir="analysis/", prefix="analysis")

save_cutoff_robustness(report, output_dir=".", prefix="analysis") writes analysis_cutoff_robustness.json, analysis_cutoff_robustness_pairs.csv and analysis_cutoff_robustness_coordination.csv with the default prefix.

amorphgen.analysis.format_cutoff_robustness(report)[source]

Format pair shares and directional coordination across the window.

amorphgen.analysis.save_cutoff_robustness(report, output_dir='.', prefix='analysis')[source]

Save the full JSON report and pooled pair/coordination CSV tables.

Returns a mapping of json, pairs and coordination to paths. Aligned per-structure contact counts and CN curves are retained in JSON. Fractions use the 0..1 scale; missing fractions are empty CSV cells.

Measured scattering and XRD

StructureAnalyser.compare_experiment() loads measured S(q) or T(r), calculates the matching ensemble curve and reports residuals and metrics with pointwise ensemble bands. The separate loader and comparison functions accept already calculated curves. StructureAnalyser.xrd_pattern() returns coherent X-ray intensity per atom on a physical 2θ axis. See Analysis for file formats, metric definitions, CLI examples and intensity conventions.

amorphgen.analysis.load_experiment(path, kind='sq', *, columns=None, delimiter=None, skiprows=0, comments='#')[source]

Load a measured S(q) or T(r) curve from a numeric text or CSV file.

Default input has exactly two columns (coordinate, observation), or three (coordinate, observation, positive one-sigma measurement uncertainty). columns=(x_column, y_column[, sigma_column]) selects zero-based columns explicitly when other columns are present. delimiter=None detects a comma in the first data line, otherwise whitespace is used. Headers must be commented or skipped explicitly using skiprows; malformed selected data are errors, never silently discarded. Blank/comment lines are ignored.

Data are sorted by coordinate with observations and uncertainties kept together. Duplicate or nonfinite coordinates, nonfinite observations, and nonpositive/nonfinite uncertainties are rejected. No units, normalisation, density, or S(q)/T(r) convention conversion is performed: q is in inverse angstrom, r in angstrom, S(q) dimensionless, and T(r)=4*pi*r*rho*g(r) in inverse square angstrom. The return value is JSON-native.

Return type:

dict

amorphgen.analysis.compare_experiment(calculated, experiment, *, kind=None, x_range=None, confidence=0.95, n_bootstrap=1000, seed=0)[source]

Compare an ensemble S(q) or T(r) result to measured observations.

Interpolate each structure onto the measured coordinates, then recompute the equal-weight ensemble mean and pointwise Student t/bootstrap intervals. Interpolation never bridges missing bins or extrapolates. At each point, only structures with supported finite values contribute. Measurements with no model support, or outside optional inclusive x_range=(lo, hi), are excluded and counted. Indices refer to the full sorted measurement grid. plot_break_before marks breaks needed before a retained point, either because measured points were excluded or because the measured grid skips an unsupported interval of the native model grid entirely. Plotting lines and bands must respect these breaks.

Residuals and bias are calculated minus observed. RMSE, MAE, and bias use equal point weights. Rw = sqrt(sum(w*residual**2)/sum(w*observed**2)), with w=1/sigma**2 if measured one-sigma uncertainties exist, otherwise w=1. Rw is a fraction, not a percentage, and is undefined for zero denominator. Its definition follows https://dictionary.iucr.org/R_factor.

Chi-square = sum((residual/sigma)**2) is provided only with measurement sigma; reduced chi-square divides by N (no fitted parameters). These are descriptive diagonal-error metrics: transformed T(r) and sampled curves can have correlated errors. They do not establish statistical significance or a p-value. Ensemble SEM is never substituted for measurement sigma. All metrics share exactly the returned grid. No fit or unit conversion is performed; experimental and model conventions must already agree.

Return type:

dict

amorphgen.analysis.format_experiment_report(result)[source]

Return a concise readable comparison report, including exclusions.

Return type:

str

amorphgen.analysis.save_experiment_comparison(result, output_dir='.', prefix='analysis', dpi=300, save_pdf=False)[source]

Save an S(q) or T(r) comparison as JSON, CSV, text and PNG/PDF.

result is returned by compare_experiment. CSV rows share the exact support used for the fit statistics and include source indices, counts, spread and confidence bounds. The plot separates measurement uncertainties from the ensemble mean confidence interval, with a residual panel below. Returns paths under keys json, csv, report, png and (when requested) pdf. Figures remain outside pyplot’s global registry.

amorphgen.analysis.compute_xrd_pattern(atoms_list, wavelength=1.5406, qmax=None, nq=300, method='direct', rmax=None, sigma_q=0.0, *, q_batch=4096, confidence=0.95, n_bootstrap=1000, seed=0)[source]

Calculate a coherent X-ray intensity profile per atom with ensemble bands.

wavelength is in Angstrom (default Cu K-alpha, 1.5406 A). The returned q grid is in inverse Angstrom, two_theta is in degrees, and intensity is in electron squared per atom. This is an equal-structure mean, without peak-height normalization. For each structure independently,

I_coh(q)/N = <f(q)>**2 * (S(q) - 1) + <f(q)**2>

uses that structure’s composition and neutral-atom Waasmaier-Kirfel form factors. S(q) follows the Faber-Ziman convention. method='direct' uses reciprocal-lattice shell means; 'ft' transforms the raw RDF on a shared rmax grid. Direct conversion evaluates form factors at shell centres, an approximation that improves as q bins narrow. FT termination ripples, including any negative intensities, are retained.

qmax=None selects min(15, 4*pi/wavelength). An explicit qmax must be physically accessible, and at most 24*pi inverse Angstrom (the form-factor fit’s validity limit). nq is the number of q bins/points, not uniformly spaced angular points. The angular mapping is 2*asin(q*wavelength/4/pi); intensities are sampled values, so no integration Jacobian is applied. rmax applies only to method='ft'. All input structures must be nonempty, fully periodic, and have a finite nonsingular three-dimensional cell. Compositions, volumes and atom counts may differ between structures.

sigma_q optionally Gaussian-rebins each coherent intensity curve in q before ensemble averaging, with reciprocal-vector weights for the direct method and equal point weights for FT. This is presentation smoothing, not an instrument-resolution model. Missing direct shells remain missing; intensity_raw and per_structure_raw retain the unsmoothed values.

per_structure and uncertainty retain the independent structure curves and pointwise uncertainty of their mean. Whole structures, not reciprocal vectors, are bootstrapped. With fewer than two contributing structures, interval bounds are unavailable. Missing values use JSON null. Correlated trajectory frames require blocking before use. Ensemble bands exclude finite-cell, form-factor and experimental systematic errors.

The output includes coherent self-scattering, but no anomalous/incoherent scattering, polarization, Lorentz, absorption, background or instrument response. Apply only corrections appropriate to the actual measurement.

amorphgen.analysis.save_xrd_pattern(result, output_dir='.', prefix='analysis', dpi=300, save_pdf=False)[source]

Save coherent X-ray intensity vs 2θ as strict JSON, CSV and PNG/PDF.

result is returned by xrd_pattern. CSV includes q (Å⁻¹), 2θ (degrees), coherent intensity per atom (electron²), the number of contributing structures, and pointwise uncertainty. Missing reciprocal shells remain empty in the table and gaps in the plot. Return a mapping of json, csv, png and optional pdf paths.

Optional material descriptors

The following functions are also exported from amorphgen.analysis. Their StructureAnalyser counterparts are bond_order(), void_distribution(), oxygen_speciation(), elastic_moduli() and vibrational_dos(). All return ensemble and per-structure results; the geometry functions need no calculator. Elastic and vibrational calculations require a live ASE calculator that can evaluate strained or displaced configurations.

from ase.io import read
from amorphgen.analysis import compute_void_distribution, compute_oxygen_speciation

frames = read("silica.extxyz", index=":")
voids = compute_void_distribution(frames, n_samples=20000, seed=42)
oxygen = compute_oxygen_speciation(frames, network_formers=["Si"],
                                  cutoff={"Si-O": 2.0})

See Analysis for physical conventions, calculator cost, CLI examples and exports. amorphgen.analysis.descriptors.save_descriptor accepts a result with the name bond_order, voids, oxygen_speciation, elastic or vdos and writes JSON, CSV and PNG, plus PDF with save_pdf=True.

amorphgen.analysis.compute_bond_order(atoms_list, cutoff='auto-rdf', qbar6_threshold=0.3, min_neighbors=4)[source]

Compute q6, qbar6 and connected crystal-like clusters for each frame.

The Steinhardt vector is the mean of Y_6m over all neighbours inside the cutoff. Lechner–Dellago averaging then averages those complex vectors over the central atom and its neighbour images, before taking the rotational invariant sqrt(4*pi/13*sum(|q6m|^2)). Both definitions include periodic self-images and repeated images in small cells, using ASE’s full periodic neighbour list, including for triclinic cells. An isolated atom has q6 = qbar6 = 0.

An atom is labelled ordered when qbar6 >= qbar6_threshold and its neighbour count is at least min_neighbors. The threshold is a material- and cutoff-dependent diagnostic, not a universal crystal classifier: calibrate it against the relevant crystal and liquid. Ordered atoms sharing a neighbour-list edge form a cluster. Cluster sizes count unique atoms in the simulation cell, even when a periodic cluster connects to itself across the boundary.

Parameters:
  • atoms_list (iterable of ase.Atoms) – Structures in input order. Empty frames and an empty list produce zero counts, fractions and means.

  • cutoff (float, dict, or str) – Neighbour cutoff in Angstrom; accepts the same scalar, pair-table and automatic forms as StructureAnalyser. All numeric cutoffs must be positive and finite. Resolved cutoffs are shared across frames so their results can be compared.

  • qbar6_threshold (float) – Inclusive ordered-atom threshold, between 0 and 1 (default 0.3).

  • min_neighbors (int) – Minimum neighbour-image count to classify an atom (default 4).

Returns:

parameters stores the resolved cutoff, qbar6_threshold and min_neighbors. per_structure contains atom arrays q6, qbar6, neighbor_counts, ordered and cluster_ids (disordered atoms have cluster ID -1), plus index, n_atoms, ordered_count, ordered_fraction, largest_cluster_size, largest_cluster_fraction, q6_mean and qbar6_mean. Fractions use all atoms in that frame as the denominator. Top-level scalar summaries are arithmetic means over frames, including largest_cluster_size; n_structures is the count. uncertainty contains standard errors and intervals of these structure means, excluding empty frames. fraction_of_sites pools ordered atoms, while fraction_of_structures counts frames with at least one ordered atom, including empty frames in its denominator.

Return type:

dict

amorphgen.analysis.compute_void_distribution(atoms_list, n_samples=10000, probe_radius=0.0, radii=None, nbins=50, seed=0, *, probe_radii=None)[source]

Sample the volume-weighted point-clearance distribution in periodic cells.

At each uniformly sampled cell point x, the clearance is min_i(min_image_distance(x, atom_i) - radius_i) in angstrom. Points with clearance >= probe_radius admit the centre of a spherical probe. This geometric accessibility does not imply a connected path to a pore.

radius and bin_edges describe clearance, not diameter or clearance minus probe radius. probability_density is normalized over the accessible points (integral one when any are found), whereas bin_volume_fraction sums to accessible_fraction. Frames are weighted by cell volume. accessible_volume is the arithmetic mean accessible volume per frame in angstrom cubed, not their sum.

probe_curve evaluates the fraction and volume accessible to each spherical probe using the same unfiltered clearance draws, including radii below probe_radius. probe_radii is a nonempty 1D sequence of finite nonnegative numbers, sorted and deduplicated in the result; when omitted it defaults to bin_edges. Curve standard errors are pointwise Monte Carlo errors, with correlated estimates across radii. clearance_quantiles gives p10, p50 and p90 of clearance conditional on accessibility at probe_radius, or None if none were observed. These are inverse empirical-CDF quantiles (the smallest sampled value reaching the target cumulative weight), with each draw weighted by its cell volume; per-frame quantiles give all draws equal weight.

n_samples independent points are drawn per frame with a local NumPy random generator. seed=None requests nondeterministic sampling. Frames receive random draws in a canonical geometry order, so a fixed seed gives the same ensemble observations when the input order changes. Per-structure results are still returned in the caller’s input order. Existing *_stderr fields quantify Monte Carlo sampling only. uncertainty separately estimates standard errors and intervals of equal-weight structure means; their observed spread includes Monte Carlo noise and between-structure variation. The per-frame Wilson intervals remain informative when no accessible/occupied points were sampled. Sample maxima are lower bounds on the true maximum clearance; no maximal-cavity search is done.

Radii are in angstrom. A mapping overrides ASE covalent radii for the specified element symbols; unspecified elements keep the ASE defaults. All frames must contain atoms and have a finite full-rank 3D periodic cell. Returns only JSON-compatible objects; inputs are not modified.

amorphgen.analysis.compute_oxygen_speciation(atoms_list, network_formers=None, cutoff='auto-rdf')[source]

Classify O by the number of bonded network-former neighbours.

Zero, one, two, three and four-or-more neighbours are labelled free, non_bridging, bridging, tricluster and higher_coordinated. Legacy fractions and fraction_of_sites are in [0, 1] and pooled by oxygen count. fraction_of_structures counts structures containing any oxygen of each species. Uncertainty uses equal-weight per-structure site fractions; oxygen-free structures have missing site fractions. Periodic images are distinct neighbours, including in small unit cells.

network_formers is an iterable of element symbols (or a comma/plus separated string). By default, the Al, B, Ge, P and Si present in the ensemble are selected. Other oxides require an explicit selection; alkali/alkaline-earth modifiers must not be counted as network formers. cutoff accepts the same forms as StructureAnalyser. All frames must have the same element set; atom counts may differ.

This is a geometric connectivity descriptor, not an assignment of charge, bond order, or hydroxyl speciation. “Free” means no selected network-former neighbour, even if other atoms are nearby.

amorphgen.analysis.compute_elastic_moduli(atoms_list, calculator=None, strain=0.005, relax=False, fmax=0.01, steps=200)[source]

Compute static elastic descriptors using a stress-capable ASE calculator.

Parameters:
  • atoms_list (iterable of ase.Atoms) – Nonempty, fully periodic 3D structures. Input cells, positions, constraints, and calculator attachments are left unchanged.

  • calculator (ASE calculator, optional) – A live calculator (for example an MLIP). If omitted, use each input’s attached calculator. Stored single-point stresses cannot supply the response to strain and are rejected. Calculator result caches may be updated during evaluation.

  • strain (float) – Positive engineering strain amplitude below one. Small amplitudes (typically 0.001–0.01) approximate linear response; check convergence against this setting for quantitative use.

  • relax (bool) – If true, relax atoms with BFGS at fixed cell before evaluating the reference and each strained configuration. Existing atomic constraints are respected. Failure to converge raises an error.

  • fmax (float) – Maximum force tolerance in eV/Angstrom for atomic relaxation.

  • steps (int) – Maximum optimization steps per configuration.

Returns:

JSON-compatible per-structure tensors, residual stresses, stability diagnostics, and Voigt/Reuss/Hill moduli in GPa (Poisson ratio is dimensionless). Ensemble means and population standard deviations use equal structure weights; unavailable values are excluded with counts. uncertainty reports standard errors and intervals of ensemble means, preserving missing per-structure values. Tensor components are flattened in row-major Voigt order, with shape metadata.

Return type:

dict

Notes

Columns differentiate stress against engineering strains in ASE Voigt order xx, yy, zz, yz, xz, xy. Shear deformation matrix entries are half the engineering strain. Each column uses independent positive and negative deformations of the reference. The raw tensor is symmetrized before computing moduli. There are 13 stress evaluations per structure.

These are zero-temperature, static, tangent stress-strain coefficients. No finite-pressure correction is applied. Relax the cell separately to near-zero stress for conventional equilibrium elastic moduli. Positive definiteness is a diagnostic of the symmetrized tensor, not a general finite-pressure stability criterion. Reuss/Hill estimates are unavailable for nonpositive or ill-conditioned tensors. Ensemble component averages assume structures are expressed in comparable Cartesian orientations.

amorphgen.analysis.compute_vibrational_dos(atoms_list, calculator=None, displacement=0.01, npoints=400, sigma=0.1)[source]

Compute a Gaussian-broadened harmonic vibrational density of states.

Parameters:
  • atoms_list (sequence of ase.Atoms) – Preferably relaxed structures. Each unconstrained structure contributes all 3N modes, including rigid translations/rotations and instabilities. Constraints are rejected; remove them explicitly to analyse all atoms.

  • calculator (ASE calculator, optional) – Force calculator shared by the calculations. If omitted, use each structure’s attached calculator. Stored single-point forces cannot describe displaced configurations and are rejected. This supports an MLIP or any other ASE calculator that can recalculate forces.

  • displacement (float) – Central finite-difference displacement in angstrom. Requires 6N force evaluations per structure and dense 3N-by-3N matrix diagonalisation; this calculation is deliberately opt-in.

  • npoints (int) – Number of frequency-grid points (at least two). Use enough points to resolve the Gaussian width over the returned frequency range.

  • sigma (float) – Gaussian standard deviation in THz.

Returns:

frequencies_thz and dos give a signed-frequency DOS with unit integral. projected_dos partitions this density by element using squared mass-weighted eigenvector components; the projections sum to dos. mode_frequencies_thz holds all individual frequencies: negative values represent imaginary modes, not negative real frequencies. Structures are pooled with equal weight per mode. per_structure provides frequencies, instability counts and Hessian asymmetry diagnostics, plus individually normalized DOS curves on the shared frequency grid. uncertainty gives equal-weight structure means, standard errors, t intervals and pointwise bootstrap bands for DOS and its element projections. Numerically zero eigenvalues are clipped only within roundoff (relative to the largest eigenvalue and matrix size).

Return type:

dict

Notes

This is a Gamma-point spectrum of each supplied cell, not a sampled Brillouin-zone phonon spectrum. No relaxation or acoustic sum rule is applied. A stationary reference structure is the caller’s responsibility. Original structures and their calculator attachments are preserved; calculators may update their ordinary internal result caches.

Bond order result

from amorphgen.analysis import StructureAnalyser

sa = StructureAnalyser("structures/", cutoff="auto-rdf")
order = sa.bond_order(qbar6_threshold=0.3, min_neighbors=4,
                      cutoff=3.5)  # example order shell; choose for the material
frame = order["per_structure"][0]
print(frame["ordered_fraction"], frame["largest_cluster_size"])
print(frame["q6"], frame["qbar6"], frame["cluster_ids"])

bond_order(..., cutoff=None) uses the analyser’s resolved cutoff unless overridden. The standalone compute_bond_order(atoms_list, cutoff="auto-rdf", qbar6_threshold=0.3, min_neighbors=4) accepts ASE frames directly. Neither method loads a calculator.

parameters records the resolved cutoff, qbar6_threshold and min_neighbors. Each per_structure item contains index, n_atoms, per-atom arrays q6, qbar6, neighbor_counts, ordered, cluster_ids, and scalar ordered_count, ordered_fraction, largest_cluster_size, largest_cluster_fraction, q6_mean and qbar6_mean. Cluster IDs are -1 for disordered atoms. Fractions divide by all atoms in the structure; cluster sizes count unique cell atoms. Top-level scalar summaries are arithmetic means across structures, including largest_cluster_size.

See Analysis for the order equations and threshold calibration, and MQ-ensemble workflow for the automatic initial-crystal retention report.

Initial-crystal retention

compute_melt_memory(initial_atoms, frames, cutoff="auto-rdf", qbar6_threshold=0.3, min_neighbors=4) takes an original ASE structure and a mapping from checkpoint labels to ASE frames. Its initial and comparisons entries contain order summaries; survival_fraction divides the retained initially ordered atom count by the original ordered count. parameters stores the fixed cutoff resolved from the input.

No initially ordered atoms gives an unavailable (None) survival fraction. Atom-count or element-sequence mismatches raise ValueError in this Python helper; the automatic MQ file report instead records unavailable checkpoint rows with reasons. Both require preserved atom indices and measure endpoint retention, so neither establishes uninterrupted crystal survival.

amorphgen.analysis.compute_melt_memory(initial_atoms, frames, cutoff='auto-rdf', qbar6_threshold=0.3, min_neighbors=4)[source]

Compare a mapping of named ASE frames to the original crystal.

Cutoffs are resolved using only initial_atoms and reused for every frame, which retains its own cell and PBC. survival_fraction is the number of initially ordered atom indices ordered at the endpoint divided by the initial ordered count. It is None if no initial atoms qualify. Count or element-sequence mismatches raise ValueError; this function never truncates or remaps atoms to manufacture correspondence.

amorphgen.analysis.format_melt_memory(report)[source]

Return a compact text report with the denominator and interpretation.

Ensemble convergence

StructureAnalyser.convergence_report(tolerances=None, *, descriptors=None, confidence=0.95, sizes=None, max_structures=1000000) gathers density, coordination, total coordination, bond-distance and bond-angle summaries. descriptors adds named aligned per-structure arrays or uncertainty summaries containing per_structure. Tolerances are positive absolute Student-t mean interval half-widths, with curve descriptors assessed at every point.

The standalone convergence_report accepts a mapping of descriptor names to scalar or structure-by-component observations, with missing entries as None or nonfinite numbers. The report includes descriptors, observed counts and half-widths, curve data, declared tolerance statuses and estimated total and additional structure counts. Arrays and missing values are JSON-compatible. Canonical reductions make the curves and forecasts independent of input order.

from amorphgen.analysis import (
    convergence_report, format_convergence_report, save_convergence_report,
)

report = convergence_report(
    {"density": [2.20, 2.25, 2.18, 2.23]}, {"density": 0.02},
    confidence=0.95, max_structures=10000,
)
print(format_convergence_report(report))
paths = save_convergence_report(report, "analysis/", save_pdf=True)

Complete-data planning curves are exact componentwise all-subset RMS Student-t half-widths for sizes 2 through the observed ensemble size; larger sizes extrapolate the variance model. Vector summaries take the maximum of those componentwise values. Missing data use an observed-availability approximation. These are conditional precision forecasts under independent sampling, not empirical histories or simultaneous confidence bands. See Declared tolerances and ensemble convergence for the formulas, CLI/YAML examples, status meanings and exported files.

amorphgen.analysis.convergence_report(descriptors, tolerances=None, confidence=0.95, sizes=None, max_structures=1000000)[source]

Report absolute precision tolerances and ensemble-size planning curves.

descriptors maps names to per-structure scalars, aligned structure-by- component matrices, or uncertainty summaries containing per_structure. Every descriptor must have the same number of structure rows. None, NaN and infinity are missing observations. Whole missing curve rows may be None. Empty curve components and fewer than two finite observations in any component cannot satisfy a tolerance or produce a sample-size forecast.

tolerances maps any subset of descriptor names to positive, finite, absolute CI half-widths in that descriptor’s own units. Vector descriptors pass only when every pointwise half-width meets their tolerance; this is not a simultaneous-coverage claim. Undeclared descriptors are reported but do not participate in the overall tolerance status. Omit tolerances for an exploratory report with no declared stopping criterion.

Curves use the full ensemble sample variance, canonically sorted to ensure exact invariance to row permutations. With complete data each component equals the RMS Student-t half-width over all subsets of size n, for 2 <= n <= observed ensemble size; larger sizes extrapolate the same model. Vector descriptors use the maximum componentwise RMS. With missing data, contributing counts are floor(n * observed_count / observed_total), a plug-in availability approximation. Forecasts keep both the observed variance and availability fixed and find the smallest total size at least as large as the current ensemble satisfying all declared bounds.

sizes optionally specifies positive integer planning sizes (including future sizes); duplicates are removed, sizes sorted, and the observed endpoint always included. By default, at most 64 observed sizes are used. max_structures caps sample-size forecasts, not the observed ensemble. Results contain only JSON-native values, with undefined numbers as null. This planning report is explicitly not sequentially valid: repeated looks do not preserve its nominal confidence coverage.

Return type:

dict

amorphgen.analysis.format_convergence_report(report)[source]

Return declared tolerances, observed uncertainty and sampling forecasts.

Uncertainty is a Student-t interval half-width for the ensemble mean. Curve-valued descriptors use their largest pointwise half-width; this is not a simultaneous confidence band. Forecasts are conditional on the observed variance and availability remaining representative.

amorphgen.analysis.save_convergence_report(report, output_dir, dpi=300, save_pdf=False)[source]

Save a report as text, strict JSON, two CSV tables and descriptor plots.

Files are named analysis_convergence*. The curves CSV has one row per descriptor, ensemble size and component, retaining unavailable values as empty cells. It includes dashed planning projections through each estimated target. The summary CSV stores component counts as a JSON scalar or list. Independent PNG figures (and optional PDFs) identify the observed endpoint, declared tolerance and projected required size. Return a mapping of paths.

Reference-validation helpers

For comparing computed metrics against literature ranges (used by amorphgen --analyse --reference REF.yaml):

Validate computed structural metrics against literature reference ranges.

A reference YAML lists expected ranges for density, bond distances, mean coordination numbers, and bond angle means. Each metric is compared to the analyser’s structure-weighted mean and its confidence interval. Intervals crossing a reference bound are “inconclusive”; otherwise a metric is labelled “match” / “concern” / “fail”. Legacy results without uncertainty metadata are still compared using their point estimates.

amorphgen.analysis.validate.validate_against_reference(analyser, reference)[source]

Compare analyser output to a reference dict (loaded from YAML).

Parameters:
  • analyser (StructureAnalyser)

  • reference (dict) – Parsed YAML with optional keys: density, bond_distances, coordination, bond_angles (see examples/reference_*.yaml).

Returns:

{"system": str, "sources": list[str], "rows": list[tuple], "intervals": dict} where each row is (descriptor, computed, expected_lo, expected_hi, units, verdict). A metric the structures do not have (an element pair absent, or with no contact inside its cutoff) keeps its row, with None and the verdict “n/a”. intervals maps descriptor names to the analyser’s uncertainty metadata. Means and t intervals use structures as independent sampling units, not pooled atoms. With uncertainty metadata but no finite interval (e.g. one structure), a finite descriptor is “inconclusive”.

Return type:

dict

amorphgen.analysis.validate.format_validation_report(result)[source]

Render the dict from validate_against_reference() as a printable table.

Energy ranking helpers

For parsing random_gen.log and ranking generated structures by total energy (used by amorphgen --rank-from-log LOG):

Energy ranking for multiple structures.

amorphgen.analysis.energy.compute_energy_ranking(atoms_list)[source]

Rank structures by potential energy.

Reads energy from atoms.info or calculator. Uncertainty of the mean is estimated over structures, with missing energies retained in input order.

amorphgen.analysis.energy.rank_from_log(logfile)[source]

Parse a random-gen log file and rank structures by total energy.

The relax loop in batch_random() prints final energy on the last optimizer step row. This function reads those rows directly, so energy ranking works for VASP/CIF outputs that don’t store energy.

Parameters:

logfile (str) – Path to random_gen.log.

Returns:

{"rows": [(idx, energy, e_per_atom, fmax, n_steps, status), ...] sorted by e_per_atom ascending, "n_atoms": int, "best": idx, "worst": idx, "spread_meV_per_atom": float}.

Return type:

dict

amorphgen.analysis.energy.format_log_ranking(result, logfile=None)[source]

Render the dict from rank_from_log() as a printable table.

Comparing ensembles

Use EnsembleSpec to describe each structure set and compare_ensembles to save RDF, coordination, bond-angle, and density comparisons.

class amorphgen.analysis.EnsembleSpec(label, files, color=None, cutoff='auto-rdf')[source]

Bases: object

Specification of one ensemble to include in a comparison.

Parameters:
  • label (str) – Display name (used in legend and axis labels).

  • files (list | str) – Either a list of structure-file paths, or a single glob string like "hybrid_runs/run_*/final_amorphous.xyz".

  • color (str | None) – Matplotlib colour. If None, one is drawn from DEFAULT_COLORS in registration order.

  • cutoff (str) – Cutoff mode for StructureAnalyser ("auto", "auto-rdf", or a numeric value). Defaults to "auto-rdf" (the first RDF minimum) — the correct neighbour cutoff for coordination counting. The plain "auto" mode can land near the bond peak and undercount CN (spurious CN 0/1), so it is not the default here.

label: str
files: list | str
color: str | None = None
cutoff: str = 'auto-rdf'
resolve_files()[source]

Expand a glob string to a sorted list of file paths.

Return type:

list[str]

analyser()[source]

Lazy-built StructureAnalyser for this ensemble.

Return type:

StructureAnalyser

classmethod from_analyser(label, analyser, color=None)[source]

Wrap an already-built StructureAnalyser.

Useful when you have an analyser object already (e.g., from StructureAnalyser.plot()) and want to reuse its loaded atoms/cutoff rather than re-reading files.

Parameters:
Return type:

EnsembleSpec

amorphgen.analysis.compare_ensembles(ensembles, rdf_pairs=None, cn_top_key=None, cn_bot_key=None, angle_keys=None, exp_density=None, output_dir='comparison_plots', prefix='', exp_label='Expt.', save_pdf=True)[source]

Run all four panel functions for the supplied ensembles.

Each descriptor produces three files: {prefix}_{descriptor}.png, {prefix}_{descriptor}.pdf, {prefix}_{descriptor}.csv.

Skip any descriptor by passing None for its key arguments.

Parameters:
Return type:

None

Submodule reference

StructureAnalyser delegates to focused submodules; advanced users can import these directly:

Submodule

Provides

analysis.rdf

Pair distribution function g(r), partial RDFs, S(q), T(r)

analysis.structure

Coordination numbers, bond distances, bond angles

analysis.bond_order

Steinhardt \(q_6\), Lechner–Dellago \(\bar q_6\) and periodic ordered clusters

analysis.melt_memory

Initial ordered-atom retention at MQ melt endpoints and snapshots

analysis.rings

Shortest-path ring statistics with periodic-image closure

analysis.voronoi

Voronoi cell volumes and connectivity

analysis.voids

Periodic Monte Carlo point clearance and accessible volume

analysis.oxygen

Oxygen classes by network-former coordination

analysis.elasticity

Stress-derived stiffness and Voigt/Reuss/Hill moduli

analysis.vibrations

Harmonic cell-mode DOS and element projections

analysis.descriptors

Optional descriptor summaries and JSON/CSV/figure export

analysis.convergence

Order-independent uncertainty curves, tolerances and sample-size planning

analysis.convergence_output

Convergence text, JSON, CSV and per-descriptor figures

analysis.energy

Total-energy parsing and ranking

analysis.cutoff

Bond-cutoff selection from g(r) first minimum

analysis.plotting

Publication-quality matplotlib helpers

analysis.validate

Reference-YAML validation

analysis.comparison_plots

Multi-ensemble comparison plots and CSV output