gemdat.rdf

RDFCollection

Bases: list[RDFData]

Collection to store group of radial distribution data.

plot(*, module, **kwargs)

See gemdat.plots.radial_distribution for more info.

Source code in src/gemdat/rdf.py

@plot_backend
def plot(self, *, module, **kwargs):
    """See [gemdat.plots.radial_distribution][] for more info."""
    return module.radial_distribution(rdfs=self, **kwargs)

RDFData(x, y, label, state) dataclass

Data class for storing radial distribution data.

Parameters:

• x (ndarray) –

  1D array with x data (bins)

• y (ndarray) –

  1D array with y data (counts)

• label (str) –

  Distance to species with this symbol label

• state (str) –

  State that the floating species is in, e.g. the jump that it is making.

plot(*, module, **kwargs)

See gemdat.plots.radial_distribution for more info.

Source code in src/gemdat/rdf.py

@plot_backend
def plot(self, *, module, **kwargs):
    """See [gemdat.plots.radial_distribution][] for more info."""
    return module.radial_distribution(rdfs=[self], **kwargs)

radial_distribution(*, transitions, floating_specie, max_dist=5.0, resolution=0.1)

Calculate and sum RDFs for the floating species in the given sites data.

Parameters:

• transitions (Transitions) –

  Input transitions data

• floating_specie (str) –

  Name of the floating specie

• max_dist (float, default: 5.0 ) –

  Max distance for rdf calculation

• resolution (float, default: 0.1 ) –

  Width of the bins

Returns:

• rdfs ( dict[str, RDFCollection] ) –

  Dictionary with rdf arrays per symbol

Source code in src/gemdat/rdf.py

def radial_distribution(
    *,
    transitions: Transitions,
    floating_specie: str,
    max_dist: float = 5.0,
    resolution: float = 0.1,
) -> dict[str, RDFCollection]:
    """Calculate and sum RDFs for the floating species in the given sites data.

    Parameters
    ----------
    transitions: Transitions
        Input transitions data
    floating_specie : str
        Name of the floating specie
    max_dist : float, optional
        Max distance for rdf calculation
    resolution : float, optional
        Width of the bins

    Returns
    -------
    rdfs : dict[str, RDFCollection]
        Dictionary with rdf arrays per symbol
    """
    # note: needs trajectory with ALL species
    trajectory = transitions.trajectory
    sites = transitions.sites
    base_structure = trajectory.get_structure(0)
    assert isinstance(base_structure, Structure)
    lattice = trajectory.get_lattice()

    coords = trajectory.positions
    sp_coords = trajectory.filter(floating_specie).positions

    states2str = _get_states(sites.labels)  # type: ignore
    states_array = _get_states_array(transitions, sites.labels)  # type: ignore
    symbol_indices = _get_symbol_indices(base_structure)

    bins = np.arange(0, max_dist + resolution, resolution)
    length = len(bins) + 1

    rdfs: dict[tuple[str, str], np.ndarray] = defaultdict(lambda: np.zeros(length, dtype=int))

    n_steps = len(trajectory)

    for i in track(range(n_steps), transient=True):
        t_coords = coords[i]
        t_sp_coords = sp_coords[i]

        dists = lattice.get_all_distances(t_sp_coords, t_coords)

        rdf = np.digitize(dists, bins, right=True)

        states = np.unique(states_array[i], axis=0)

        t_states = states_array[i]

        for state in states:
            k_idx = np.argwhere(t_states == state)
            state_str = states2str[state]

            for symbol, symbol_idx in symbol_indices.items():
                rdf_state = rdf[k_idx, symbol_idx].flatten()
                rdfs[state_str, symbol] += np.bincount(rdf_state, minlength=length)

    ret: dict[str, RDFCollection] = {}

    for (state, symbol), values in rdfs.items():
        rdf_data = RDFData(
            x=bins,
            # Drop last element with distance > max_dist
            y=values[:-1],
            label=symbol,
            state=state,
        )
        ret.setdefault(state, RDFCollection())
        ret[state].append(rdf_data)

    return ret

radial_distribution_between_species(*, trajectory, specie_1, specie_2, max_dist=5.0, resolution=0.1)

Calculate RDFs from specie_1 to specie_2.

Parameters:

• trajectory (Trajectory) –

  Input trajectory.

• specie_1 (str | Collection[str]) –

  Name of specie or list of species

• specie_2 (str | Collection[str]) –

  Name of specie or list of species

• max_dist (float, default: 5.0 ) –

  Max distance for rdf calculation

• resolution (float, default: 0.1 ) –

  Width of the bins

Returns:

• rdf ( RDFData ) –

  RDF data for the given species.

Source code in src/gemdat/rdf.py

def radial_distribution_between_species(
    *,
    trajectory: Trajectory,
    specie_1: str | Collection[str],
    specie_2: str | Collection[str],
    max_dist: float = 5.0,
    resolution: float = 0.1,
) -> RDFData:
    """Calculate RDFs from specie_1 to specie_2.

    Parameters
    ----------
    trajectory: Trajectory
        Input trajectory.
    specie_1: str | list[str]
        Name of specie or list of species
    specie_2: str | list[str]
        Name of specie or list of species
    max_dist: float, optional
        Max distance for rdf calculation
    resolution: float, optional
        Width of the bins

    Returns
    -------
    rdf : RDFData
        RDF data for the given species.
    """
    coords_1 = trajectory.filter(specie_1).coords
    coords_2 = trajectory.filter(specie_2).coords

    if coords_1.ndim == 2:
        coords_1 = coords_1[None, :, :]
    if coords_2.ndim == 2:
        coords_2 = coords_2[None, :, :]

    num_time_steps, num_atoms_2, _ = coords_2.shape
    _, num_atoms_1, _ = coords_1.shape

    lattices = trajectory.lattice
    if lattices is None:
        raise RuntimeError('trajectory.lattice unexpectedly None')

    if trajectory.constant_lattice:
        lattices = np.array([lattices] * num_time_steps)

    bins = np.arange(0, max_dist + resolution, resolution)

    def normalize(radius: np.ndarray, particle_vol: float) -> np.ndarray:
        """Normalize bin to volume."""
        shell = (radius + resolution) ** 3 - radius**3
        return particle_vol * (4 / 3) * np.pi * shell

    counts = np.zeros(len(bins) - 1, dtype=float)

    for t, lat in enumerate(lattices):
        lat = Lattice(lat)
        dist_mat = lat.get_all_distances(coords_1[t], coords_2[t])
        if specie_1 == specie_2:
            mask = ~np.eye(num_atoms_1, dtype=bool)
            dists_t = dist_mat[mask]  # removes i == j for self pairs
        else:
            dists_t = dist_mat.ravel()
        hist_t, _ = np.histogram(dists_t, bins=bins, density=False)

        particle_vol_t = num_atoms_2 * num_atoms_1 / lat.volume
        norm_t = normalize(bins, particle_vol_t)[:-1]

        counts += hist_t / norm_t

    counts /= num_time_steps

    str1 = specie_1 if isinstance(specie_1, str) else '/'.join(specie_1)
    str2 = specie_2 if isinstance(specie_2, str) else '/'.join(specie_2)

    return RDFData(x=bins[:-1], y=counts, label=f'{str1}-{str2}', state='')