atlas_q.mpo_ops#

Matrix Product Operator (MPO) Operations

MPOs represent operators on quantum states in tensor network form: - Hamiltonian evolution - Observable expectation values - Noise channels - Time evolution

Author: ATLAS-Q Contributors Date: October 2025 License: MIT

class atlas_q.mpo_ops.MPO(tensors, n_sites)[source]#

Bases: object

Matrix Product Operator

Represents an operator as a chain of 4-tensors: O = Σ W[0]_{s₀s₀’} W[1]_{s₁s₁’} … W[n-1]_{sₙ₋₁sₙ₋₁’}

Each tensor W[i] has shape [χ_L, d, d, χ_R] where: - χ_L, χ_R: left and right bond dimensions - d: physical dimension (2 for qubits)

Methods

`from_local_ops`(ops[, device])	Create MPO from list of local operators (one per site)
`from_operators`(ops[, device])	Alias for from_local_ops
`identity`(n_sites[, device, dtype])	Create identity MPO

tensors: List[Tensor]#

n_sites: int#

static identity(n_sites, device='cuda', dtype=torch.complex64)[source]#

Create identity MPO

static from_local_ops(ops, device='cuda')[source]#

Create MPO from list of local operators (one per site)

Args:: ops: List of 2×2 operators for each site

static from_operators(ops, device='cuda')[source]#

Alias for from_local_ops

__add__(other)[source]#

Add two MPOs by expanding bond dimensions

For A + B, we create a new MPO with bond dimension χ_A + χ_B that represents the sum of the two operators.

class atlas_q.mpo_ops.MPOBuilder[source]#

Bases: object

Helper class to build common MPOs

Methods

`from_local_terms`(n_sites, local_terms[, ...])	Build Hamiltonian from local interaction terms
`heisenberg_hamiltonian`(n_sites[, Jx, Jy, ...])	Heisenberg Hamiltonian: H = Σᵢ (Jₓ XᵢXᵢ₊₁ + Jᵧ YᵢYᵢ₊₁ + Jᵧ ZᵢZᵢ₊₁)
`identity_mpo`(n_sites[, device, dtype])	Create identity MPO (wrapper for MPO.identity)
`ising_hamiltonian`(n_sites[, J, h, device, dtype])	Transverse-field Ising Hamiltonian: H = -J Σᵢ ZᵢZᵢ₊₁ - h Σᵢ Xᵢ
`local_operator`(op, site, n_sites[, device, ...])	Create MPO for a single-site operator at specified site
`maxcut_hamiltonian`(edges[, weights, ...])	MaxCut QAOA Hamiltonian for graph optimization: H = Σ_{(i,j)∈E} w_{ij} (1 - ZᵢZⱼ) / 2
`molecular_hamiltonian_from_specs`([molecule, ...])	Build molecular Hamiltonian from molecular specifications using PySCF.
`sum_local_operators`(n_sites, local_ops[, ...])	Create MPO for sum of local operators: Σᵢ Oᵢ

static identity_mpo(n_sites, device='cuda', dtype=torch.complex64)[source]#

Create identity MPO (wrapper for MPO.identity)

static local_operator(op, site, n_sites, device='cuda', dtype=torch.complex64)[source]#

Create MPO for a single-site operator at specified site

Args:: op: 2x2 operator matrix site: Site index (0-indexed) n_sites: Total number of sites device: torch device dtype: data type
Returns:: MPO representing I ⊗ … ⊗ op ⊗ … ⊗ I

static sum_local_operators(n_sites, local_ops, device='cuda', dtype=torch.complex64)[source]#

Create MPO for sum of local operators: Σᵢ Oᵢ

Args:: n_sites: Total number of sites local_ops: List of (site, operator) tuples device: torch device dtype: data type
Returns:: MPO representing sum of all operators
Example:: # Build total magnetization Mz = Σ Zᵢ Z = torch.tensor([[1, 0], [0, -1]]) ops = [(i, Z) for i in range(n_sites)] Mz = MPOBuilder.sum_local_operators(n_sites, ops)

static ising_hamiltonian(n_sites, J=1.0, h=0.5, device='cuda', dtype=torch.complex64)[source]#

Transverse-field Ising Hamiltonian: H = -J Σᵢ ZᵢZᵢ₊₁ - h Σᵢ Xᵢ

Args:: n_sites: Number of sites J: Coupling strength h: Transverse field

Virtual bond structure (D=3): - 0→0: identity track - 0→1: emit Z (open ZZ term) - 1→2: close with -J Z - 2→2: identity tail - 0→2: local field -h X

static heisenberg_hamiltonian(n_sites, Jx=1.0, Jy=1.0, Jz=1.0, device='cuda', dtype=torch.complex64)[source]#

Heisenberg Hamiltonian: H = Σᵢ (Jₓ XᵢXᵢ₊₁ + Jᵧ YᵢYᵢ₊₁ + Jᵧ ZᵢZᵢ₊₁)

static maxcut_hamiltonian(edges, weights=None, n_sites=None, device='cuda', dtype=torch.complex64)[source]#

MaxCut QAOA Hamiltonian for graph optimization: H = Σ_{(i,j)∈E} w_{ij} (1 - ZᵢZⱼ) / 2

This Hamiltonian encodes the MaxCut problem where we want to maximize the number of edges between two sets (minimize edges within sets).

Args:

edges: List of (i, j) tuples representing graph edges weights: Optional edge weights (default: all 1.0) n_sites: Number of nodes (inferred from edges if not provided) device: ‘cuda’ or ‘cpu’ dtype: torch dtype

Returns:

MPO representation of MaxCut Hamiltonian

Example:

```python # Triangle graph (nodes 0, 1, 2) edges = [(0,1), (1,2), (0,2)] H = MPOBuilder.maxcut_hamiltonian(edges, device=’cuda’)

# Use with QAOA from atlas_q import get_vqe_qaoa qaoa = get_vqe_qaoa() config = qaoa[‘QAOAConfig’](p=2, max_iter=100) solver = qaoa[‘QAOA’](H, config) energy, params = solver.run() ```

static from_local_terms(n_sites, local_terms, device='cuda', dtype=torch.complex64)[source]#

Build Hamiltonian from local interaction terms

Args:

n_sites: Number of sites local_terms: List of (site1, site2, operator) tuples

operator should be a 4x4 tensor for two-site interactions

device: torch device dtype: data type

Returns:

MPO representing sum of local terms

Example:

# Build ZZ interaction Hamiltonian Z = torch.tensor([[1, 0], [0, -1]], dtype=torch.complex64) local_terms = [(i, i+1, torch.kron(Z, Z)) for i in range(n_sites-1)] H = MPOBuilder.from_local_terms(n_sites, local_terms)

static molecular_hamiltonian_from_specs(molecule='H2', basis='sto-3g', charge=0, spin=0, mapping='jordan_wigner', device='cuda', dtype=torch.complex64)[source]#

Build molecular Hamiltonian from molecular specifications using PySCF.

This function computes the electronic Hamiltonian for a molecule and converts it to an MPO suitable for VQE or other quantum algorithms.

Args:

molecule: Molecular formula or geometry string: Examples: ‘H2’, ‘LiH’, ‘H2O’, or geometry string

basis: Gaussian basis set (sto-3g, 6-31g, cc-pvdz, etc.) charge: Total molecular charge spin: Spin multiplicity (2S, where S is total spin) mapping: Fermion-to-qubit mapping (‘jordan_wigner’, ‘bravyi_kitaev’, ‘parity’) device: ‘cuda’ or ‘cpu’ dtype: torch dtype

Returns:

MPO representation of molecular Hamiltonian

Example:

```python from atlas_q import get_mpo_ops, get_vqe_qaoa

# H2 molecule mpo_mod = get_mpo_ops() H = mpo_mod[‘MPOBuilder’].molecular_hamiltonian_from_specs(

molecule=’H2’, basis=’sto-3g’, device=’cuda’

)

# Run VQE vqe_mod = get_vqe_qaoa() config = vqe_mod[‘VQEConfig’](n_layers=2, max_iter=100) vqe = vqe_mod[‘VQE’](H, config) energy, params = vqe.run() print(f”Ground state energy: {energy:.6f} Ha”) ```

Note:

Requires pyscf package: pip install pyscf

atlas_q.mpo_ops.apply_mpo_to_mps(mpo, mps, chi_max=128, eps=1e-08)[source]#

Apply MPO to MPS: |ψ’⟩ = O |ψ⟩

Uses zipper/zip-up algorithm for efficient contraction

Args:: mpo: Matrix product operator mps: Matrix product state (AdaptiveMPS) chi_max: Maximum bond dimension after compression eps: Truncation tolerance
Returns:: New MPS after applying MPO

atlas_q.mpo_ops.expectation_value(mpo, mps, use_gpu_optimized=True)[source]#

Compute ⟨ψ|O|ψ⟩ where O is an MPO and |ψ⟩ is an MPS

Args:: mpo: Matrix Product Operator mps: Matrix Product State use_gpu_optimized: Use GPU-optimized contractions (torch.compile)
Returns:: Complex expectation value

atlas_q.mpo_ops.correlation_function(op1, site1, op2, site2, mps)[source]#

Compute two-point correlation function: ⟨ψ| O₁(site1) O₂(site2) |ψ⟩

Args:: op1: First operator (2×2) site1: First site op2: Second operator (2×2) site2: Second site mps: MPS state
Returns:: Correlation ⟨O₁ O₂⟩

atlas_q.mpo_ops.pauli_string_to_mpo(pauli_string, device='cuda', dtype=torch.complex128)[source]#

Convert a Pauli string to an MPO.

Args:

pauli_string: String like “IXYZ” representing I⊗X⊗Y⊗Z device: torch device dtype: torch dtype

Returns:

MPO representation of the Pauli operator

Example:

>>> mpo = pauli_string_to_mpo("ZZII", device="cuda")  # Z⊗Z⊗I⊗I

atlas_q.mpo_ops.apply_pauli_exp_to_mps(mps, pauli_string, coeff, theta, chi_max=128)[source]#

Apply exp(theta * coeff * P) to MPS in-place, where P is a Pauli string.

IMPORTANT: This implements U(θ) = exp(θ * coeff * P) - If coeff = i·a (purely imaginary from anti-Hermitian UCCSD), this gives exp(i*(aθ)*P) → unitary - If coeff = a (real), this gives exp(aθ*P) → must include i factor for unitarity

For unitary Pauli exponentials: exp(i * λ * P) = cos(λ) I + i sin(λ) P (P² = I)

Implementation: Builds a single MPO M = cos(λ)I + i·sin(λ)P and applies it once. This avoids MPS summation errors that would break unitarity.

Args:: mps: AdaptiveMPS to modify in-place pauli_string: Pauli string like “IXYZ” coeff: Complex coefficient from generator (often purely imaginary for UCCSD) theta: Variational parameter (real) chi_max: Maximum bond dimension after compression

Overview#

The mpo_ops module provides Matrix Product Operator (MPO) operations for representing quantum operators as tensor networks. MPOs enable efficient computation of expectation values, operator application, and correlation functions for systems where operators have local structure.

Key Features#

Efficient operator representation: O(n D²) vs. O(4ⁿ) for dense matrices
Hamiltonian builders: Ising, Heisenberg, molecular chemistry, MaxCut
Expectation values: \(\langle\psi|\hat{O}|\psi\rangle\) in O(n χ² D²) time
Operator application: \(|\psi'\rangle = \hat{O}|\psi\rangle\)
Correlation functions: Two-point and multi-point correlators
GPU acceleration: CUDA-enabled tensor contractions

Why MPOs?#

Dense operator matrices scale as \(d^{2n}\) (16 GB for 15 qubits). MPOs exploit locality:

Local Hamiltonians: Many-body Hamiltonians with nearest-neighbor or short-range interactions can be represented exactly with small bond dimension D (typically D ≤ 10).

Examples:

Ising model: D = 3
Heisenberg model: D = 5
Molecular Hamiltonians: D = O(N_orbitals)

Mathematical Background#

Matrix Product Operator Representation#

An n-site operator \(\hat{O}\) is decomposed as:

\[\hat{O} = \sum_{s_1,s_1',\ldots,s_n,s_n'} W^{[1]}_{s_1 s_1'} W^{[2]}_{s_2 s_2'} \cdots W^{[n]}_{s_n s_n'} |s_1 \ldots s_n\rangle \langle s_1' \ldots s_n'|\]

where each tensor \(W^{[i]}\) has shape \([\chi_{i-1}, d, d, \chi_i]\):

\(d = 2\) for qubits
\(\chi_i\) is the MPO bond dimension (controls accuracy)

Storage: O(n D² d²) vs. O(d^{2n}) for dense operator

Hamiltonian as MPO#

Consider the transverse-field Ising model:

\[\hat{H} = -J \sum_{i=1}^{n-1} Z_i Z_{i+1} - h \sum_{i=1}^n X_i\]

This can be exactly represented as an MPO with bond dimension D = 3:

\[\begin{split}W^{[i]} = \begin{pmatrix} I & 0 & 0 \\ Z & 0 & 0 \\ -h X & -J Z & I \end{pmatrix}\end{split}\]

for internal sites, with boundary conditions for first and last sites.

Expectation Value#

Computing \(\langle\psi|\hat{O}|\psi\rangle\) for MPS \(|\psi\rangle\) and MPO \(\hat{O}\):

\[\langle\psi|\hat{O}|\psi\rangle = \text{Contract}(A^*, W, A)\]

Complexity: O(n χ² D²) where χ is MPS bond dimension, D is MPO bond dimension

Classes#

`MPO`	Matrix Product Operator
`MPOBuilder`	Helper class to build common MPOs

Functions#

`expectation_value`	Compute ⟨ψ\|O\|ψ⟩ where O is an MPO and \|ψ⟩ is an MPS
`apply_mpo_to_mps`	Apply MPO to MPS: \|ψ'⟩ = O \|ψ⟩
`correlation_function`	Compute two-point correlation function: ⟨ψ\| O₁(site1) O₂(site2) \|ψ⟩

MPO#

class atlas_q.mpo_ops.MPO(tensors, n_sites)[source]#

Bases: object

Matrix Product Operator

Represents an operator as a chain of 4-tensors: O = Σ W[0]_{s₀s₀’} W[1]_{s₁s₁’} … W[n-1]_{sₙ₋₁sₙ₋₁’}

Each tensor W[i] has shape [χ_L, d, d, χ_R] where: - χ_L, χ_R: left and right bond dimensions - d: physical dimension (2 for qubits)

Methods

`from_local_ops`(ops[, device])	Create MPO from list of local operators (one per site)
`from_operators`(ops[, device])	Alias for from_local_ops
`identity`(n_sites[, device, dtype])	Create identity MPO

Matrix Product Operator representation for quantum operators.

Stores operator as chain of 4-tensors with efficient contraction algorithms. Supports operator addition, scaling, and composition.

Constructor:

from atlas_q.mpo_ops import MPO
import torch

# Create MPO from tensor list
tensors = [...]  # List of shape [chi_L, d, d, chi_R]
mpo = MPO(tensors, device='cuda')

Parameters:

tensors (list[torch.Tensor]): List of MPO tensors
device (str): ‘cuda’ or ‘cpu’

Storage: Approximately \(16 n D^2 d^2\) bytes for complex64

Methods

`identity`(n_sites[, device, dtype])	Create identity MPO
`from_local_ops`(ops[, device])	Create MPO from list of local operators (one per site)
`from_operators`(ops[, device])	Alias for from_local_ops
`__add__`(other)	Add two MPOs by expanding bond dimensions

Class Methods:

classmethod identity(n_sites, device='cuda')[source]#

Create identity operator MPO.

Parameters:

n_sites (int) – Number of sites
device (str) – Device

Returns:

Identity MPO

Return type:

MPO

Example:

from atlas_q.mpo_ops import MPO

I = MPO.identity(n_sites=10, device='cuda')

classmethod from_local_ops(operators, device='cuda')[source]#

Create MPO from list of local (single-site) operators.

Parameters:

operators (list) – List of 2×2 operator matrices
device (str) – Device

Returns:

MPO representing tensor product

Return type:

MPO

Example:

import torch
from atlas_q.mpo_ops import MPO

X = torch.tensor([[0, 1], [1, 0]], dtype=torch.complex64)
Z = torch.tensor([[1, 0], [0, -1]], dtype=torch.complex64)
I = torch.eye(2, dtype=torch.complex64)

# Create X ⊗ Z ⊗ I
mpo = MPO.from_local_ops([X, Z, I], device='cuda')

classmethod from_operators(op_strings, coeffs, n_sites, device='cuda')[source]#

Create MPO from sum of operator strings.

Parameters:

op_strings (list) – List of operator strings (e.g., [‘XXI’, ‘ZZI’])
coeffs (list) – Corresponding coefficients
n_sites (int) – Number of sites
device (str) – Device

Returns:

MPO

Return type:

MPO

Example:

# H = 0.5 X₀X₁ + 0.3 Z₀Z₁
mpo = MPO.from_operators(
    op_strings=['XXI', 'ZZI', 'IIX'],
    coeffs=[0.5, 0.3, 0.2],
    n_sites=3,
    device='cuda'
)

Operator Arithmetic:

__add__(other)[source]#

Add two MPOs: \(\hat{O}_1 + \hat{O}_2\)

Parameters:: other (MPO) – Another MPO
Returns:: Sum MPO
Return type:: MPO

Bond dimension grows: \(D_{\text{sum}} = D_1 + D_2\)

Example:

H1 = MPOBuilder.ising_hamiltonian(10, J=1.0, h=0.0)
H2 = MPOBuilder.ising_hamiltonian(10, J=0.0, h=0.5)
H_total = H1 + H2  # Combined Hamiltonian

__mul__(scalar)#

Scalar multiplication: \(c \hat{O}\)

Parameters:: scalar (float) – Scalar coefficient
Returns:: Scaled MPO
Return type:: MPO

tensors: List[Tensor]#

n_sites: int#

static identity(n_sites, device='cuda', dtype=torch.complex64)[source]#

Create identity MPO

static from_local_ops(ops, device='cuda')[source]#

Create MPO from list of local operators (one per site)

Args:: ops: List of 2×2 operators for each site

static from_operators(ops, device='cuda')[source]#

Alias for from_local_ops

__add__(other)[source]#

Add two MPOs by expanding bond dimensions

For A + B, we create a new MPO with bond dimension χ_A + χ_B that represents the sum of the two operators.

MPOBuilder#

class atlas_q.mpo_ops.MPOBuilder[source]#

Bases: object

Helper class to build common MPOs

Methods

`from_local_terms`(n_sites, local_terms[, ...])	Build Hamiltonian from local interaction terms
`heisenberg_hamiltonian`(n_sites[, Jx, Jy, ...])	Heisenberg Hamiltonian: H = Σᵢ (Jₓ XᵢXᵢ₊₁ + Jᵧ YᵢYᵢ₊₁ + Jᵧ ZᵢZᵢ₊₁)
`identity_mpo`(n_sites[, device, dtype])	Create identity MPO (wrapper for MPO.identity)
`ising_hamiltonian`(n_sites[, J, h, device, dtype])	Transverse-field Ising Hamiltonian: H = -J Σᵢ ZᵢZᵢ₊₁ - h Σᵢ Xᵢ
`local_operator`(op, site, n_sites[, device, ...])	Create MPO for a single-site operator at specified site
`maxcut_hamiltonian`(edges[, weights, ...])	MaxCut QAOA Hamiltonian for graph optimization: H = Σ_{(i,j)∈E} w_{ij} (1 - ZᵢZⱼ) / 2
`molecular_hamiltonian_from_specs`([molecule, ...])	Build molecular Hamiltonian from molecular specifications using PySCF.
`sum_local_operators`(n_sites, local_ops[, ...])	Create MPO for sum of local operators: Σᵢ Oᵢ

Factory class for constructing common quantum Hamiltonians as MPOs.

Provides pre-built Hamiltonians for:

Spin models (Ising, Heisenberg, XY)
Molecular chemistry (via PySCF)
Combinatorial optimization (MaxCut, QAOA)

Methods

`identity_mpo`(n_sites[, device, dtype])	Create identity MPO (wrapper for MPO.identity)
`ising_hamiltonian`(n_sites[, J, h, device, dtype])	Transverse-field Ising Hamiltonian: H = -J Σᵢ ZᵢZᵢ₊₁ - h Σᵢ Xᵢ
`heisenberg_hamiltonian`(n_sites[, Jx, Jy, ...])	Heisenberg Hamiltonian: H = Σᵢ (Jₓ XᵢXᵢ₊₁ + Jᵧ YᵢYᵢ₊₁ + Jᵧ ZᵢZᵢ₊₁)
`molecular_hamiltonian_from_specs`([molecule, ...])	Build molecular Hamiltonian from molecular specifications using PySCF.
`maxcut_hamiltonian`(edges[, weights, ...])	MaxCut QAOA Hamiltonian for graph optimization: H = Σ_{(i,j)∈E} w_{ij} (1 - ZᵢZⱼ) / 2

static identity_mpo(n_sites, device='cuda')[source]#

Create identity operator.

Parameters:

n_sites (int) – Number of sites
device (str) – Device

Returns:

Identity MPO

Return type:

MPO

static ising_hamiltonian(n_sites, J=1.0, h=0.5, periodic=False, device='cuda')[source]#

Create transverse-field Ising Hamiltonian.

\[\hat{H} = -J \sum_{i=1}^{n-1} Z_i Z_{i+1} - h \sum_{i=1}^n X_i\]

Parameters:

n_sites (int) – Number of sites
J (float) – Coupling strength
h (float) – Transverse field strength
periodic (bool) – Periodic boundary conditions
device (str) – Device

Returns:

Ising MPO (bond dimension D=3)

Return type:

MPO

Example:

from atlas_q.mpo_ops import MPOBuilder

# Ferromagnetic (J > 0)
H = MPOBuilder.ising_hamiltonian(n_sites=10, J=1.0, h=0.5, device='cuda')

# Antiferromagnetic (J < 0)
H = MPOBuilder.ising_hamiltonian(n_sites=10, J=-1.0, h=0.5, device='cuda')

static heisenberg_hamiltonian(n_sites, Jx=1.0, Jy=1.0, Jz=1.0, device='cuda')[source]#

Create Heisenberg Hamiltonian.

\[\hat{H} = \sum_{i=1}^{n-1} (J_x X_i X_{i+1} + J_y Y_i Y_{i+1} + J_z Z_i Z_{i+1})\]

Parameters:

n_sites (int) – Number of sites
Jx (float) – X coupling
Jy (float) – Y coupling
Jz (float) – Z coupling
device (str) – Device

Returns:

Heisenberg MPO (bond dimension D=5)

Return type:

MPO

Special cases:

Jx = Jy = Jz: Isotropic Heisenberg
Jx = Jy, Jz = 0: XY model
Jx = Jy = 0: Ising model

Example:

# Isotropic Heisenberg (SU(2) symmetric)
H = MPOBuilder.heisenberg_hamiltonian(n_sites=10, Jx=1.0, Jy=1.0, Jz=1.0)

# XXZ model
H = MPOBuilder.heisenberg_hamiltonian(n_sites=10, Jx=1.0, Jy=1.0, Jz=2.0)

static xy_hamiltonian(n_sites, Jx=1.0, Jy=1.0, device='cuda')#

Create XY model Hamiltonian.

\[\hat{H} = \sum_{i=1}^{n-1} (J_x X_i X_{i+1} + J_y Y_i Y_{i+1})\]

Parameters:

n_sites (int) – Number of sites
Jx (float) – X coupling
Jy (float) – Y coupling
device (str) – Device

Returns:

XY MPO (bond dimension D=3)

Return type:

MPO

static molecular_hamiltonian_from_specs(molecule, basis='sto-3g', geometry=None, device='cuda')[source]#

Create molecular Hamiltonian using PySCF.

Parameters:

molecule (str) – Molecule name (‘H2’, ‘LiH’, ‘H2O’, etc.)
basis (str) – Basis set (‘sto-3g’, ‘6-31g’, etc.)
geometry (list) – Optional custom geometry as [(atom, coords), …]
device (str) – Device

Returns:

Molecular Hamiltonian MPO

Return type:

MPO

Requires: PySCF installed (pip install pyscf)

Example:

# H2 molecule at equilibrium
H = MPOBuilder.molecular_hamiltonian_from_specs(
    molecule='H2',
    basis='sto-3g',
    device='cuda'
)

# LiH with custom geometry
H = MPOBuilder.molecular_hamiltonian_from_specs(
    molecule='LiH',
    basis='6-31g',
    geometry=[('Li', (0, 0, 0)), ('H', (1.6, 0, 0))],  # Angstroms
    device='cuda'
)

static maxcut_hamiltonian(graph, device='cuda')[source]#

Create MaxCut Hamiltonian for graph.

\[\hat{H} = -\frac{1}{2} \sum_{(i,j) \in E} (I - Z_i Z_j)\]

Parameters:

graph (networkx.Graph) – Graph for MaxCut problem
device (str) – Device

Returns:

MaxCut MPO

Return type:

MPO

Example:

import networkx as nx
from atlas_q.mpo_ops import MPOBuilder

# 4-node cycle graph
G = nx.Graph()
G.add_edges_from([(0, 1), (1, 2), (2, 3), (3, 0)])

H = MPOBuilder.maxcut_hamiltonian(G, device='cuda')

# Use with QAOA
from atlas_q.vqe_qaoa import QAOA
qaoa = QAOA(hamiltonian=H, p=3)
result = qaoa.run()

static identity_mpo(n_sites, device='cuda', dtype=torch.complex64)[source]#

Create identity MPO (wrapper for MPO.identity)

static local_operator(op, site, n_sites, device='cuda', dtype=torch.complex64)[source]#

Create MPO for a single-site operator at specified site

Args:: op: 2x2 operator matrix site: Site index (0-indexed) n_sites: Total number of sites device: torch device dtype: data type
Returns:: MPO representing I ⊗ … ⊗ op ⊗ … ⊗ I

static sum_local_operators(n_sites, local_ops, device='cuda', dtype=torch.complex64)[source]#

Create MPO for sum of local operators: Σᵢ Oᵢ

Args:: n_sites: Total number of sites local_ops: List of (site, operator) tuples device: torch device dtype: data type
Returns:: MPO representing sum of all operators
Example:: # Build total magnetization Mz = Σ Zᵢ Z = torch.tensor([[1, 0], [0, -1]]) ops = [(i, Z) for i in range(n_sites)] Mz = MPOBuilder.sum_local_operators(n_sites, ops)

static ising_hamiltonian(n_sites, J=1.0, h=0.5, device='cuda', dtype=torch.complex64)[source]#

Transverse-field Ising Hamiltonian: H = -J Σᵢ ZᵢZᵢ₊₁ - h Σᵢ Xᵢ

Args:: n_sites: Number of sites J: Coupling strength h: Transverse field

Virtual bond structure (D=3): - 0→0: identity track - 0→1: emit Z (open ZZ term) - 1→2: close with -J Z - 2→2: identity tail - 0→2: local field -h X

static heisenberg_hamiltonian(n_sites, Jx=1.0, Jy=1.0, Jz=1.0, device='cuda', dtype=torch.complex64)[source]#

Heisenberg Hamiltonian: H = Σᵢ (Jₓ XᵢXᵢ₊₁ + Jᵧ YᵢYᵢ₊₁ + Jᵧ ZᵢZᵢ₊₁)

static maxcut_hamiltonian(edges, weights=None, n_sites=None, device='cuda', dtype=torch.complex64)[source]#

MaxCut QAOA Hamiltonian for graph optimization: H = Σ_{(i,j)∈E} w_{ij} (1 - ZᵢZⱼ) / 2

This Hamiltonian encodes the MaxCut problem where we want to maximize the number of edges between two sets (minimize edges within sets).

Args:

Returns:

MPO representation of MaxCut Hamiltonian

Example:

```python # Triangle graph (nodes 0, 1, 2) edges = [(0,1), (1,2), (0,2)] H = MPOBuilder.maxcut_hamiltonian(edges, device=’cuda’)

# Use with QAOA from atlas_q import get_vqe_qaoa qaoa = get_vqe_qaoa() config = qaoa[‘QAOAConfig’](p=2, max_iter=100) solver = qaoa[‘QAOA’](H, config) energy, params = solver.run() ```

static from_local_terms(n_sites, local_terms, device='cuda', dtype=torch.complex64)[source]#

Build Hamiltonian from local interaction terms

Args:

n_sites: Number of sites local_terms: List of (site1, site2, operator) tuples

operator should be a 4x4 tensor for two-site interactions

device: torch device dtype: data type

Returns:

MPO representing sum of local terms

Example:

static molecular_hamiltonian_from_specs(molecule='H2', basis='sto-3g', charge=0, spin=0, mapping='jordan_wigner', device='cuda', dtype=torch.complex64)[source]#

Build molecular Hamiltonian from molecular specifications using PySCF.

This function computes the electronic Hamiltonian for a molecule and converts it to an MPO suitable for VQE or other quantum algorithms.

Args:

molecule: Molecular formula or geometry string: Examples: ‘H2’, ‘LiH’, ‘H2O’, or geometry string

Returns:

MPO representation of molecular Hamiltonian

Example:

```python from atlas_q import get_mpo_ops, get_vqe_qaoa

# H2 molecule mpo_mod = get_mpo_ops() H = mpo_mod[‘MPOBuilder’].molecular_hamiltonian_from_specs(

molecule=’H2’, basis=’sto-3g’, device=’cuda’

)

Note:

Requires pyscf package: pip install pyscf

Functions#

atlas_q.mpo_ops.expectation_value(mpo, mps, use_gpu_optimized=True)[source]#

Compute ⟨ψ|O|ψ⟩ where O is an MPO and |ψ⟩ is an MPS

Args:: mpo: Matrix Product Operator mps: Matrix Product State use_gpu_optimized: Use GPU-optimized contractions (torch.compile)
Returns:: Complex expectation value

Compute expectation value \(\langle\psi|\hat{O}|\psi\rangle\).

Signature:

expectation_value(mps, mpo) -> complex

Parameters:

mps (AdaptiveMPS): State as Matrix Product State
mpo (MPO): Operator as Matrix Product Operator

Returns:

complex: Expectation value

Complexity: O(n χ² D²) where n is number of sites, χ is MPS bond dimension, D is MPO bond dimension

Example:

from atlas_q.mpo_ops import MPOBuilder, expectation_value
from atlas_q.adaptive_mps import AdaptiveMPS

# Create Hamiltonian
H = MPOBuilder.ising_hamiltonian(n_sites=10, J=1.0, h=0.5, device='cuda')

# Create state
mps = AdaptiveMPS(num_qubits=10, bond_dim=16, device='cuda')
# ... prepare state ...

# Compute energy
energy = expectation_value(mps, H)
print(f"Energy: {energy.real:.6f} Ha")

atlas_q.mpo_ops.apply_mpo_to_mps(mpo, mps, chi_max=128, eps=1e-08)[source]#

Apply MPO to MPS: |ψ’⟩ = O |ψ⟩

Uses zipper/zip-up algorithm for efficient contraction

Args:: mpo: Matrix product operator mps: Matrix product state (AdaptiveMPS) chi_max: Maximum bond dimension after compression eps: Truncation tolerance
Returns:: New MPS after applying MPO

Apply MPO to MPS: \(|\psi'\rangle = \hat{O}|\psi\rangle\).

Signature:

apply_mpo_to_mps(mps, mpo, chi_max=None) -> AdaptiveMPS

Parameters:

mps (AdaptiveMPS): Input state
mpo (MPO): Operator to apply
chi_max (int, optional): Maximum bond dimension for result

Returns:

AdaptiveMPS: Resulting state

Complexity: O(n χ² D²) with bond dimension growth χ’ ≤ χ·D

Example:

from atlas_q.mpo_ops import MPO, apply_mpo_to_mps
from atlas_q.adaptive_mps import AdaptiveMPS

# Create state |ψ⟩
mps = AdaptiveMPS(num_qubits=10, bond_dim=16, device='cuda')

# Create operator (e.g., Hadamard on all qubits)
H_tensor = torch.tensor([[1, 1], [1, -1]], dtype=torch.complex64) / (2**0.5)
H_mpo = MPO.from_local_ops([H_tensor] * 10, device='cuda')

# Apply: |ψ'⟩ = H|ψ⟩
mps_prime = apply_mpo_to_mps(mps, H_mpo, chi_max=64)

atlas_q.mpo_ops.correlation_function(op1, site1, op2, site2, mps)[source]#

Compute two-point correlation function: ⟨ψ| O₁(site1) O₂(site2) |ψ⟩

Args:: op1: First operator (2×2) site1: First site op2: Second operator (2×2) site2: Second site mps: MPS state
Returns:: Correlation ⟨O₁ O₂⟩

Compute two-point correlation function \(\langle\psi|\hat{O}_i \hat{O}_j|\psi\rangle\).

Signature:

correlation_function(mps, op_i, op_j, site_i, site_j) -> complex

Parameters:

mps (AdaptiveMPS): State
op_i (torch.Tensor): Operator at site i (2×2 matrix)
op_j (torch.Tensor): Operator at site j (2×2 matrix)
site_i (int): First site index
site_j (int): Second site index

Returns:

complex: Correlation value

Complexity: O(|j-i| χ²)

Example:

from atlas_q.mpo_ops import correlation_function
from atlas_q.adaptive_mps import AdaptiveMPS
import torch

mps = AdaptiveMPS(num_qubits=20, bond_dim=32, device='cuda')
# ... prepare state ...

# Compute <Z_5 Z_10>
Z = torch.tensor([[1, 0], [0, -1]], dtype=torch.complex64, device='cuda')
corr = correlation_function(mps, Z, Z, site_i=5, site_j=10)
print(f"<Z_5 Z_10> = {corr.real:.6f}")

# Spin-spin correlation function
for j in range(20):
    corr_j = correlation_function(mps, Z, Z, site_i=0, site_j=j)
    print(f"<Z_0 Z_{j}> = {corr_j.real:.4f}")

Performance Characteristics#

Computational Complexity#

Operation	Complexity
Expectation value	O(n χ² D²)
Apply MPO to MPS	O(n χ² D²)
Correlation function	O(\|i-j\| χ²)
MPO addition	O(n D₁ D₂)

where:

n: number of sites
χ: MPS bond dimension
D: MPO bond dimension

Memory Usage#

Structure	Memory
MPO storage	O(n D² d²)
MPS storage	O(n χ² d)
Expectation value work	O(χ² D)

Example: 50 sites, χ=128, D=5, d=2

MPO: 50 × 5² × 2² × 8 bytes = 40 KB
MPS: 50 × 128² × 2 × 8 bytes = 13 MB
Working memory: 128² × 5 × 8 bytes = 0.6 MB

Total: ~14 MB (dominated by MPS)

Benchmark Results#

From scripts/benchmarks/validate_all_features.py:

# Expectation value: 50 sites, χ=128, D=3 (Ising)
Time: 0.08 sec (GPU), 0.35 sec (CPU)
Speedup: 4.4×

# Apply MPO: 50 sites, χ=64, D=3
Time: 0.15 sec (GPU), 0.62 sec (CPU)
Speedup: 4.1×

# Correlation function: 50 sites, distance=25, χ=128
Time: 0.05 sec (GPU), 0.18 sec (CPU)
Speedup: 3.6×

Examples#

Ising Hamiltonian#

from atlas_q.mpo_ops import MPOBuilder, expectation_value
from atlas_q.adaptive_mps import AdaptiveMPS

# Build Hamiltonian: H = -J Σ Z_i Z_{i+1} - h Σ X_i
H = MPOBuilder.ising_hamiltonian(n_sites=10, J=1.0, h=0.5, device='cuda')

# Create ground state (all |0⟩ for J > 0, h = 0)
mps = AdaptiveMPS(num_qubits=10, bond_dim=8, device='cuda')

# Compute energy
energy = expectation_value(mps, H)
print(f"Energy: {energy.real:.6f}")

Heisenberg Hamiltonian#

# H = Σ (Jx X_i X_{i+1} + Jy Y_i Y_{i+1} + Jz Z_i Z_{i+1})
H = MPOBuilder.heisenberg_hamiltonian(
    n_sites=10,
    Jx=1.0,
    Jy=1.0,
    Jz=1.0,
    device='cuda'
)

# Use with TDVP for ground state
from atlas_q.tdvp import TDVP
mps = AdaptiveMPS(num_qubits=10, bond_dim=32, device='cuda')
tdvp = TDVP(mps=mps, hamiltonian=H, dt=0.1j, method='two_site')
energy = tdvp.run(n_steps=100)

Molecular Hamiltonian (H2)#

# Build H2 Hamiltonian
H = MPOBuilder.molecular_hamiltonian_from_specs(
    molecule='H2',
    basis='sto-3g',
    device='cuda'
)

# VQE for ground state energy
from atlas_q.vqe_qaoa import VQE
mps = AdaptiveMPS(num_qubits=4, bond_dim=16, device='cuda')  # H2 has 4 qubits
vqe = VQE(mps=mps, hamiltonian=H, optimizer='COBYLA')
result = vqe.run()
print(f"Ground state energy: {result['energy']:.6f} Ha")

MaxCut Hamiltonian#

import networkx as nx
from atlas_q.mpo_ops import MPOBuilder
from atlas_q.vqe_qaoa import QAOA

# Create graph
G = nx.cycle_graph(4)  # 4-node cycle

# Build MaxCut Hamiltonian
H = MPOBuilder.maxcut_hamiltonian(G, device='cuda')

# Solve with QAOA
from atlas_q.adaptive_mps import AdaptiveMPS
mps = AdaptiveMPS(num_qubits=4, bond_dim=8, device='cuda')
qaoa = QAOA(mps=mps, hamiltonian=H, p=3)
result = qaoa.run()

print(f"MaxCut value: {result['max_cut']}")
print(f"Best partition: {result['partition']}")

Custom Hamiltonian#

import torch
from atlas_q.mpo_ops import MPO

# Define local operators
X = torch.tensor([[0, 1], [1, 0]], dtype=torch.complex64, device='cuda')
Y = torch.tensor([[0, -1j], [1j, 0]], dtype=torch.complex64, device='cuda')
Z = torch.tensor([[1, 0], [0, -1]], dtype=torch.complex64, device='cuda')
I = torch.eye(2, dtype=torch.complex64, device='cuda')

# Build H = 0.5 X₀X₁ + 0.3 Y₀Y₁ + 0.2 Z₀Z₁
H = MPO.from_operators(
    op_strings=['XXI', 'YYI', 'ZZI'],
    coeffs=[0.5, 0.3, 0.2],
    n_sites=3,
    device='cuda'
)

MPO Addition#

# H_total = H_kinetic + H_potential
H_kinetic = MPOBuilder.ising_hamiltonian(10, J=1.0, h=0.0)
H_potential = MPOBuilder.ising_hamiltonian(10, J=0.0, h=0.5)
H_total = H_kinetic + H_potential

# Bond dimension increases: D_total = D_kin + D_pot = 3 + 3 = 6

Correlation Function Analysis#

from atlas_q.mpo_ops import correlation_function
import matplotlib.pyplot as plt

n_sites = 50
mps = AdaptiveMPS(num_qubits=n_sites, bond_dim=32, device='cuda')
# ... prepare state (e.g., ground state of Heisenberg) ...

Z = torch.tensor([[1, 0], [0, -1]], dtype=torch.complex64, device='cuda')

# Compute spin-spin correlation
correlations = []
for j in range(n_sites):
    corr = correlation_function(mps, Z, Z, site_i=0, site_j=j)
    correlations.append(corr.real)

# Plot
plt.plot(range(n_sites), correlations)
plt.xlabel('Distance')
plt.ylabel('<Z_0 Z_j>')
plt.title('Spin-Spin Correlation Function')
plt.show()

Use Cases#

When to Use MPOs#

Local Hamiltonians: Nearest-neighbor or short-range interactions
Expectation values: Need frequent energy calculations (VQE, TDVP)
Large systems: >20 qubits where dense matrices infeasible
Quantum chemistry: Molecular Hamiltonians with PySCF
Combinatorial optimization: MaxCut, graph problems with QAOA

MPO Bond Dimensions#

Common models and their bond dimensions:

Hamiltonian	Bond Dimension
Ising (transverse field)	D = 3
Heisenberg (isotropic)	D = 5
XY model	D = 3
Molecular (N orbitals)	D ~ 4N
MaxCut (graph)	D ~ max_degree

Rule of thumb: Sparse graphs and local interactions → small D

Cross-References#

References#

Key papers on MPOs:

Schollwöck, U. (2011). “The density-matrix renormalization group in the age of matrix product states.” Annals of Physics, 326(1), 96-192.
McCulloch, I. P. (2007). “From density-matrix renormalization group to matrix product states.” Journal of Statistical Mechanics: Theory and Experiment, 2007(10), P10014.
Crosswhite, G. M. & Bacon, D. (2008). “Finite automata for caching in matrix product algorithms.” Physical Review A, 78(1), 012356.

atlas_q.mpo_ops#

Overview#

Key Features#

Why MPOs?#

Mathematical Background#

Matrix Product Operator Representation#

Hamiltonian as MPO#

Expectation Value#

Classes#

Functions#

MPO#

MPOBuilder#

Functions#

Performance Characteristics#

Computational Complexity#

Memory Usage#

Benchmark Results#

Examples#

Ising Hamiltonian#

Heisenberg Hamiltonian#

Molecular Hamiltonian (H2)#

MaxCut Hamiltonian#

Custom Hamiltonian#

MPO Addition#

Correlation Function Analysis#

Use Cases#

When to Use MPOs#

MPO Bond Dimensions#

Cross-References#

See Also#

References#