Twistor-Encoded Quantum State Stabilization with IBM's 127-Qubit Quantum Computer

The Mutual Information Between Qubits (MIBQ) above (full code on Qwork) shows Classical correlation between the measurement outcomes c_0​ and c_1​, plotted over increasing idle delay. At 0 - 900 ns, the mutual information stays nearly zero. This confirms that the entangled state’s correlation is non-classical and the outcomes of the two qubits remain statistically independent, exactly what we expect from teleportation and coherent Bell states. At ~1800 ns, MI begins rising, peaking at ~0.025 bits at 3640 ns. The twistor encoding sustains non-classical independence, a signature of quantum teleportation and coherent entanglement. The eventual rise in MI indicates coherence decay into classical correlation, the state starts to 'bleed' quantum information into observable statistical structure. Twistor encoding delays the quantum-to-classical transition. Mutual information is a quantitative marker of decoherence onset, and this curve lets us pinpoint the exact breakdown threshold.
The Entropy Growth Rate (EGR) (full code on Qwork) shows the first derivative of the Shannon entropy over time: (dS)/(dt)​, showing how quickly the system is losing information and becoming mixed. Nearly zero growth until ~450 ns. Then, a rapid increase, peaking near 0.0003 bits/ns around 900 ns. After that, the system enters a plateaued high-entropy-growth regime. This curve identifies a quantum phase transition, from a structured, geometry-guided manifold into a chaotic, maximally mixed space. The abrupt jump mirrors the fidelity derivative’s behavior, both identifying ~900 ns as the Twistor-protection decay point. Twistor encoding slows down entropy production. Entropy rate spiking corresponds to the state escaping from the geometric basin of attraction, entering standard noisy Hilbert diffusion.
The Heatmap: Bitstring Drift From Ideal (HBDFI) above (full code on Qwork) shows population drift across Hamming shells from the ideal state '00', plotted as a heatmap. Distance 0: starts at 85%, gently decays, then collapses hard to 49% at 3640 ns. Distance 1 ('01', '10') rises steadily. Distance 2 ('11') is nearly zero for early delays, then surges at the end. Twistor encoding traps error inside a low-dimensional geodesic shell early on, limiting error to Hamming Distance 1. Only after prolonged delay do we see spread into higher Hamming distance zones, indicating true multi-error decoherence. This validates the idea that Twistor geometry compresses the error space into structured shells. The heatmap is a direct error-topology fingerprint. It confirms the presence of error basin formation, where only certain transitions are geometrically allowed until that structure fails.
The Hilbert Phase Drift of State Projection (HPDoSP) above (full code on Qwork) shows The quantum state, projected as a complex number using bitstring mappings, with the phase (arg) of the projection plotted over time. Smooth linear phase drift beginning near zero and rising to over 1 radian. No oscillations, bifurcations, or chaotic jumps, indicating deterministic evolution. Twistor states evolve along a geodesic phase path, and this path is not disrupted by noise until very late. The linearity of this drift suggests a holomorphic flow, consistent with twistor theory, where complex projective evolution governs time. This confirms the twistor-phase-locking hypothesis: noise does not induce random phase collapse, but a structured rotation through the Hilbert sphere until geometric constraints decay.
In the end, this experiment implemented a novel twistor-inspired quantum error correction (QEC) protocol on IBM's 127-qubit 'ibm_sherbrooke' to test how long quantum coherence can be passively preserved in a geometrically encoded Bell state. By applying a structured Twistor entanglement and introducing idle delays of increasing duration, this circuit measured the degradation of fidelity, entropy growth, mutual information, and phase evolution. The results show that the Twistor encoding significantly suppressed decoherence early on, maintaining high fidelity (0.85 at 0 ns) and near-zero mutual information, while constraining errors within low Hamming shells and minimizing entropy dispersion. Only after ~1800 ns did fidelity rapidly decay, entropy and classical correlations rise, and the state exhibit topological collapse. These findings confirm that Twistor geometry acts as a passive QEC mechanism, redirecting noise along geodesic attractors in Hilbert space, delaying entropic drift, and preserving quantum information via holomorphic evolution, offering a hint into the possibility of geometry-native fault tolerance in quantum computation.

Code:

# Main Circuit
# Imports
import numpy as np
import json
import logging
import pandas as pd
from qiskit import QuantumCircuit, transpile, QuantumRegister, ClassicalRegister
from qiskit_ibm_runtime import QiskitRuntimeService, Session, SamplerV2
from qiskit.circuit.library import RYGate, CXGate
from qiskit.visualization import plot_histogram
import matplotlib.pyplot as plt

# Setup logging
logging.basicConfig(level=logging. INFO)
logger = logging.getLogger(__name__)

# Load IBMQ account
service = QiskitRuntimeService(
    channel='ibm_quantum',
    instance='ibm-q/open/main',
    token='YOUR_IBMQ_API_O-`'  
)

backend_name = 'ibm_sherbrooke'
backend = service.backend(backend_name)

# Backend timing resolution (dt in seconds)
dt = backend.configuration().dt
logger. info(f"Backend dt: {dt:.10e} seconds")

# Load calibration data and select best qubits
def load_calibration_data(file_path):
    logger. info("Loading calibration data from %s", file_path)
    calibration_data = pd. read_csv(file_path)
    calibration_data.columns = calibration_data.columns.str.strip()
    logger. info("Calibration data loaded successfully")
    return calibration_data

def select_best_qubits(calibration_data, n_qubits):
    logger. info("Selecting best qubits based on T1, T2, and √x error")
    qubits_sorted = calibration_data.sort_values(
        by=["√x (sx) error", "T1 (us)", "T2 (us)"],
        ascending=[True, False, False]
    )
    best_qubits = qubits_sorted["Qubit"].head(n_qubits).tolist()
    logger. info("Selected qubits: %s", best_qubits)
    return best_qubits

# Load most recent calibration
data_path = '/Users/steventippeconnic/Downloads/ibm_sherbrooke_calibrations_2025-04-03T05_07_35Z.csv'
calibration_data = load_calibration_data(data_path)
selected_qubits = select_best_qubits(calibration_data, 2)

# Only use exact dt-multiples 
valid_delays_dt_units = [0, 1024, 2048, 4096, 8192, 16384]
results = {}

# Twistor encoding function
def twistor_encoding(qc, qubit_indices):
    for i, q in enumerate(qubit_indices):
        theta = np.pi / 4 * (i + 1)
        qc.append(RYGate(theta), [q])
        if i < len(qubit_indices) - 1:
            qc.append(CXGate(), [q, qubit_indices[i + 1]])
    qc.barrier()

# Run experiments for each delay duration
for delay_dt in valid_delays_dt_units:
    delay_ns = int(delay_dt * dt * 1e9)
    logger. info(f"Running circuit with delay = {delay_dt} dt units (~{delay_ns} ns)")

    qr = QuantumRegister(2, 'q')  # q0 and q1 in Bell state
    cr = ClassicalRegister(2, 'c')  # Measure both
    qc = QuantumCircuit(qr, cr)

    # Prepare Bell state
    qc.h(qr[0])
    qc .cx(qr[0], qr[1])
    qc.barrier()

    # Twistor entanglement
    twistor_encoding(qc, [qr[0], qr[1]])

    # Idle delay in dt units
    qc.delay(delay_dt, qr[0], unit='dt')
    qc.delay(delay_dt, qr[1], unit='dt')
    qc.barrier()

    # Bell-basis measurement
    qc. cx(qr[0], qr[1])
    qc.h(qr[0])
    qc.barrier()
    qc.measure(qr[0], cr[0])
    qc.measure(qr[1], cr[1])

    # Transpile and run
    transpiled_qc = transpile(qc, backend=backend, optimization_level=3, scheduling_method='alap')

    with Session(service=service, backend=backend) as session:
        sampler = SamplerV2(session=session)
        job = sampler. run([transpiled_qc], shots=8192)
        job_result = job.result()

        data_bin = job_result._pub_results[0]['__value__']['data']
        classical_register = transpiled_qc.cregs[0].name
        counts = data_bin[classical_register].get_counts() if classical_register in data_bin else {}

        # Save result 
        results[str(delay_dt)] = counts

# Save results to Json
output_path = '/Users/steventippeconnic/Documents/QC/Twistor_QEC_Stability_Results_0.json'
with open(output_path, 'w') as f:
    json.dump(results, f, indent=4)

# Visual
for delay, counts in results.items():
    plot_histogram(counts, title=f"Delay: {delay} dt")
    plt. show()

# Fidelity Calculation for delay = 0
counts = results['0']
total_shots = sum(counts.values())
correct = counts.get('00', 0) + counts.get('11', 0)  # Bell state measurement ideal
fidelity = correct / total_shots
print(f"Fidelity at 0 delay: {fidelity:.5f}")

# End. 

# /////////////////////////////////
# Full code for all visualizations using run data 
# Imports 
import json
import matplotlib.pyplot as plt
import numpy as np
from scipy.special import rel_entr
from collections import Counter
from scipy.stats import entropy
from math import log2

# Load results
file_path = '/Users/steventippeconnic/Documents/QC/Twistor_QEC_Stability_Results_0.json'
with open(file_path, 'r') as f:
    data = json.load(f)

# Define IBM Sherbrooke dt (from run and in seconds)
dt = 2.2222222222222221e-10

# Normalize counts
def normalize(counts):
    total = sum(counts.values())
    return {k: v / total for k, v in counts.items()}

# Parse delays and convert to ns for plotting
delays_dt = sorted([int(k) for k in data.keys()])
delays_ns = [dt * d * 1e9 for d in delays_dt]  

# Constants
ideal_bits = {'00', '11'}

# Fidelity vs Delay
fidelities = []
for d in delays_dt:
    counts = data[str(d)]
    total = sum(counts.values())
    correct = sum(counts.get(b, 0) for b in ideal_bits)
    fidelities.append(correct / total)

plt.figure(figsize=(8, 6))
plt.plot(delays_ns, fidelities, 'o-', label='Fidelity')
plt.xlabel('Delay (ns)')
plt.ylabel('Fidelity')
plt.title('Fidelity vs Delay')
plt.grid(True)
plt.show()

# Shannon Entropy vs Delay
def shannon_entropy(counts):
    total = sum(counts.values())
    probs = [v / total for v in counts.values()]
    return -sum(p * np.log2(p) for p in probs if p > 0)

entropies = [shannon_entropy(data[str(d)]) for d in delays_dt]

plt.figure(figsize=(8, 6))
plt.plot(delays_ns, entropies, 'o--', color='purple')
plt.xlabel('Delay (ns)')
plt.ylabel('Shannon Entropy (bits)')
plt.title('Entropy of Measurement Outcomes')
plt.grid(True)
plt.show()

# Hamming Shell Distribution
def hamming_shells(counts):
    shells = {0: 0, 1: 0, 2: 0}
    for bitstring, count in counts.items():
        hd = sum([bitstring[i] != '0' for i in range(2)])
        shells[hd] += count
    total = sum(counts.values())
    return [shells[0]/total, shells[1]/total, shells[2]/total]

shell_distributions = np.array([hamming_shells(data[str(d)]) for d in delays_dt])

plt.figure(figsize=(8, 6))
plt.plot(delays_ns, shell_distributions[:,0], label='Distance 0')
plt.plot(delays_ns, shell_distributions[:,1], label='Distance 1')
plt.plot(delays_ns, shell_distributions[:,2], label='Distance 2')
plt.xlabel('Delay (ns)')
plt.ylabel('Proportion')
plt.title('Error Shell Distribution (Hamming Distance)')
plt.legend()
plt.grid(True)
plt.show()

# Bitstring Evolution (for 00, 01, 10, 11)
labels = ['00', '01', '10', '11']
bit_counts = {label: [data[str(d)].get(label, 0)/sum(data[str(d)].values()) for d in delays_dt] for label in labels}

plt.figure(figsize=(8, 6))
for label in labels:
    plt.plot(delays_ns, bit_counts[label], label=label)
plt.xlabel('Delay (ns)')
plt.ylabel('Normalized Frequency')
plt.title('Bitstring Probability Evolution')
plt.legend()
plt.grid(True)
plt.show()

# Parity Bias over Time
even_parity = []
odd_parity = []
for d in delays_dt:
    counts = normalize(data[str(d)])
    even = sum(v for k, v in counts.items() if k.count('1') % 2 == 0)
    odd = sum(v for k, v in counts.items() if k.count('1') % 2 == 1)
    even_parity.append(even)
    odd_parity.append(odd)

plt.figure(figsize=(8, 6))
plt.plot(delays_ns, even_parity, label='Even Parity')
plt.plot(delays_ns, odd_parity, label='Odd Parity')
plt.xlabel('Delay (ns)')
plt.ylabel('Probability')
plt.title('Parity Symmetry Evolution')
plt.legend()
plt.grid(True)
plt.show()

# Fidelity Derivative (Rate of Change)
ideal_bits = {'00', '11'}
fidelities = []
for d in delays_dt:
    counts = data[str(d)]
    total = sum(counts.values())
    correct = sum(counts.get(b, 0) for b in ideal_bits)
    fidelities.append(correct / total)

fidelity_deriv = np.gradient(fidelities, delays_ns)

plt.figure(figsize=(8, 6))
plt.plot(delays_ns, fidelity_deriv, 'o-', color='red')
plt.xlabel('Delay (ns)')
plt.ylabel('dFidelity/dt')
plt.title('Fidelity Decay Rate')
plt.grid(True)
plt.axhline(0, linestyle='--', color='gray')
plt.show()

# KL Divergence from Initial Distribution
def kl_div(p, q):
    # Fill missing keys
    all_keys = set(p.keys()).union(q.keys())
    p_arr = np.array([p.get(k, 1e-10) for k in all_keys])
    q_arr = np.array([q.get(k, 1e-10) for k in all_keys])
    return np.sum(rel_entr(p_arr, q_arr))

baseline = normalize(data[str(delays_dt[0])])
kl_values = []
for d in delays_dt:
    dist = normalize(data[str(d)])
    kl = kl_div(dist, baseline)
    kl_values.append(kl)

plt.figure(figsize=(8, 6))
plt.plot(delays_ns, kl_values, 'o-', color='darkgreen')
plt.xlabel('Delay (ns)')
plt.ylabel('KL Divergence (bits)')
plt.title('KL Divergence from Initial Distribution')
plt.grid(True)
plt.show()

# Bitstring Projection into Complex Plane (ℂ^2)
def bitstring_to_complex(s):
    # Map bitstring to complex number: z = c0 + i·c1
    return int(s[0]) + 1j * int(s[1])

complex_points = []
weights = []

for d in delays_dt:
    counts = normalize(data[str(d)])
    total = sum(counts.values())
    z = sum(bitstring_to_complex(k) * v for k, v in counts.items())
    complex_points.append(z)
    weights.append(abs(z))

xs = [z.real for z in complex_points]
ys = [z.imag for z in complex_points]

plt.figure(figsize=(8, 6))
plt.plot(xs, ys, 'o-', linewidth=2)
for i, d in enumerate(delays_ns):
    plt.text(xs[i] + 0.02, ys[i] + 0.02, f'{int(d)}ns', fontsize=9)
plt.xlabel('Re')
plt.ylabel('Im')
plt.title('Quantum State Collapse Path in Complex Plane')
plt.grid(True)
plt.axhline(0, color='gray', linestyle='--')
plt.axvline(0, color='gray', linestyle='--')
plt.show()

# Mutual Information between c0 and c1
def mutual_information(counts):
    total = sum(counts.values())
    px = {'0': 0, '1': 0}
    py = {'0': 0, '1': 0}
    pxy = {'00': 0, '01': 0, '10': 0, '11': 0}

    for bitstring, count in counts.items():
        x = bitstring[0]
        y = bitstring[1]
        pxy[bitstring] += count / total
        px[x] += count / total
        py[y] += count / total

    mi = 0
    for xy, p in pxy.items():
        if p > 0:
            mi += p * log2(p / (px[xy[0]] * py[xy[1]]))
    return mi

mi_values = [mutual_information(data[str(d)]) for d in delays_dt]

plt.figure(figsize=(8, 6))
plt.plot(delays_ns, mi_values, marker='o')
plt.title("Mutual Information Between Qubits")
plt.xlabel("Delay (ns)")
plt.ylabel("MI(c₀ : c₁) [bits]")
plt.grid(True)
plt.show()

# Entropy Velocity (dS/dt)
def shannon_entropy(counts):
    total = sum(counts.values())
    probs = [v / total for v in counts.values()]
    return -sum(p * log2(p) for p in probs if p > 0)

entropies = [shannon_entropy(data[str(d)]) for d in delays_dt]
entropy_rate = np.gradient(entropies, delays_ns)

plt.figure(figsize=(8, 6))
plt.plot(delays_ns, entropy_rate, marker='o', color='orange')
plt.title("Entropy Growth Rate (dS/dt)")
plt.xlabel("Delay (ns)")
plt.ylabel("Entropy Rate [bits/ns]")
plt.grid(True)
plt.axhline(0, linestyle='--', color='gray')
plt.show()

# Bitstring Transition Heatmap (vs 00)
import seaborn as sns
import pandas as pd

bitstrings = ['00', '01', '10', '11']
distance_from_00 = {'00': 0, '01': 1, '10': 1, '11': 2}

df = pd.DataFrame(index=delays_ns, columns=['Distance 0', 'Distance 1', 'Distance 2'])
for i, d in enumerate(delays_dt):
    counts = normalize(data[str(d)])
    shell_counts = {0: 0, 1: 0, 2: 0}
    for b in bitstrings:
        shell_counts[distance_from_00[b]] += counts.get(b, 0)
    for k in shell_counts:
        df.loc[delays_ns[i], f"Distance {k}"] = shell_counts[k]

df = df.astype(float)
sns.heatmap(df.transpose(), annot=True, cmap='rocket', cbar_kws={'label': 'Probability'})
plt.title("Heatmap: Bitstring Drift from Ideal (00)")
plt.xlabel("Delay (ns)")
plt.ylabel("Hamming Shell")
plt.show()

# Hilbert Phase Drift Map
def bitstring_to_complex(bs):
    return int(bs[0]) + 1j * int(bs[1])

complex_points = []
for d in delays_dt:
    counts = normalize(data[str(d)])
    z = sum(bitstring_to_complex(b) * p for b, p in counts.items())
    complex_points.append(z)

phases = [np.angle(z) for z in complex_points]

plt.figure(figsize=(8, 6))
plt.plot(delays_ns, phases, marker='o', color='green')
plt.title("Hilbert Phase Drift of State Projection")
plt.xlabel("Delay (ns)")
plt.ylabel("Phase [radians]")
plt.grid(True)
plt.axhline(0, linestyle='--', color='gray')
plt.show()