Start a calculation

Install and choose an interface

Use Python 3.12 in a dedicated environment. Install the supplied wheel, for example:

python -m venv .venv
.venv\Scripts\python -m pip install path\to\hcphases-0.3.2-py3-none-any.whl

NumPy and PyYAML are the runtime dependencies. No package-index publication is assumed by this command. For exact reproduction of a saved run, install its recorded package artifact and NumPy version. Current bounded execution requires Windows; workers=1 is the implemented sampling configuration.

Read Choose and define a model to choose between the higher-degree and vertex interfaces. The units are \(k_B=\hbar=1\); T is temperature and \(\delta_\tau=1/(TM)\) is the imaginary-time step. The shared Homological-code Hamiltonians reference fixes the projector normalization of J,K.

A complete higher-degree calculation

This example measures two Wilson operators, their connected correlation and a local X expectation for a small two-dimensional edge-qubit model. It shows the full plan/run/audit/load/continue workflow. The short histories make the example quick to execute; they do not establish precision or equilibration.

Save the following function in your own project’s source directory. Pass it a fresh project directory and a nonnegative integer seed that you own.

from pathlib import Path
from hcphases import (
    ModelSpec, PeriodicLattice, BinaryFields, TrotterSpec,
    Cell, LocalX, Wilson, Connected, coordinate_cycle,
    SamplingPlan, StoragePlan, DiagnosticsPlan, ResourceBudget,
    plan, run, replay, load_result, resume,
)

def general_example(project, seed):
    project = Path(project)
    project.mkdir(parents=True, exist_ok=False)
    lattice = PeriodicLattice((2, 2))
    model = ModelSpec(d=2, p=1, N=2, lattice=lattice,
                      J=1.0, K=1.0, T=1.0,
                      sources=BinaryFields(g=0.35, h=0.2))
    first = Wilson(coordinate_cycle(lattice, axes=(0,), fixed={1: 0}),
                   lattice_lengths=lattice.lengths, model_p=1)
    second = Wilson(coordinate_cycle(lattice, axes=(0,), fixed={1: 1}),
                    lattice_lengths=lattice.lengths, model_p=1)
    observables = {
        "loop_a": first,
        "loop_b": second,
        "connected": Connected(first, second),
        "local_x": LocalX(Cell((0, 0), (0,))),
    }
    budget = ResourceBudget(wall_seconds=60, memory_bytes=2**30,
                            output_bytes=32*2**20, blas_threads=1)
    job = plan(
        model, observables, representation=TrotterSpec(slices=2),
        sampling=SamplingPlan(chains=2, starts=("ordered", "hot"),
                              warmup_sweeps=2, measurement_sweeps=8,
                              measure_every=2),
        storage=StoragePlan(snapshot_every=2, checkpoint_every=3),
        diagnostics=DiagnosticsPlan(batch_sizes=(1, 2)),
        resources=budget, seed=seed,
    )
    job.write(project / "plan.json")
    result = run(job, output=project / "runs" / "first")
    if result.execution_status != "completed":
        raise RuntimeError(result.receipt)
    audit = replay(result.run_path, output=project / "audits" / "replay",
                   budget=budget)
    if audit["status"] != "passed":
        raise RuntimeError(audit)
    loaded = load_result(result.run_path,
                         output=project / "audits" / "load", budget=budget)
    continued = resume(result.run_path, "checkpoints-index.json",
                       output=project / "runs" / "continued",
                       additional_sweeps=4)
    if continued.execution_status != "completed":
        raise RuntimeError(continued.receipt)
    return job, result, audit, loaded, continued

Run this from your own study working directory, with output outside the installed package and its source checkout. Call it, for example, with general_example("my-project", seed=12345). In a script, put that call under if __name__ == "__main__":. The seed is illustrative: select and record your own independent stream roots.

Inspect the returned objects:

job, result, audit, loaded, continued = general_example("my-project", 12345)
print(job.capability_report)
print(job.resource_estimate)
print(result.estimates["Z:connected"])
print(result.estimates["X:local_x"])
print(result.diagnostics["Z:connected"])

Connected produces pooled a,b,d means and mean=d-a*b. X: and Z: identify different measurement ensembles, not two columns of one joint quantum measurement. The continuation returns the new segment’s statistics; it does not pool the parent history automatically.

A complete vertex-code calculation

The vertex interface describes the quantum code through K,g,T and derives its classical couplings internally. The following isotropic specialization uses \(\delta_\tau=1\), \(K=2\kappa\) and \(g=-\log(\tanh\kappa)/2\). See The vertex code and its two-dimensional image.

import math
from pathlib import Path
from hcphases import (
    VertexModel, VertexObservable, VertexSamplingPlan,
    ResourceBudget, run_vertex, read_vertex_records, audit_vertex_run,
    rank_normalized_split_rhat, bulk_ess, tail_ess,
)

def vertex_example(project, seed):
    project = Path(project)
    project.mkdir(parents=True, exist_ok=False)
    n, kappa = 4, 0.44
    model = VertexModel(length=n, slices=n, K=2*kappa,
                        g=-math.log(math.tanh(kappa))/2, T=1/n)
    observables = {
        "m2": VertexObservable(kind="magnetization", power=2),
        "m4": VertexObservable(kind="magnetization", power=4),
        "pair": VertexObservable(points=((0, 0), (n//2, 0))),
    }
    sampling = VertexSamplingPlan(
        chains=2, starts=("ordered", "hot"), warmup_steps=8,
        sample_steps=32, measure_every=2, snapshot_every=3,
    )
    budget = ResourceBudget(wall_seconds=60, memory_bytes=2**30,
                            output_bytes=8*2**20, blas_threads=1)
    run_path = project / "runs" / "first"
    receipt = run_vertex(model, observables, sampling=sampling,
                         resources=budget, seed=seed, output=run_path)
    if receipt["terminal_reason"] != "completed":
        raise RuntimeError(receipt)
    audit = audit_vertex_run(run_path, output=project / "audits" / "replay",
                             resources=budget)
    if audit["terminal_reason"] != "completed":
        raise RuntimeError(audit)
    metadata, samples = read_vertex_records(run_path)
    index = {name: j for j, name in enumerate(metadata["columns"])}
    m2 = samples[:, :, index["m2"]]
    m4 = samples[:, :, index["m4"]]
    mean2 = float(m2.mean())
    binder = None if mean2 == 0 else float(1 - m4.mean() / (3*mean2**2))
    diagnostics = {"rhat": rank_normalized_split_rhat(m2),
                   "bulk_ess": bulk_ess(m2), "tail_ess": tail_ess(m2)}
    return metadata, samples, binder, diagnostics

samples has axes (chain, retained_record, observable). Always obtain the observable order from metadata["columns"]; it need not match insertion order in your original Python dictionary. The Binder ratio is formed from pooled moments, not by averaging per-record ratios. Its uncertainty needs a joint analysis of the m2 and m4 columns; the diagnostic values above alone do not supply an interval for that ratio.

Where to go next