Densities

The classes under the density module handle the processing and plotting of charge and band densities. Note that there is no explicit dependence on any of the classes used in C2a, and therefore could be used standalone. However, it only supports orthorhombic cells.

class conquest2a.density.bandden(directory: str, hkl: tuple[int, int, int], offset: float, operations: str, band: int, spin: int | None = None, kpt: int | None = None)

Class which looks for band density files in a directory.

Parameters:

directory – Path to directory containing band density files

locate_banddens_files() list[str]

Gets the paths to all band density files, sorted by band number, \(k\) point then spin.

Raises:

FileNotFoundError – If self.directory does not exist.

Returns:

List of all paths to band density files.

Return type:

list[str]

class conquest2a.density.chden(directory: str, hkl: tuple[int, int, int], offset: float, operations: str, charge_stub: str = 'chden', spin: int | None = None)

Class which looks for charge density files in a directory.

Parameters:

directory – Path to directory containing band density files

locate_chden_files() list[str]

Gets the paths to all charge density files.

For spin-unpolarised calculation, CONQUEST will output just chden.cube. For spin-polarised calculation, CONQUEST will output chden_up.cube and chden_dn.cube.

This output may change if the user has set Process.ChargeStub in Conquest_input which defaults to chden. Since reading Conquest_input directly is unsupported, just assume the user will pass in any stub.

Raises:

FileNotFoundError – If self.directory does not exist.

Returns:

List of all paths to band density files.

Return type:

list[str]

class conquest2a.density.density(hkl: tuple[int, int, int], offset: float, paths: Sequence[str], operations: str | None = None)

Process volumetric density data

CONQUEST can output any number of .cube files, e.g.:
  • chden_up.cube/chden_dn.cube for a spin polarised calculation

  • a single cube file, e.g. chden.cube for unpolarised calculations

  • Band densities

This class uses ASE’s .cube file processor to get atoms and data. It is therefore suitable for generic .cube file processing.

An arbitrary number of .cube files can be supplied via paths. If more than one path is given, operations is a string describing how each subsequent file is combined with the running total, applied left to right, one character per file after the first, much like VESTA:

  • "+" add

  • "-" subtract

  • "x" multiply

  • "/" divide

Example: for paths=[a, b, c, d] and operations="+/-", the resulting data is ((a + b) / c) - d.

Parameters:
  • hkl (INT ARRAY) – The \(hkl\) slice of the crystal to plot charge densities in.

  • offset – The \(hkl\) direction defines a family of planes.

Use offset to select which one (i.e. wherein the unit cell). :type offset: float :param paths: Path(s) to one or more charge density (.cube) files. At least

one path must be provided.

Parameters:

operations (str | None, optional) – How to combine each file in paths with the running total of those before it, as a string with one of "+", "-", "x", "/" per character. Must have exactly len(paths) - 1 characters. Not required (and ignored) when only one path is given, defaults to None.

Raises:
  • ValueError – If all Miller indices are 0: cannot slice through the origin.

  • ValueError – If no paths are provided.

  • ValueError – If the length of operations does not match len(paths) - 1.

  • ValueError – If operations contains a character other than "+-x/".

  • ValueError – If cube files being combined have mismatched grid shapes.

  • ZeroDivisionError – If a "/" operation would divide by a zero-valued voxel.

densities: list[tuple[NDArray[float64], Atoms]]

The (data, atoms) pair loaded from each cube file in paths, in order.

extract_slice(n_points: int = 1000, interp_order: int = 5) tuple[Any, NDArray[float64], NDArray[float64], NDArray[float64], NDArray[float64], NDArray[float64]]

Sample the charge density on the \([hkl]\) plane at fractional offset.

For each point on a 2D \((t_1, t_2)\) grid centred on the plane origin:

  1. Compute its Cartesian position:

    \(p = \text{origin} + t_1 v_1 + t_2 v_2\)

  2. Convert to fractional coordinates: s = p @ inv(cell)

  3. Wrap into the unit cell to enforce periodic boundaries.

  4. Map fractional -> voxel index and interpolate via scipy.ndimage.map_coordinates.

Parameters:
  • n_points (int, optional) – Number of points to sample on the slice, defaults to 1000

  • interp_order (int, optional) – The polynomial degree for interpolation, defaults to 5

Returns:

A tuple containing:

  • density (REAL ARRAY) – Charge density in \(e/a_0^3\), shape (n_points, n_points)

  • v1 (REAL ARRAY) – First unit vector spanning the \([hkl]\) plane

  • v2 (REAL ARRAY) – Second unit vector spanning the \([hkl]\) plane

  • t1 (REAL ARRAY) – Scalar grid along \(v_1\)

  • t2 (REAL ARRAY) – Scalar grid along \(v_2\)

  • origin (float) – Cartesian slice origin in Bohr

inplane_basis() tuple[NDArray[float64], NDArray[float64], NDArray[float64]]

Return two orthonormal Cartesian vectors \((v_1, v_2)\) spanning the \([hkl]\) plane and the unit plane-normal \(\hat{n}\).

1. The plane normal in Cartesian space is \(n = ha^* + kb^* + lc^*\) where \(a^*, b^*, c^*\) are the reciprocal lattice vectors (rows of inv(cell)^T). 2. Two fractional-space null vectors of \([hkl]\) are found via SVD. 3. Those are converted to Cartesian, then Gram-Schmidt orthonormalised so the axes are perpendicular in real-space.

plane_origin() NDArray[float64]

Return a Cartesian point lying on the plane \(ha + kb + lc =\) offset.

offset is a dimensionless fractional intercept (0-1 spans one interplanar period). Pick the simplest fractional coordinate satisfying the plane equation.

Returns:

Origin of the slice.

Return type:

REAL ARRAY

project_atoms(v1: NDArray[float64], v2: NDArray[float64], n_hat: NDArray[float64], origin: NDArray[float64], thickness: float = 0.5) tuple[NDArray[float64], NDArray[float64], list[str]]

Project atoms within thickness of the slice plane onto \((v_1, v_2)\) axes.

Parameters:
  • v1 (REAL ARRAY) – First unit vector spanning the \([hkl]\) plane

  • v2 (REAL ARRAY) – Second unit vector spanning the \([hkl]\) plane

  • n_hat (REAL ARRAY) – Unit vector defining the slice

  • origin (REAL ARRAY) – The origin of the slice

  • thickness (float, optional) – The Cartesian distance perpendicular to the plane to consider atoms as lying on the slice, defaults to 0.5 Bohr.

Returns:

A tuple containing:

  • t1_proj (REAL ARRAY) – Positions of label along \(v_1\)

  • t2_proj (REAL ARRAY) – Positions of label along \(v_2\)

  • syms (list[str]) – List of atom labels

class conquest2a.density.plot_densities(source: density | chden | bandden, show_atoms: bool = False, extension: str = 'png', cbar_label: str | None = None)

Helper class to slice, analyse, and plot volumetric density data.

Takes a density instance directly, or with chden or bandden that exposes a .density attribute

Parameters:
  • source (density | chden | bandden) – The density data to plot.

  • show_atoms (bool, optional) – Whether to show atom labels on the plot, defaults to False.

  • extension (str, optional) – File extension used for auto-generated filenames (leading “.” optional), defaults to “png”.

  • cbar_label (str | None, optional) – Label for the colour bar. If not given, a sensible default is chosen based on the type of source – charge density and band density are conventionally expressed in different units.

plot_slice(density_grid: NDArray[float64], t1: NDArray[float64], t2: NDArray[float64], v1: NDArray[float64], v2: NDArray[float64], atom_data: tuple[NDArray[float64], NDArray[float64], list[str]] | None = None, log_scale: bool = False, vmin: float | None = 0.0, vmax: float | None = None, output: str | None = None, figsize: tuple[float, float] | None = None, savefig_kwargs: dict[str, Any] | None = None, **imshow_kwargs: Any) None

Create the plot and save it to disk.

Parameters:
  • density_grid (REAL ARRAY) – Sliced density data, shape (n_points, n_points).

  • t1 (REAL ARRAY) – Scalar grid along \(v_1\)

  • t2 (REAL ARRAY) – Scalar grid along \(v_2\)

  • v1 (REAL ARRAY) – First unit vector spanning the \([hkl]\) plane

  • v2 (REAL ARRAY) – Second unit vector spanning the \([hkl]\) plane

  • atom_data (tuple[REAL_ARRAY, REAL_ARRAY, list[str]] | None, optional) – Data for atom labels, see density.project_atoms(), defaults to None

  • log_scale (bool, optional) – Whether to plot the density on a base-10 logarithmic scale - useful for revealing details without colour clipping, defaults to False

  • vmin (float | None, optional) – Minimum value to set the colour scale at (ignored if log_scale), defaults to 0.0

  • vmax (float | None, optional) – Maximum value to set the colour scale at (ignored if log_scale), defaults to None

  • output (str | None, optional) – Filename to save as. Will save with a useful name if not provided, defaults to None

  • figsize (tuple[float, float] | None, optional) – Override the auto-computed figure size, defaults to None

  • savefig_kwargs (dict[str, Any] | None, optional) – Extra keyword arguments forwarded to matplotlib.pyplot.savefig() (e.g. dpi, transparent, bbox_inches), defaults to None

  • imshow_kwargs – Any further keyword arguments (e.g. cmap, interpolation, alpha, aspect) are forwarded to matplotlib.pyplot.imshow(), overriding the defaults.

run(filename: str | None = None, thickness: float = 0.5, log_scale: bool = False, vmin: float | None = 0.0, vmax: float | None = None, figsize: tuple[float, float] | None = None, savefig_kwargs: dict[str, Any] | None = None, **imshow_kwargs: Any) None

Runs the full sequence of steps:

  1. Extracts a slice from the underlying density

  2. Optionally projects nearby atoms onto the slice

  3. Plots and saves the result

Parameters:
  • filename (str | None, optional) – Filename to save as. Will save with a useful name if not provided, defaults to None

  • thickness (float, optional) – The Cartesian distance perpendicular to the plane to consider atoms as lying on the slice, defaults to 0.5 Bohr.

  • log_scale (bool, optional) – Whether to plot the density on a base-10 logarithmic scale, defaults to False

  • vmin (float | None, optional) – Minimum value to set the colour scale at, defaults to 0.0

  • vmax (float | None, optional) – Maximum value to set the colour scale at, defaults to None

  • figsize (tuple[float, float] | None, optional) – Override the auto-computed figure size, defaults to None

  • savefig_kwargs (dict[str, Any] | None, optional) – Extra keyword arguments forwarded to matplotlib.pyplot.savefig(), defaults to None

  • imshow_kwargs – Any further keyword arguments (e.g. cmap, interpolation) are forwarded to matplotlib.pyplot.imshow().