Partitioned TreeTNs

tensor4all-partitionedtreetn stores TreeTN subdomains as eagerly masked patches. It replaces the retired tensor4all-partitionedtt crate and supports named chains, branched trees, and multiple site indices on one node.

This crate provides partition algebra, TreeTN-general adaptive patching, and adaptive patched interpolation of a function through any tree interpolation engine. It does not implement an engine itself; the TreeTCI engine lives in tensor4all-treetci.

Construct an eager patch

Projectors use zero-based coordinates and full index identity. Construction retains every site axis but masks values outside the selected coordinates:

use tensor4all_core::{DynIndex, IdxTensor};
use tensor4all_partitionedtreetn::{Projector, SubDomainTreeTN};
use tensor4all_treetn::TreeTN;
fn main() -> Result<(), Box<dyn std::error::Error>> {
let site = DynIndex::new_dyn(2);
let tensor = IdxTensor::from_dense(
    vec![site.clone()],
    vec![3.0_f64, 1.0e12],
)?;
let tree = TreeTN::from_tensors(vec![tensor], vec!["root".to_string()])?;
let patch = SubDomainTreeTN::new(
    tree,
    Projector::from_pairs([(site.clone(), 0)])?,
)?;

let node = patch.data().node_index(&"root".to_string()).ok_or("missing root")?;
assert_eq!(patch.data().tensor(node).ok_or("missing tensor")?.to_vec::<f64>()?,
           vec![3.0, 0.0]);
assert!((patch.norm_squared()? - 9.0).abs() < 1.0e-12);
Ok(())
}

Norms, inner products, contraction, truncation, and summation use this stored masked value directly. No projector is re-applied and no full network is densified.

Adaptive patching

Every truncating or contracting operation takes an explicit existing node name as its center. PatchingOptions::cutoff is a non-negative dimensionless discarded-weight budget shared across patches, not an SVD truncation-policy knob, and PatchingOptions::allocation selects how it is spread. BudgetAllocation::Volume (the default, arXiv:2602.22372) assigns the absolute local threshold cutoff * ||F||^2 * volume_p / total_volume, applies each whole threshold at the patch’s local SVD truncations, drops a patch whose squared norm is at most its threshold, and then splits patches that remain above the bond cap. BudgetAllocation::Norm instead applies the relative squared tail-sum cutoff cutoff at every local SVD — the same rule as an ITensor-style relative cutoff in TreeTN::truncate — so each patch is accurate relative to itself, no patch is dropped, and projected children inherit no absolute budget. The same choice applies to add_with_patching, contract_adaptive, and truncate_adaptive. cutoff is best effort for the final whole-network error; max_bond_dim is a hard cap. Inputs that share an equal projector key are summed before patching:

use tensor4all_core::{DynIndex, IdxTensor};
use tensor4all_partitionedtreetn::{
    add_with_patching, BudgetAllocation, PatchSplitStrategy, PatchingOptions, SubDomainTreeTN,
};
use tensor4all_treetn::TreeTN;
fn main() -> Result<(), Box<dyn std::error::Error>> {
let site0 = DynIndex::new_dyn(2);
let bond = DynIndex::new_dyn(2);
let site1 = DynIndex::new_dyn(2);
let left = IdxTensor::from_dense(
    vec![site0.clone(), bond.clone()],
    vec![1.0_f64, 0.0, 0.0, 1.0],
)?;
let right = IdxTensor::from_dense(
    vec![bond, site1],
    vec![1.0_f64, 0.0, 0.0, 1.0],
)?;
let patch = SubDomainTreeTN::from_treetn(
    TreeTN::from_tensors(vec![left, right], vec![0usize, 1])?,
)?;
let result = add_with_patching(
    vec![patch],
    &0,
    &PatchingOptions {
        cutoff: 0.0,
        allocation: BudgetAllocation::Volume,
        max_bond_dim: Some(1),
        patch_order: vec![site0],
        split_strategy: PatchSplitStrategy::Sequential,
    },
)?;

assert_eq!(result.len(), 2);
assert!(result.values().all(|patch| patch.max_bond_dim() <= 1));
Ok(())
}

PatchSplitStrategy::Sequential follows patch_order. The default ExactParameterGain forms and budget-truncates every candidate’s children, then compares checked sums of logical local tensor element counts. Structured storage payload length and AD state are not used as the metric.

Adaptive patched interpolation

adaptive_interpolation::patched_interpolate interpolates a function on an arbitrary named tree directly into a partition. It runs an engine implementing tensor4all_treetn::interpolation::TreeInterpolator, such as tensor4all_treetci::TreeTciInterpolator, on the whole domain. A patch that cannot be accepted is split: the next unfixed site of patch_order is fixed and every child is interpolated in turn. For quantics grids, list the most significant bits first.

The evaluator receives a column-major [n_sites, n_points] batch of zero-based full-domain points in the derived site order (nodes in ascending name order, each node’s sites in the given order) and the run’s EvalContext. Every numerical stage of the run - the engine, re-embedding, measurement, and the norm readout - runs on the engine behind that context; nothing enters the process-global session, so a run completes on the engine you pass in (for example one lane’s engine) and needs no process-global warm-up. The returned partition and report are published as compact host storage and stay usable after the session closes. Each patch caches its samples and hands them to its children, so no point is evaluated twice, and a patch with at most one unfixed site is evaluated exactly without the engine.

The driver has three entries that share one wave loop over the BFS generations of the patch queue. patched_interpolate runs every patch on the caller’s session, one after another. adaptive_interpolation::patched_interpolate_on_lanes fans each wave out over the lanes of a caller-owned hataori_tenferro::TenferroLanes, so every patch runs in the child session of the lane that pulled it; the root wave and waves smaller than the lane count run one patch at a time, so a single item keeps the backend’s whole team. adaptive_interpolation::patched_interpolate_mpi (feature hataori-mpi or hataori-mpi-rsmpi-rt) distributes each wave over MPI ranks through hataori::pmap_lanes, with the root owning the queue and the result. The schedule does not change the report: seeds come from the seed and the patch path, counters are per patch, and the wave merge is in canonical path order.

The accuracy requirement is an ErrorNorm with an ErrorTolerance { rtol, atol }. The default, ErrorNorm::L2, bounds the L2 error over the whole domain by delta = max(atol, rtol * S), where S is an L2 norm of the function, usually given as L2Reference::Given(S). The driver measures every accepted and zero patch itself.

use std::collections::BTreeMap;
use tensor4all_core::{ColMajorArray, ColMajorArrayRef, DynIndex, IdxTensor};
use tensor4all_partitionedtreetn::adaptive_interpolation::{
    patched_interpolate, GlobalL2Error, PatchedInterpolationOptions,
};
use tensor4all_partitionedtreetn::{ErrorNorm, ErrorTolerance, L2Reference};
use tensor4all_tensorbackend::{explicit::Session, CpuBackend};
use tensor4all_treetci::TreeTciInterpolator;
use tensor4all_treetn::interpolation::{EvalContext, Parallelism};
use tensor4all_treetn::NodeNameNetwork;
fn main() -> Result<(), Box<dyn std::error::Error>> {
// A junction "c" of degree three without a site, and three leaves.
let (x, y, z) = (DynIndex::new_dyn(2), DynIndex::new_dyn(3), DynIndex::new_dyn(3));
let mut topology = NodeNameNetwork::new();
for node in ["c", "x", "y", "z"] {
    topology.add_node(node.to_string())?;
}
for leaf in ["x", "y", "z"] {
    topology.add_edge(&"c".to_string(), &leaf.to_string())?;
}
let node_sites = BTreeMap::from([
    ("c".to_string(), vec![]),
    ("x".to_string(), vec![x.clone()]),
    ("y".to_string(), vec![y.clone()]),
    ("z".to_string(), vec![z.clone()]),
]);

// Site order [x, y, z]. f vanishes for x = 0 and has rank three otherwise.
let f = |p: &[usize]| (p[0] * (1 + p[1] + p[2]).pow(2)) as f64;
let values: Vec<f64> = (0..18).map(|k| f(&[k % 2, (k / 2) % 3, k / 6])).collect();
let norm = values.iter().map(|v| v * v).sum::<f64>().sqrt();
let options = PatchedInterpolationOptions::new(3)
    .with_error_norm(ErrorNorm::l2(L2Reference::Given(norm)))
    .with_tolerance(ErrorTolerance { rtol: 1e-10, atol: 0.0 })
    .with_patch_order(vec![x.clone(), y.clone()]);
let mut backend = CpuBackend::with_threads(1)?;
let mut ctx = EvalContext::new(Session::new(backend.engine_mut()), Parallelism::Serial);
let result = patched_interpolate(
    &TreeTciInterpolator::default(),
    topology,
    node_sites,
    ColMajorArray::new(vec![1, 2, 2], vec![3, 1]).expect("2D pivots"),
    &mut ctx,
    |batch: ColMajorArrayRef<'_, usize>,
     _ctx: &mut EvalContext<'_>|
     -> anyhow::Result<Vec<f64>> { Ok(batch.data().chunks(3).map(f).collect()) },
    &options,
)?;

// The x = 0 half is a measured zero patch; the rest is split at y.
assert_eq!(result.report.zero_patches.len(), 1);
assert_eq!(result.partition.len(), 3);
// Every patch was measured exhaustively: the global error is certified.
let error = result.report.norm.l2_error().unwrap();
assert!(matches!(error.global, GlobalL2Error::Certified { .. }));

let reference = IdxTensor::from_dense(vec![x, y, z], values)?;
let dense = result.partition.to_treetn()?.contract_to_tensor()?;
assert!(dense.sub(&reference)?.norm()? <= result.report.norm.delta().unwrap() + 1e-12 * norm);
Ok(())
}

Under ErrorNorm::L2 the allowance is split by patch volume: every accepted or zero patch must have a root-mean-square residual of at most tau = delta / sqrt(|X|), which is also the engine’s absolute tolerance. A patch with at most max(max_exhaustive_points, samples) points (see VerificationOptions) is measured at every point; a larger patch on fresh uniform samples, followed by an independent audit sample. A failed measurement of a converged run reruns the engine with the worst measured points added to its pivots, then splits the patch. What a run can claim is the report’s GlobalL2Error:

• Certified: every contribution is exact or exhaustive, and the absolute error is at most delta up to a small relative margin and a rounding term from a calibrated, not proven, model (rounding_limited says when that term is not below tau); relative_error_bound bounds E / ||f|| when ||f~|| exceeds E.
• Audited: some contribution was sampled, and every sampled one was audited. The audited mean square is an unbiased estimate (its square root is not) with a standard error, never a bound and not a guarantee: a residual concentrated on unsampled points is missed by the acceptance sample and the audit alike, and the standard error, computed from the same points, does not reveal it. A localized feature that enters a patch only through a corner or an edge is such a case; the audited error can then be orders of magnitude too small. More samples make a miss less likely but do not exclude it. The opt-in with_cache_candidates(true) (child patches also start from the largest values the parent already evaluated in them) can prevent a miss where the parent sampled the feature, but not one where the parent never sampled it or where the engine saw the feature without resolving it and the acceptance sample misses the residual.
• AcceptanceOnly: audits were disabled. The combined acceptance statistics are neither a bound nor an estimate.
• ToleranceNotMet: some patch was retained by the minimum patch size without meeting its allowance (below). It takes precedence over the other three, is never certified, and its basis says which of them describes its measured value; an exact or exhaustive basis still bounds the error.

Patch-size bounds

Two options bound the patch size in generalized bits: every unfixed site of a patch counts as one bit, whatever its dimension, so a fused quantics site of dimension 4 is one bit.

• with_min_patch_bits(m) stops splitting where the children would have fewer than m bits. Such a patch is retained with its last engine run and recorded as PatchStatus::WithinTolerance or PatchStatus::ToleranceNotMet; report.tolerance_met() says whether every patch met its allowance, and the error of a run with a ToleranceNotMet patch can exceed delta without limit.
• with_capped_patches(CappedPatches::AcceptUpTo { bits }) accepts a patch that reaches the bond cap when it has at most bits unfixed sites and passes its error check. Larger capped patches still split, and a patch that converges below the cap is never split for its size. The bound must be at least the minimum.

By default neither applies, and the driver may split down to exact patches in the worst case. On binary layouts, a capped bound of at most log2(max(max_exhaustive_points, samples)) bits keeps every capped acceptance exhaustively measured.

use std::collections::BTreeMap;
use tensor4all_core::{ColMajorArray, ColMajorArrayRef, DynIndex, IdxTensor};
use tensor4all_partitionedtreetn::adaptive_interpolation::{
    patched_interpolate, GlobalL2Error, PatchStatus, PatchedInterpolationOptions,
    ToleranceNotMetBasis,
};
use tensor4all_partitionedtreetn::{ErrorNorm, ErrorTolerance, L2Reference};
use tensor4all_tensorbackend::{explicit::Session, CpuBackend};
use tensor4all_treetci::TreeTciInterpolator;
use tensor4all_treetn::interpolation::{EvalContext, Parallelism};
use tensor4all_treetn::NodeNameNetwork;
fn main() -> Result<(), Box<dyn std::error::Error>> {
// f(x, y) = 1 + [x == y] on two sites of dimension 3 has rank three, which
// a cap of two cannot represent. A minimum of two bits retains the root.
let (x, y) = (DynIndex::new_dyn(3), DynIndex::new_dyn(3));
let mut topology = NodeNameNetwork::new();
topology.add_node(0usize)?;
topology.add_node(1usize)?;
topology.add_edge(&0, &1)?;
let node_sites = BTreeMap::from([(0usize, vec![x.clone()]), (1, vec![y.clone()])]);
let f = |p: &[usize]| if p[0] == p[1] { 2.0 } else { 1.0 };
let values: Vec<f64> = (0..9).map(|k| f(&[k % 3, k / 3])).collect();
let norm = values.iter().map(|v| v * v).sum::<f64>().sqrt();
let options = PatchedInterpolationOptions::new(2)
    .with_error_norm(ErrorNorm::l2(L2Reference::Given(norm)))
    .with_tolerance(ErrorTolerance { rtol: 1e-10, atol: 0.0 })
    .with_min_patch_bits(2);
let mut backend = CpuBackend::with_threads(1)?;
let mut ctx = EvalContext::new(Session::new(backend.engine_mut()), Parallelism::Serial);
let result = patched_interpolate(
    &TreeTciInterpolator::default(),
    topology,
    node_sites,
    ColMajorArray::new(vec![0, 0], vec![2, 1]).expect("2D pivots"),
    &mut ctx,
    |batch: ColMajorArrayRef<'_, usize>,
     _ctx: &mut EvalContext<'_>|
     -> anyhow::Result<Vec<f64>> { Ok(batch.data().chunks(2).map(f).collect()) },
    &options,
)?;

assert_eq!(result.report.splits, 0);
assert_eq!(result.report.accepted[0].status, PatchStatus::ToleranceNotMet);
assert!(!result.report.tolerance_met());
let error = result.report.norm.l2_error().unwrap();
let GlobalL2Error::ToleranceNotMet { measured_rms, basis, .. } = &error.global else {
    panic!("expected ToleranceNotMet");
};
// The nine points were measured exhaustively, so the measured error is the
// dense one.
assert!(matches!(basis, ToleranceNotMetBasis::ExactOrExhaustive { .. }));
let reference = IdxTensor::from_dense(vec![x, y], values)?;
let dense = result.partition.to_treetn()?.contract_to_tensor()?;
let residual_rms = dense.sub(&reference)?.norm()? / 3.0;
assert!((measured_rms - residual_rms).abs() <= 1e-12 * norm);
Ok(())
}

The reference norm defaults to L2Reference::Required, which fails before any evaluation unless rtol = 0 (use atol alone) or the root has at most one site (the reference is then exact). L2Reference::MonteCarlo estimates it from uniform root samples; the estimate is heavy-tailed for localized functions and can make the allowance looser than requested. With a tiny or zero tolerance the driver may split down to exact patches, so set max_patches.

ErrorNorm::sampled_max() keeps the M2 criterion: the engine’s sampled error estimate against max(atol, rtol * max_reference), where max_reference is a function value (sampled_max_with_reference(max_abs), or the largest root sample). It runs no measurement, is neither a certified bound nor a measured error, and makes no L2 claim. Under every norm, a patch whose candidate samples are all exactly zero is first screened (and under L2 measured); zero patches are reported in report.zero_patches and left out of the partition, which treats an absent patch as zero.

Execution is sequential. For a fixed seed, a deterministic evaluator, and a deterministic engine, the report and every stored node tensor (values and positional axis order) are identical across runs on fresh threads within one process as long as the measured network values are reproducible. This is tested for f64 on trees with exactly one site per node. On trees with a site-free node or a node with several sites, or with f32/c32 data, the cached evaluator can round differently between threads and processes, an open limitation tracked in tenferro-rs#1963. Issue #795 closed the tracking/test gap through #816; it did not fix the planner. The intended scope is bitwise identical results on the same machine and build across threads, thread counts, and processes, with no cross-machine promise. Generic-path measurements still need deterministic contraction paths; reproducibility across processes is not claimed until a two-process test passes. Issue #791 was fixed by #793: TreeTN site/edge construction and dense-axis order no longer inherit that hash-iteration nondeterminism. site_space is still a set with unspecified iteration order; cached generic evaluation retains the separate planner limitation above.

Reconstruction with a fixed global L2 tolerance

Use reconstruction::reconstruct when approximation must be measured against one immutable target, rather than the local discarded-weight cutoff used above. ReconstructionTarget::from_partition validates and snapshots disjoint patches and pins their combined L2 norm. ReconstructionTolerance { rtol, atol } sets the fixed allowance max(atol, rtol * reference_scale).

The rank goal is soft: a split must improve rank, and a pairwise merge must reduce the sum of operand ranks. Unprofitable sums remain as superposition terms. The output’s regions are disjoint, but terms within a region can overlap. Use regions() to consume them. into_partition() rejects a region containing multiple terms; it never implicitly sums them.

This example is included directly from the checked executable source:

use tensor4all_core::{DynIndex, IdxTensor};
use tensor4all_partitionedtreetn::{reconstruction::*, PartitionedTreeTN, PatchSplitStrategy, SubDomainTreeTN, TreeTN};
fn main() -> Result<(), Box<dyn std::error::Error>> {
    let x = DynIndex::new_dyn(2);
    let y = DynIndex::new_dyn(2);
    let bond = DynIndex::new_dyn(2);
    // Column-major cores for T(x,y) = delta(x,y).
    let tree = TreeTN::from_tensors(
        vec![
            IdxTensor::from_dense(vec![x.clone(), bond.clone()], vec![1.0, 0.0, 0.0, 1.0])?,
            IdxTensor::from_dense(vec![bond, y], vec![1.0, 0.0, 0.0, 1.0])?,
        ],
        vec![0usize, 1],
    )?;
    let reference = tree.clone().to_dense()?; // Small reference example only.
    let partition = PartitionedTreeTN::from_subdomain(SubDomainTreeTN::from_treetn(tree)?)?;
    let target = ReconstructionTarget::from_partition(&partition)?;
    let output = reconstruct(
        &target,
        &0,
        ReconstructionTolerance {
            rtol: 1e-8,
            atol: 0.0,
        },
        &ReconstructionOptions {
            target_bond_dim: Some(1),
            patch_order: vec![x],
            split_strategy: PatchSplitStrategy::Sequential,
            ..Default::default()
        },
    )?;
    assert!((output.report().reference_scale - 2.0_f64.sqrt()).abs() < 1e-12);
    assert_eq!(output.report().region_count, 2);
    assert_eq!(output.report().max_bond_dim, 1);
    assert!(output.report().error_bound <= output.report().absolute_tolerance);
    // Conversion succeeds because each region has one retained term.
    let reconstructed = output.into_partition()?.to_treetn()?.to_dense()?;
    assert!(reconstructed.sub(&reference)?.maxabs()? < 1e-12);
Ok(())
}

ReconstructionTarget::from_tensor_products accepts pairs of patches on the same named topology and independent site spaces. Their output node owns both factors’ external indices. Each product norm is the product of its factor norms; orthogonal product-patch norm squares are then added. Products remain factorized until reconstruction begins. This is a tensor product, not an elementwise product or an induced operator norm.

Each accepted compression is checked against its uncompressed local input by an explicit difference-network norm. Residuals add within a region and combine in quadrature across disjoint regions. The report is a numerical a posteriori bound; it excludes floating-point roundoff. rtol = atol = 0 disables approximate compression. Reaching max_regions retains higher rank without relaxing the error allowance. A partial patch_order list constrains the search independently of any future QFT-selected index subset.

Reconstruction reuses PatchSplitStrategy. Sequential tries only the first unprojected nontrivial index in patch_order and stops that region on no gain or insufficient region capacity, without trying later indices. The default ExactParameterGain compares all permitted candidates by logical parameter count. For contiguous QTT intervals, supply the bits MSB first and select Sequential, even when the TT stores those bits in reverse order. An empty order uses all external indices in deterministic identity order, not numeric bit significance.

Applying a QFT to a subset of sites

ReconstructionTarget::from_subset_operator(&preimage, &center, &operator, &selection, &options) prepares the images of an existing linear operator acting on an ordered subset of the target’s site indices. selection holds one full site index per operator node, in the operator’s own node order; for a quantics Fourier transform its node 0 is the most significant input bit. Build that operator with tensor4all_quanticstransform::quantics_fourier_operator and the crate stays free of a simplett-stack runtime dependency.

The selection may skip sites and spectators keep their identity, dimension, and node assignment. A spectator may share its node with a selected index. Several selected indices may share one node as well: the operator MPO nodes carrying them are fused into one multi-site node, which is exact and keeps the preimage node name. Only a group whose operator nodes are threaded through another owner’s node cannot be fused locally, and that selection is rejected with repair guidance.

The transformed output norm is never measured. SubsetOperatorOptions::unitary selects the amplification factor: false (default) uses the selected-space operator’s Frobenius norm, an upper bound on its induced amplification; true is a caller guarantee that the operator preserves the L2 norm, giving factor one. The preimage’s reference scale is multiplied by that factor, so successive applications propagate the scale instead of recomputing an output norm. For a general operator rtol is therefore relative to that scale, not to the actual ||A x||_2; pass unitary = true for a Fourier transform, whose construction error is accounted separately from the reconstruction bound.

The operator is applied exactly and the images stay separate, so a global direct sum is never formed.

The transform stores frequency bit t at selected position t without an output bit-reversal permutation. For contiguous output patches, supply patch_order = [r1, ..., rR] with PatchSplitStrategy::Sequential.

Level-coupled merge-refine scheduling

reconstruction::schedule_merge_refine(&preimage, &center, &operator, &selection, &subset, tolerance, &MergeRefineOptions) instead runs the complementary input-merge/output-refine trajectory of the patched Fourier algorithm. The preimage must already be the dyadic input leaves of the d selected binary indices, sharing identical spectator constraints. Each leaf fixes a contiguous prefix of those indices, so leaf depths may differ as long as the leaves form a prefix code; all 2^d coordinate assignments must be covered by default, or a sparse subset with CoverageContract::ZeroForMissingLeaves, where every omitted assignment contributes exactly zero. Level zero applies the complete transform once per leaf. Level t merges the input siblings by removing the constraint on k_(d-t+1) and refines the output by fixing r_t, restricting every contribution to its output region before adding it, so no sum over the whole output domain is ever assembled and each computed object is reused by its descendants.

MergeRefineOptions::output_depth stops after that many levels; None (the default) refines every selected bit. MergeRefineOptions::coordinate_groups supplies several coordinate axes explicitly: each CoordinateGroup gives that axis’ selected indices in input order and, separately, in output order, one level advances every non-exhausted axis by one bit, and a two-axis level therefore combines four input children and produces four output children. The default None treats selection as one axis with the documented reversed placement. MergeRefineOptions::target_bond_dim is a soft rank goal: None keeps the trajectory uniform and exact, while Some(goal) makes it adaptive. An output region is then refined only when its retained terms exceed the goal and refining lowers its maximum retained rank. A region that stops refining keeps its output key and still merges its input siblings on every later level, so it consolidates the input prefixes it held instead of freezing them. A merged pair is combined into one term when it fits the goal or lowers the bond dimension (otherwise both operands stay as separate terms of the region’s superposition), and every truncation candidate must fit the item’s share of the global allowance. An unaffordable or gainless candidate is retained exactly rather than violating the accuracy contract, and exceeding MergeRefineOptions::max_terms returns a resource-limit error.

MergeRefineOptions::apply_options opts into a truncating operator application. The schedule then applies the operator with those options and exactly, measures each leaf’s deviation, and charges that sum before any scheduling compression, so an approximate application is never a silent error; a measured application error already above the allowance is rejected. Every other entry point applies exactly. MergeRefineReport records the schedule shape (level_count, applied_operator_count, additions, projections, compression_attempts, compressions, work_items_per_level, peak_work_items, refined_regions, stopped_regions) next to the pinned reference scale, the allowance, and the measured error_bound. Every measured component — an application error, a truncation residual, or a dropped norm — is recorded once, at the level and region where it was measured: components of one level live in disjoint regions and combine by the Euclidean norm, while components of different levels can be nested and add by the triangle inequality. The bound therefore never exceeds the allowance and does not grow with the number of refined regions. It excludes operator-construction error and any error of a caller-approximated operator outside the charged application deviation. Per-input-branch refinement depths with lazy reconciliation and benchmark evidence in the intended patched-input regime remain follow-up work. The transform is never padded and its length never changes: the operator’s input and output are the same selected indices.

use std::collections::HashMap;
use tensor4all_core::{DynIndex, IdxTensor};
use tensor4all_partitionedtreetn::{reconstruction::*, PartitionedTreeTN, Projector, SubDomainTreeTN, TreeTN};
use tensor4all_treetn::{IndexMapping, LinearOperator};
fn main() -> Result<(), Box<dyn std::error::Error>> {
    let site = DynIndex::new_dyn(2);
    let tree = TreeTN::from_tensors(
        vec![IdxTensor::from_dense(vec![site.clone()], vec![3.0, 4.0])?],
        vec![0usize],
    )?;
    // The schedule consumes the 2^d dyadic input leaves of the selected sites.
    let full = SubDomainTreeTN::from_treetn(tree)?;
    let leaves = (0..2)
        .map(
            |value| -> Result<SubDomainTreeTN, Box<dyn std::error::Error>> {
                let projector = Projector::from_pairs([(site.clone(), value)])?;
                full.project(&projector)?.ok_or_else(|| "zero leaf".into())
            },
        )
        .collect::<Result<Vec<_>, _>>()?;
    let preimage =
        ReconstructionTarget::from_partition(&PartitionedTreeTN::from_subdomains(leaves)?)?;

    // A normalized Hadamard transform written as a one-node MPO.
    let (internal_input, internal_output) = (DynIndex::new_dyn(2), DynIndex::new_dyn(2));
    let scale = 1.0 / 2.0_f64.sqrt();
    let mpo = TreeTN::from_tensors(
        vec![IdxTensor::from_dense(
            vec![internal_input.clone(), internal_output.clone()],
            vec![scale, scale, scale, -scale],
        )?],
        vec![0usize],
    )?;
    let mut input_mapping = HashMap::new();
    input_mapping.insert(
        0usize,
        IndexMapping {
            true_index: site.clone(),
            internal_index: internal_input,
        },
    );
    let mut output_mapping = HashMap::new();
    output_mapping.insert(
        0usize,
        IndexMapping {
            true_index: site.clone(),
            internal_index: internal_output,
        },
    );

    let result = schedule_merge_refine(
        &preimage,
        &0,
        &LinearOperator::new(mpo, input_mapping, output_mapping),
        std::slice::from_ref(&site),
        &SubsetOperatorOptions { unitary: true },
        ReconstructionTolerance {
            rtol: 1e-12,
            atol: 0.0,
        },
        &MergeRefineOptions::default(),
    )?;

    // One complete transform per input leaf, one level, one addition per child.
    let report = result.report();
    assert_eq!(report.applied_operator_count, 2);
    assert_eq!(report.additions, 2);
    assert_eq!(report.work_items_per_level, vec![2, 2]);
    assert!((report.reference_scale - 5.0).abs() < 1e-12);

    // The fully refined result is one patch per output coordinate.
    let partition = result.into_partition()?;
    let region = |value: usize| -> Result<f64, Box<dyn std::error::Error>> {
        let patch = partition
            .values()
            .find(|patch| patch.projector().get(&site) == Some(value))
            .ok_or("missing output region")?;
        Ok(patch.norm()?)
    };
    assert!((region(0)? - 7.0 / 2.0_f64.sqrt()).abs() < 1e-12);
    assert!((region(1)? - 1.0 / 2.0_f64.sqrt()).abs() < 1e-12);
Ok(())
}

Dtype and topology

A partition is homogeneous: all patches must use the same IdxTensor scalar dtype and the same named topology and site-index assignment. Both f64 and Complex64 are supported. Topology is not restricted to a chain; a TreeTN with a central named node and three named leaves is a valid partition input. See the Tree Tensor Networks guide for constructing branched networks and selecting contraction/truncation options.

Migration from tensor4all-partitionedtt

The retired tensor4all-partitionedtt crate (removed in #915 D3, settling #825) worked on the positional SimpleTensorTrain representation. This crate carries the same partition algebra on named TreeTN values:

The adaptive interpolation driver needs an engine that implements TreeInterpolator; tensor4all-treetci::TreeTciInterpolator is the TreeTCI one. benchmarks/rust/benchmark_patching.rs is the ported patching benchmark.