Analyzing AKMC output

Use soda_analyzer to inspect states, elementary processes, and the sequence of transitions taken by AKMC. These functions apply to the EAM, GAP, and MACE examples.

Run the snippets after the calculation has finished, from its run directory, with the matching SODA environment activated. For a copied database, use the standalone analyzer setup and replace output/output.h5 with its path. Loading these records does not require the original force potential or model.

What the database contains

Item	Meaning
State	Registered atomic structure, identified by stateID. A state can be visited more than once.
Process	Discovered transition rooted at a state, with reactant, saddle, product, barriers, and prefactors. Its procID is local to that state.
Search summary	Optional record of a search attempt. Searches can fail or rediscover a process, so summary counts differ from process counts.
KMC transition	One selected process and its destination, recorded in dynamics_table with a waiting time and cumulative simulated time.

The process catalogue lists available transitions. The dynamics table lists transitions actually taken. Neither the number of stored states nor the number of completed searches is the KMC step count. Unlike APE analysis, AKMC does not require descriptor-cluster queries.

Results at each phase

Function signatures show key arguments; ... marks additional options in the linked API reference.

AKMC samples the system’s time evolution, and the trajectory ends wherever the run was stopped. Use the catalogue and trajectory queries directly to inspect the saved exploration.

Result to inspect	Helper or property
Registered minima	db.state_list with state_load=True, or db.read_state(stateID, proc_load, structure_load) for one state; use state.atoms for its structure.
Search attempts, when summaries were saved	state.search_stats and state.termination_list.
Discovered processes, barriers, and prefactors	state.get_proc_list(structure_load=None), then process.get_barrier_list() and process.get_prefactor_list().
Recorded process-to-state connections	Read procid and product_stnum from each state’s process_table; see catalogue connections.
Transitions actually taken	db.read_dynamics() for transition triples.
KMC history with times, barriers, and rates	db.h5["dynamics_table"][()].reshape(-1); see the history example below.

Catalogue connections can include processes never taken in the KMC trajectory.

Read states, barriers, and rates

This opens the database read-only and lists its processes:

import soda_analyzer as soda

temperature = 600.0  # Match [Main] temperature in this run.
db = soda.SodaLibrary(
    "output/output.h5",
    workflow="akmc",
    state_load=True,
    proc_load=True,
    structure_load=False,
    storesummariesinfile=False,
)
try:
    print("States:", db.state_num)
    for state in db.state_list:
        print("State:", state.stateID, "Processes:", state.proc_num)
        for process in state.proc_list:
            print("  Process:", process.procID,
                  "Barriers [eV]:", process.get_barrier_list(),
                  "Prefactors [1/s]:", process.get_prefactor_list(),
                  "Rates [1/s]:", process.get_rates(T=temperature))
finally:
    db.close()

Arrays above are ordered forward, reverse. Barriers are measured from the corresponding minimum to the saddle. get_rates() defaults to 300 K; pass the run’s temperature explicitly. Recalculating rates at another temperature does not produce a new KMC trajectory or establish sufficient sampling at that temperature. The KMC rate table also applies its thermal window, so a catalogue process need not have been eligible for selection.

storesummariesinfile=False skips optional per-search summaries, making this snippet suitable for runs launched with either summary setting. Set it to True if you saved summaries and need to inspect them. For large databases, use proc_load=False and call state.get_proc_list() for selected states. Keep the database open while querying its objects.

For the single-adatom EAM benchmark, the expected catalogue contains 36 states and four processes per state: the adatom can hop to neighboring sites in four surface directions. The supplied run takes 20,000 KMC transitions, repeatedly visiting those states. Fewer than four recorded processes at a state indicates incomplete discovery there. These expected counts apply to the fixed-surface adatom example, not the island or oxygen examples.

State and process functions

The summaries below come from the source docstrings. Here, state belongs to db.state_list, and process belongs to that state’s proc_list.

get_proc_list(structure_load=None)	Load and cache process objects belonging to this state.
get_atoms()	Build an ASE structure from stored coordinates and fixed-atom flags.
get_barrier_list()	Read forward and reverse process activation barriers.
get_prefactor_list()	Read forward and reverse process attempt frequencies.
get_rates(T=300)	Calculate forward and reverse rates from this process summary.
get_stat_repeats()	Read the number of times this process was rediscovered.
get_rsp_atoms()	Build ASE structures for the reactant, saddle, and product.
read_dynamics()	Read recorded AKMC state-to-state transitions.

Read catalogue connections

AKMC stores resolved destinations in each state’s process_table. The read_process_product_states() helper reads APE’s environment connection table, which ordinary AKMC output does not populate. An empty result from that helper therefore does not mean that AKMC has no recorded connections.

import soda_analyzer as soda

db = soda.SodaLibrary("output/output.h5", workflow="akmc",
                      state_load=False, structure_load=False,
                      storesummariesinfile=False)
try:
    edges = []
    for state_id in db.state_id_list:
        table = db.h5.get(f"States/{state_id}/process_table")
        if table is None:
            continue
        for row in table[()].reshape(-1):
            process_id = int(row["procid"])
            product_id = int(row["product_stnum"])
            if process_id >= 0 and product_id >= 0:
                edges.append((int(state_id), process_id, product_id))
    print("Resolved catalogue connections:", len(edges))
    print("First (state, process, product state) triples:", edges[:5])
finally:
    db.close()

A negative product_stnum means the destination has not been resolved to a registered state; the example omits those rows. These connections describe the known network. Use the KMC history below for transitions actually taken.

Read the KMC history and simulated time

db.read_dynamics() returns an (N, 3) NumPy array in recorded order, with one row per KMC transition: [reactant_state_id, process_id, product_state_id]. The example below also extracts the visited-state ID sequence and reads the full structured table for time, rate, and energy fields:

import soda_analyzer as soda

db = soda.SodaLibrary("output/output.h5", workflow="akmc",
                      state_load=False, structure_load=False,
                      storesummariesinfile=False)
try:
    if "dynamics_table" not in db.h5:
        print("No KMC history is present.")
    else:
        history = db.h5["dynamics_table"][()].reshape(-1)
        print("Recorded KMC steps:", len(history))
        steps = db.read_dynamics()
        print("First transition triples:", steps[:5])
        state_ids = (
            [int(steps[0, 0])] + steps[:, 2].astype(int).tolist()
            if len(steps) else []
        )
        print("First visited state IDs:", state_ids[:10])
        for row in history[:5]:
            print("Step:", int(row["step_number"]),
                  "From:", int(row["reactant_stnum"]),
                  "Process:", int(row["process_procid"]),
                  "To:", int(row["product_stnum"]),
                  "Waiting time [s]:", float(row["step_time"]),
                  "Simulated time [s]:", float(row["total_time"]))
        if len(history):
            print("Recorded simulated time [s]:", float(history[-1]["total_time"]))
        else:
            print("No KMC transition was recorded.")
finally:
    db.close()

state_ids includes the starting state and preserves repeated visits. For example, 3 → 7 → 3 gives [3, 7, 3]. With at least one recorded transition, N transitions give N + 1 state IDs; an empty history gives [].

Field	Meaning
step_number	Zero-based index in the recorded transition history.
reactant_stnum, process_procid, product_stnum	Starting state, selected process within that state, and destination state.
step_time, total_time	This step’s waiting time and accumulated simulated time, in seconds. These are not wall-clock timings.
barrier, rate	Selected forward process barrier in eV and its stored rate in s⁻¹. The rate is not the sum of all escape rates.
energy	Reactant-state energy in eV before this transition.

Normally the waiting time is sampled from an exponential distribution using the sum of eligible escape rates. [Debug] use_mean_time = true instead uses its mean. An empty history can coexist with discovered processes if the run ended before a KMC step.

Export a selected transition

This selects the first recorded KMC transition and writes its reactant, saddle, and product as three frames. It does not assume state or process 0 was selected:

from ase.io import write
import soda_analyzer as soda

db = soda.SodaLibrary("output/output.h5", workflow="akmc",
                      state_load=False, structure_load=False,
                      storesummariesinfile=False)
try:
    if "dynamics_table" not in db.h5 or db.h5["dynamics_table"].size == 0:
        print("No KMC transition is available to export.")
    else:
        first = db.h5["dynamics_table"][()].reshape(-1)[0]
        state_id = int(first["reactant_stnum"])
        process_id = int(first["process_procid"])
        images = [
            soda.get_atoms(db.h5, stateID=state_id,
                           procID=process_id, rsp=rsp, delta=True)
            for rsp in (0, 1, 2)
        ]
        write("first_kmc_process.xyz", images)
        print("Wrote reactant, saddle, product to first_kmc_process.xyz")
finally:
    db.close()

Open first_kmc_process.xyz in an atomistic viewer to inspect the moving atoms. rsp=0, 1, and 2 select reactant, saddle, and product; delta=True reconstructs the process structures from their stored changes. These three frames illustrate a transition path, not a time-resolved molecular-dynamics trajectory.

When comparing runs, compare structures, barriers, rate settings, and statistics rather than expecting state IDs or event sequences to match. Longer kinetic analyses also need sufficient sampling and correct handling of periodic crossings; a single step cannot determine a diffusion coefficient.

See the states and processes API reference and SodaLibrary reference for full signatures, loading options, and return values.