Skip to content

Latest commit

 

History

History
530 lines (412 loc) · 18.4 KB

File metadata and controls

530 lines (412 loc) · 18.4 KB

Command-line interface

The structure-factor command computes a connected two-dimensional structure factor from particle coordinates in a periodic rectangular box. It exposes all estimators in the Python library and can:

  • print a short calculation summary;
  • save the numerical arrays in a compressed NumPy archive;
  • save publication-style PNG figures;
  • display the figures interactively.

The command performs one calculation per invocation. It is therefore suitable both for exploratory use in a terminal and for shell scripts that sweep over files or parameters.

Installation and invocation

Install the project to make the console command available:

python3 -m pip install .

The general syntax is:

structure-factor POSITION_FILE Lx Ly [OPTIONS]

The equivalent module invocation is:

python3 -m structure_factor POSITION_FILE Lx Ly [OPTIONS]

From a source checkout where the package has not been installed, prepend the source directory to PYTHONPATH:

PYTHONPATH=src python3 -m structure_factor positions.dat 128 128 --no-show

Use the installed command in production scripts when possible. The PYTHONPATH=src form is mainly convenient while developing the project.

To inspect the syntax implemented by the installed version:

structure-factor --help
structure-factor --version

Required positional arguments

structure-factor POSITION_FILE Lx Ly
Argument Meaning
POSITION_FILE Whitespace-delimited text file containing one x y coordinate pair per row.
Lx Positive finite periodic-box length along the x direction.
Ly Positive finite periodic-box length along the y direction.

All lengths in a command must use the same physical unit. For example, if the coordinates and box lengths are measured in particle diameters, --grid-size, --particle-diameter, and --kernel-length must use particle diameters too. Reported wave numbers are angular wave numbers in the inverse of that unit.

The box is half-open:

$$ 0 \leq x < L_x, \qquad 0 \leq y < L_y. $$

Coordinates exactly equal to Lx or Ly are invalid. The command does not silently wrap, clip, or otherwise repair coordinates outside the box.

Input-file format and validation

A valid position file has at least one row and exactly two numeric columns:

0.50 0.50
1.25 3.75
7.90 2.00

The loader uses NumPy's whitespace-delimited text reader. Integer notation, decimal notation, and scientific notation are accepted. Blank lines and lines beginning with # follow NumPy's usual text-file behavior. Commas are not field separators.

The command rejects a file if it:

  • cannot be opened or read;
  • is empty;
  • is not a numeric table;
  • has a number of columns other than two;
  • contains NaN, positive infinity, or negative infinity;
  • contains a coordinate outside the half-open box.

The loaded coordinates are converted to contiguous float64 values. Input order is preserved. Periodicity enters through the reciprocal modes and the grid geometry; it is not used to normalize invalid coordinates.

For reproducible analyses, write the exact box lengths alongside the data. They cannot be inferred unambiguously from the largest observed coordinates.

Complete option reference

Option Default Used by Meaning
--method {particle,binary,density,occupancy,kernel} particle all calculations Select the estimator.
--grid-size FLOAT 2.0 all grid methods Requested maximum cell spacing. The effective spacing can be slightly smaller so an integer number of cells tiles the box exactly.
--density-threshold FLOAT 0.7 binary Nonnegative local packing-fraction cutoff. A cell exactly at the cutoff is dense.
--particle-diameter FLOAT 1.0 binary Positive particle diameter used to convert occupancy to local packing fraction.
--max-mode INTEGER 8 particle Positive reciprocal-lattice cutoff in both directions.
--direct-chunk-size INTEGER 4096 particle Positive number of particles processed in each direct-sum chunk. It changes temporary memory use, not the result.
--kernel {gaussian,top-hat,disk} gaussian kernel Smoothing-window shape.
--kernel-length FLOAT --grid-size kernel Positive Gaussian standard deviation, square top-hat half-width, or disk radius.
--no-radial-average radial average enabled all methods Skip radial binning and omit the radial arrays and radial figure.
--output-data PATH no data file all methods Save arrays and selected metadata in a compressed .npz archive.
--output-dir DIRECTORY no figure directory all methods Save the figures available for the result as PNG files.
--no-show interactive display enabled plotting Do not open interactive plot windows. This does not prevent saving figures.
--version command Print the package version and exit successfully.
-h, --help command Print command help and exit successfully.

Options that do not apply to the selected estimator do not alter its numerical result. For example, --grid-size does not change a particle calculation, and --max-mode does not change an occupancy calculation. The command still constructs and validates one complete configuration, however, so explicitly supplied sizes and integer counts must be valid even when another method is selected.

Grid size

For binary, density, occupancy, and kernel, the requested maximum spacing $\Delta$ produces

$$ M_x = \left\lceil \frac{L_x}{\Delta} \right\rceil, \qquad M_y = \left\lceil \frac{L_y}{\Delta} \right\rceil. $$

The actual cell sizes are

$$ h_x = \frac{L_x}{M_x}, \qquad h_y = \frac{L_y}{M_y}. $$

Consequently, the --grid-size value is a resolution target, not necessarily the exact width of every cell. A smaller value creates more cells, extends the FFT to higher wave numbers, and uses more memory. Scientific results should be checked for convergence under grid refinement.

Particle mode cutoff and chunk size

The particle estimator samples

$$ \mathbf{k}_{mn}

\left( \frac{2\pi m}{L_x}, \frac{2\pi n}{L_y} \right), \qquad -M \leq m,n \leq M, $$

where M is --max-mode. The resulting spectrum has shape (2*M + 1, 2*M + 1). Its cost grows in proportion to N * (2*M + 1)^2, so increasing the mode cutoff can be substantially more expensive for large snapshots.

--direct-chunk-size bounds the temporary phase matrices used during this sum. Lower it when memory is constrained; raise it only if benchmarking shows that larger batches help on the current machine. Changing it does not change the wave-number grid, normalization, or mathematical result.

The CLI accepts one common integer cutoff. The Python API additionally permits different x and y cutoffs.

Binary-field parameters

For a cell occupancy $n_{ij}$, the binary estimator first computes the local packing fraction

$$ \phi_{ij}

\frac{n_{ij}\pi\sigma^2}{4h_xh_y}, $$

where $\sigma$ is --particle-diameter. It then assigns -1 when $\phi_{ij}$ is below --density-threshold and +1 when it is equal to or above the threshold.

The resulting spectrum describes phase morphology rather than microscopic particle-density fluctuations. Its absolute amplitude is therefore not directly comparable to the particle-, occupancy-, density-, or kernel-normalized spectra.

Kernel parameters

--kernel-length has a different geometric meaning for each kernel:

  • gaussian: standard deviation of a radial Gaussian;
  • top-hat: half-width of an axis-aligned square;
  • disk: radius of a circular window.

The discrete periodic kernel is normalized to conserve particle mass. No deconvolution is applied, so the result is a deliberately filtered spectrum. If the kernel length is smaller than an effective grid cell, the command emits a runtime warning because the sampled window can be under-resolved and may collapse to the occupancy result.

Estimator examples

The following examples assume positions.dat belongs to a 128 x 128 periodic box. --no-show makes each command safe on a headless machine.

Direct particle definition

structure-factor positions.dat 128 128 \
  --method particle \
  --max-mode 16 \
  --direct-chunk-size 4096 \
  --output-data results/particle.npz \
  --no-show

At each nonzero sampled reciprocal mode, this computes

$$ S(\mathbf{k})

\frac{1}{N} \left| \sum_{j=1}^{N} e^{-i\mathbf{k}\cdot\mathbf{r}_j} \right|^2. $$

It does not discretize particle positions and is the appropriate reference when a low-wave-number microscopic structure factor is required.

Occupancy FFT

structure-factor positions.dat 128 128 \
  --method occupancy \
  --grid-size 1.0 \
  --output-data results/occupancy.npz \
  --no-show

This bins particles into counts $n_{ij}$, subtracts the mean count, and divides the Fourier power by the particle number.

Number-density FFT

structure-factor positions.dat 128 128 \
  --method density \
  --grid-size 1.0 \
  --output-data results/density.npz \
  --no-show

This exposes $n_{ij}/(h_xh_y)$ as its real-space field. Its Fourier amplitude is multiplied by the cell area, so its dimensionless spectrum agrees with the occupancy spectrum up to floating-point round-off. The two methods differ in the units of source_field, not in their reported spectral normalization.

Binary coarse-grained field

structure-factor positions.dat 128 128 \
  --method binary \
  --grid-size 2.0 \
  --density-threshold 0.7 \
  --particle-diameter 1.0 \
  --output-data results/binary.npz \
  --output-dir figures/binary \
  --no-show

This is useful for domain morphology and phase-separation analyses where a thresholded dilute/dense field is the observable of interest.

Gaussian kernel

structure-factor positions.dat 128 128 \
  --method kernel \
  --grid-size 1.0 \
  --kernel gaussian \
  --kernel-length 2.0 \
  --output-data results/kernel-gaussian.npz \
  --no-show

Square top-hat kernel

structure-factor positions.dat 128 128 \
  --method kernel \
  --grid-size 1.0 \
  --kernel top-hat \
  --kernel-length 2.0 \
  --output-data results/kernel-top-hat.npz \
  --no-show

Disk kernel

structure-factor positions.dat 128 128 \
  --method kernel \
  --grid-size 1.0 \
  --kernel disk \
  --kernel-length 2.0 \
  --output-data results/kernel-disk.npz \
  --no-show

For all three kernels, the source field contains smoothed particle mass per cell and sums to the particle count to floating-point precision.

Connected spectrum and radial average

Every method reports a connected spectrum. The mean field is removed, and the zero-wave-number cell is set exactly to zero. For the unconnected particle definition, the corresponding zero-mode value would be $N$; neither that value nor the connected zero alone estimates compressibility from one fixed-particle-number snapshot.

Unless --no-radial-average is present, the command also groups the actual wave-number magnitudes into circular bins. The default bin width is

$$ \min\left(\frac{2\pi}{L_x},\frac{2\pi}{L_y}\right). $$

The zero mode, unmatched negative Nyquist modes on even FFT grids, empty bins, and modes outside the largest symmetrically supported reciprocal-space circle are omitted. The saved radial arrays contain the mean actual wave number, mean power, and contributing mode count for every retained nonempty bin.

Numerical output

Use --output-data to write a compressed NumPy archive:

structure-factor positions.dat 128 128 \
  --method occupancy \
  --grid-size 1.0 \
  --output-data results/occupancy.npz \
  --no-show

Missing parent directories are created automatically. If the path does not end in .npz, NumPy appends that suffix. Text metadata is stored as NumPy Unicode arrays rather than pickled Python objects, so the file can be loaded with pickle disabled:

import numpy as np

with np.load("results/occupancy.npz", allow_pickle=False) as result:
    print(result.files)
    print(result["method"].item())
    print(result["spectrum"].shape)

NPZ schema

Two-dimensional arrays use x-first ordering: array[x_index, y_index]. In particular, spectrum[i, j] belongs to (kx[i], ky[j]).

The following keys are always present:

Key Shape Meaning
method scalar Estimator name as a NumPy Unicode value.
particle_count scalar Number of input particles, $N$.
box_size (2,) (Lx, Ly) in coordinate units.
grid_size scalar Requested target grid size. It is saved even for the particle method and is not the derived effective cell width.
spectrum (len(kx), len(ky)) Connected two-dimensional float64 structure factor.
kx (len(kx),) Shifted angular wave numbers for the first spectrum dimension.
ky (len(ky),) Shifted angular wave numbers for the second spectrum dimension.

The following keys are conditional:

Key Present when Meaning
radial_wave_numbers radial averaging enabled Mean actual wave-number magnitude in each retained bin.
radial_spectrum radial averaging enabled Mean spectral power in each retained bin.
radial_mode_counts radial averaging enabled Number of two-dimensional modes contributing to each bin.
source_field any grid method Real-space field transformed by the selected estimator.
occupancy any grid method Raw particle count in every grid cell.
max_mode particle Direct reciprocal-lattice cutoff.
kernel kernel Selected kernel name as a NumPy Unicode value.
kernel_length kernel Resolved smoothing length; this equals grid_size when no explicit length was supplied.

When --no-radial-average is used, all three radial keys are absent rather than being stored as null values. For the particle method, source_field and occupancy are absent. For grid methods, source_field, occupancy, and spectrum share the configured grid shape.

The archive stores the arrays needed for analysis and selected estimator-specific metadata; it is not a complete copy of every command-line token. In particular, the current schema does not store density_threshold, particle_diameter, or direct_chunk_size. Keep the command, a script, or a separate provenance record when exact reconstruction of a run matters.

Figure output and interactive display

The plotting behavior is controlled independently by --output-dir and --no-show:

Options Save PNG files Open interactive windows
neither no yes
--no-show no no
--output-dir DIR yes yes
--output-dir DIR --no-show yes no

The last form is recommended for servers, schedulers, containers, and automated pipelines:

MPLBACKEND=Agg structure-factor positions.dat 128 128 \
  --method binary \
  --grid-size 2.0 \
  --output-dir figures/binary \
  --no-show

The output directory and missing parents are created automatically. PNG files are written at 160 DPI with a tight bounding box:

Filename Present when Contents
spectrum.png always Two-dimensional connected spectrum on physical kx and ky axes, with positive power displayed logarithmically.
source_field.png grid methods Occupancy, density, binary phase field, or kernel-smoothed real-space field in the physical box.
radial_spectrum.png radial averaging enabled Log-log radial spectrum with the strongest displayed mode and its wavelength $2\pi/k_\mathrm{peak}$ marked.

The binary source-field plot uses categorical colors for -1 and +1. Continuous source fields use a quantitative color scale. The spectrum plot masks the connected DC mode and can include a low-wave-number inset when the available resolution makes it useful.

If --no-radial-average is supplied, radial_spectrum.png is not generated. The particle method has no gridded real-space field, so it does not generate source_field.png.

In batch environments, setting MPLBACKEND=Agg makes the noninteractive backend explicit. With --output-dir --no-show, figures are saved and then closed. With --no-show and no output directory, Matplotlib is not imported at all.

Standard output and exit behavior

After a successful calculation, the command prints:

Computed occupancy structure factor on 128 x 128 modes.
Radial spectrum contains 63 nonempty bins.

The dimensions and bin count depend on the chosen parameters. The radial line is omitted when --no-radial-average is used.

Normal calculations return exit status 0. --help and --version also exit successfully without running a calculation.

Syntax errors, unreadable input, invalid coordinates, invalid configuration, and caught output errors produce a concise message on standard error and exit status 2, for example:

structure-factor: error: max_mode must be a positive integer; got 0

This behavior makes failures detectable in shell automation:

if ! structure-factor positions.dat 128 128 --no-show; then
  echo "structure-factor calculation failed" >&2
  exit 1
fi

Practical recommendations

  • Use --no-show in scripts, even when no figures are requested, so a process never waits for a graphical window.
  • Combine --output-data with --output-dir when both machine-readable data and inspection figures are useful; one calculation produces both.
  • Quote paths containing spaces.
  • Start particle calculations with the default low mode cutoff and estimate the cost before increasing it.
  • Check grid estimators at more than one --grid-size. High-wave-number power is affected by grid assignment, and this release applies no assignment-window correction.
  • Treat occupancy and density as two real-space representations of the same normalized gridded spectrum.
  • Do not compare the absolute binary-field amplitude directly with particle-normalized spectra.
  • Treat kernel results as filtered spectra. Top-hat and disk transforms have zeros, so the release intentionally performs no unstable deconvolution.
  • Record the full command and software version with scientific results. The NPZ archive is convenient and pickle-free, but it does not contain every CLI parameter or the original particle coordinates.

For the estimator definitions and normalizations, see Scientific methods. For programmatic use, see the Python API guide. Runnable workflows for the bundled xy.dat snapshot are collected in Examples.