Code Examples

#!/usr/bin/env python3
"""Complete CV cluster state design pipeline."""

import numpy as np
from cv_cluster import ClusterCompiler
from cv_cluster.graph import GraphState
from cv_cluster.analysis import MonteCarloAnalyzer
from cv_cluster.loss import LossSimulator
from cv_cluster.export import GDSExporter

# ============== STEP 1: DESIGN GRAPH ==============
print("Step 1: Defining graph topology...")

# Choose topology based on application
graph = GraphState.hexagonal(7)  # 7-mode hexagonal for MBQC

print(f"  Modes: {graph.num_modes}")
print(f"  Edges: {graph.num_edges}")
print(f"  Connectivity: {graph.avg_degree:.2f}")

# ============== STEP 2: CREATE COMPILER ==============
print("\nStep 2: Initializing compiler...")

compiler = ClusterCompiler(
    graph=graph,
    squeezing_dB=6.0,           # Input squeezed states
    mesh_architecture="clements"  # Interferometer type
)

print(f"  Architecture: Clements mesh")
print(f"  Phase parameters: {compiler.num_phases}")

# ============== STEP 3: OPTIMIZE PHASES ==============
print("\nStep 3: Optimizing phases...")

result = compiler.optimize_phases(
    method="L-BFGS-B",
    tol=1e-8,
    maxiter=1000
)

print(f"  Converged in {result.iterations} iterations")
print(f"  Avg nullifier variance: {result.avg_variance:.4f} SNU")
print(f"  Max nullifier variance: {result.max_variance:.4f} SNU")

# ============== STEP 4: MONTE CARLO ANALYSIS ==============
print("\nStep 4: Running Monte Carlo tolerance analysis...")

mc = MonteCarloAnalyzer(compiler)
mc_result = mc.run(
    num_trials=2000,
    phase_error_std=0.5,     # degrees
    bs_imbalance_std=0.02,   # 2% splitting error
    threshold=1.0             # SNU criterion
)

print(f"  Success rate: {mc_result.success_rate:.1%}")
print(f"  95th percentile variance: {mc_result.percentile_95:.4f} SNU")

# ============== STEP 5: LOSS ANALYSIS ==============
print("\nStep 5: Analyzing loss tolerance...")

loss_sim = LossSimulator(
    compiler,
    loss_per_layer_dB=0.3,
    num_layers=compiler.mesh_depth
)
loss_result = loss_sim.analyze()

print(f"  Total loss: {loss_result.total_loss_dB:.2f} dB")
print(f"  Entangled modes: {loss_result.entangled_count}/{graph.num_modes}")

# ============== STEP 6: EXPORT FOR FABRICATION ==============
print("\nStep 6: Exporting design...")

exporter = GDSExporter(compiler, result.optimal_phases)
exporter.export("hexagonal_cluster.gds")
exporter.export_phases("phase_settings.json")

print("  Generated: hexagonal_cluster.gds")
print("  Generated: phase_settings.json")
print("\n✓ Pipeline complete!")

"""Nullifier-based entanglement verification."""

import numpy as np
from cv_cluster import ClusterCompiler
from cv_cluster.graph import GraphState

def verify_cluster_state(compiler, threshold=1.0):
    """
    Verify entanglement via nullifier variances.

    For an ideal cluster state: δ_j = p_j - Σ_k A_jk x_k → 0
    Entanglement witness: ⟨δ_j²⟩ < 1 SNU for all j
    """
    n = compiler.graph.num_modes
    adj = compiler.graph.adjacency_matrix
    V = compiler.get_covariance_matrix()

    variances = []
    for j in range(n):
        # Build nullifier coefficient vector
        # c = [0...0, 0...0] with c[n+j] = 1 (p_j)
        # and c[k] = -A_jk for all k (x terms)
        c = np.zeros(2 * n)
        c[n + j] = 1  # p_j coefficient
        for k in range(n):
            c[k] = -adj[j, k]  # -A_jk * x_k

        # Variance: ⟨δ_j²⟩ = c^T V c
        var = c @ V @ c
        variances.append(var)

    variances = np.array(variances)

    # Check entanglement criterion
    all_entangled = np.all(variances < threshold)
    avg_var = np.mean(variances)
    max_var = np.max(variances)

    return {
        'variances': variances,
        'average': avg_var,
        'maximum': max_var,
        'all_entangled': all_entangled,
        'entangled_count': np.sum(variances < threshold)
    }

# Example usage
graph = GraphState.square(4)
compiler = ClusterCompiler(graph, squeezing_dB=6.0)
compiler.optimize_phases()

result = verify_cluster_state(compiler)
print(f"Mode variances (SNU):")
for i, v in enumerate(result['variances']):
    status = "✓ entangled" if v < 1 else "✗ separable"
    print(f"  δ_{i}²: {v:.4f} {status}")
print(f"\nCluster valid: {result['all_entangled']}")

"""Symplectic transformation of Gaussian states."""

import numpy as np

def beamsplitter_symplectic(theta, modes, n_total):
    """2N×2N symplectic matrix for a beamsplitter."""
    S = np.eye(2 * n_total)
    c, s = np.cos(theta), np.sin(theta)
    i, j = modes

    # x quadratures
    S[i, i] = c
    S[i, j] = s
    S[j, i] = -s
    S[j, j] = c

    # p quadratures
    S[n_total+i, n_total+i] = c
    S[n_total+i, n_total+j] = s
    S[n_total+j, n_total+i] = -s
    S[n_total+j, n_total+j] = c

    return S

def phase_shifter_symplectic(phi, mode, n_total):
    """2N×2N symplectic matrix for a phase shifter."""
    S = np.eye(2 * n_total)
    c, s = np.cos(phi), np.sin(phi)

    S[mode, mode] = c
    S[mode, n_total+mode] = s
    S[n_total+mode, mode] = -s
    S[n_total+mode, n_total+mode] = c

    return S

def clements_mesh(phases, n_modes):
    """Build full Clements mesh symplectic matrix."""
    S = np.eye(2 * n_modes)
    idx = 0

    # Triangular arrangement of beamsplitters
    for layer in range(n_modes - 1):
        start = layer % 2
        for i in range(start, n_modes - 1, 2):
            # Phase shifter
            S = phase_shifter_symplectic(phases[idx], i, n_modes) @ S
            idx += 1
            # Beamsplitter
            S = beamsplitter_symplectic(phases[idx], (i, i+1), n_modes) @ S
            idx += 1

    return S

def propagate_covariance(V_in, S):
    """Propagate covariance: V_out = S V_in S^T"""
    return S @ V_in @ S.T

# Example: 4-mode squeezed states through mesh
n = 4
r = 6.0 / (10 * np.log10(np.e) * 2)  # 6 dB squeezing

# Initial squeezed vacuum covariance
V0 = np.diag([np.exp(-2*r)/2]*n + [np.exp(2*r)/2]*n)

# Random phases
num_phases = 2 * (n * (n - 1))
phases = np.random.uniform(0, 2*np.pi, num_phases)

S = clements_mesh(phases, n)
V_out = propagate_covariance(V0, S)

print(f"Output covariance matrix shape: {V_out.shape}")
print(f"Symplecticity check: {np.allclose(np.linalg.det(S), 1.0)}")

"""CV GHZ-type state generation via cluster rail."""

from cv_cluster import ClusterCompiler
from cv_cluster.graph import GraphState
import numpy as np

def create_ghz_cluster(num_modes, squeezing_dB):
    """
    Create GHZ-rail cluster state.

    Linear graph: 0 - 1 - 2 - ... - (N-1)
    Nullifiers: p_0 - x_1, p_1 - x_0 - x_2, ..., p_{N-1} - x_{N-2}
    """
    graph = GraphState.ghz_rail(num_modes)

    compiler = ClusterCompiler(
        graph=graph,
        squeezing_dB=squeezing_dB,
        mesh_architecture="clements"
    )

    # Optimize for minimal nullifier variance
    result = compiler.optimize_phases(tol=1e-10)

    return compiler, result

def measure_ghz_correlations(compiler):
    """Measure X-X and P-P correlations characteristic of GHZ states."""
    V = compiler.get_covariance_matrix()
    n = compiler.graph.num_modes

    # Extract X-X correlations (upper-left block)
    V_xx = V[:n, :n]

    # Extract P-P correlations (lower-right block)
    V_pp = V[n:, n:]

    # Sum variance: Var(X_total) for GHZ
    x_sum_var = np.sum(V_xx)
    p_diff_var = V_pp[0,0] + V_pp[-1,-1] - 2*V_pp[0,-1]

    return {
        'x_correlation_matrix': V_xx,
        'p_correlation_matrix': V_pp,
        'x_sum_variance': x_sum_var,
        'p_diff_variance': p_diff_var
    }

# Generate 8-mode GHZ cluster
compiler, opt_result = create_ghz_cluster(8, squeezing_dB=8.0)

print("8-Mode CV GHZ Cluster State")
print("=" * 40)
print(f"Average nullifier variance: {opt_result.avg_variance:.4f} SNU")
print(f"Optimization iterations: {opt_result.iterations}")

correlations = measure_ghz_correlations(compiler)
print(f"\nGHZ correlations:")
print(f"  Var(ΣX): {correlations['x_sum_variance']:.4f}")
print(f"  Var(P_0-P_7): {correlations['p_diff_variance']:.4f}")