Delayed-Choice Entanglement Swapping on IBM's 127-Qubit Quantum Computer

The 3D Bar Plot of State Probabilities above (code below) shows the probabilities of all states measured in the experiment. A small set of states dominates the distribution, with significantly higher probabilities.
The 3D Delayed-Choice Qubit Correlations above (code below) shows the relationship between the delayed-choice qubit (q_4​) and the entangled states (q_0​ and q_3​). Clusters in the 3D space show strong correlations between the delayed-choice qubit and the entangled states. These clusters validate the role of q_4​ in retroactively influencing the outcomes of q_0​ and q_3​. The delayed-choice mechanism successfully demonstrates its influence, as evidenced by clear separations in the state probabilities conditioned on q_4.
The 3D Cumulative Probability Surface above (code below) shows the cumulative probabilities of the measured states. The surface rises steeply initially, with the first 10 - 15 states contributing the bulk of the cumulative probability. The surface flattens quickly, indicating that the majority of states contribute negligibly. The steep rise in cumulative probabilities shows the dominance of a small set of states, confirming that the system behaves as expected under ideal delayed-choice entanglement swapping conditions. The flattening means effective suppression of noise or unexpected transitions.
The 3D State Transition Path above (code below) shows transitions between dominant states and their nearest neighbors, with transition weights proportional to the probabilities of the connected states. Strong connections are observed between dominant states and their closest neighbors, showing coherence and connectivity. Sparse and weak connections to other states show limited noise interference. The transitions show that the dominant states form a coherent and strong subspace in the larger Hilbert space. This reinforces the conclusion that the delayed-choice mechanism effectively projects the system into expected quantum states with minimal leakage to unintended states.
In the end, this experiment successfully showed delayed-choice entanglement swapping on IBM's ibm_kyiv. By entangling qubits through intermediary measurements and conditioning outcomes on the delayed-choice qubit, this circuit observed strong correlations, validating the influence of later measurements on prior entanglement outcomes. The result confirms that the entanglement between q_0​ and q_3​ is established only after the measurement of q_4​, showing how future decisions can influence past quantum correlations.

Code:

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

# Logging Setup
logging.basicConfig(level=logging. INFO, format="%(asctime)s - %(levelname)s - %(message)s")
logger = logging.getLogger(__name__)

# Load Backend Calibration Data
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

# Parse Calibration Data
def select_best_qubits(calibration_data, n_qubits):
    logger. info("Selecting the best qubits based on T1, T2, and error rates")
    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 Calibration Data
calibration_file = '/Users/Downloads/ibm_kyiv_calibrations_2025-01-15T19_06_40Z.csv'
calibration_data = load_calibration_data(calibration_file)

# Select Best Qubits for Experiment
num_qubits = 6  # 4 qubits + 2 control qubits for delayed choice
best_qubits = select_best_qubits(calibration_data, num_qubits)

# IBMQ 
logger. info("Setting up IBM Q service")
service = QiskitRuntimeService(
    channel="ibm_quantum",
    instance="ibm-q/open/main",
    token='YOUR_IBMQ_API_KEY_O-`'
)

backend_name = "ibm_kyiv"
backend = service.backend(backend_name)
logger. info("Backend selected: %s", backend_name)

# Quantum and Classical Registers
qr = QuantumRegister(num_qubits, "q")
cr = ClassicalRegister(num_qubits, "c")
qc = QuantumCircuit(qr, cr)

# Create Two EPR Pairs (Qubit Pairs: (1,2) and (3,4))
qc.h(qr[0])
qc. cx(qr[0], qr[1])
qc.h(qr[2])
qc. cx(qr[2], qr[3])

# Prepare a Delayed-Choice Bell-State Measurement (Victor's Measurement)
qc.h(qr[4])  # Delayed-choice control qubit
qc.h(qr[5])  # Auxiliary qubit for measurement
qc. cx(qr[2], qr[3])  # Entangle qubits 2 and 3
qc.ccx(qr[4], qr[2], qr[5])  # Conditional entanglement based on delayed-choice qubit

# Perform Entanglement Swapping
qc.h(qr[4])  # Finalize delayed-choice decision
qc. cx(qr[4], qr[2])  # Entangle qubit 1 with Qubit 4 via intermediate qubits

# Measurement for Verification
qc.measure(qr, cr)

# Transpile 
logger. info("Transpiling the quantum circuit for the backend")
qc_transpiled = transpile(qc, backend=backend, optimization_level=3)
logger. info("Circuit transpilation complete")

# Execute 
shots = 16384
with Session(service=service, backend=backend) as session:
    sampler = SamplerV2(session=session)
    logger. info("Executing the circuit on the backend")
    job = sampler. run([qc_transpiled], shots=shots)
    job_result = job.result()

    # Extract counts
    data_bin = job_result._pub_results[0]['__value__']['data']
    if 'c' in data_bin:
        counts = data_bin['c'].get_counts()
    else:
        raise KeyError("No valid key found in data_bin to extract counts.")

    # Save Json
    results_data = {
        "experiment_name": "Delayed-Choice Entanglement Swapping",
        "raw_counts": counts
    }
    output_file = '/Users/Documents/Delayed_Choice_Entanglement_Swapping_0.json'
    with open(output_file, "w") as f:
        json.dump(results_data, f, indent=4)
    logger. info(f"Results saved to {output_file}")

# Visual
plot_histogram(counts)
plt.title("Delayed-Choice Entanglement Swapping")
plt. show()

////////////////////////
Code for All Visuals from Run Data
# Imports
import matplotlib.pyplot as plt
import json
import numpy as np
import seaborn as sns
from mpl_toolkits.mplot3d import Axes3D

# Load Run Results
file_path = '/Users/Documents/Delayed_Choice_Entanglement_Swapping_0.json'
with open(file_path, 'r') as f:
    data = json.load(f)

raw_counts = data["raw_counts"]
total_counts_v = sum(raw_counts.values())

# Convert counts to arrays
states = list(raw_counts.keys())
counts = np.array(list(raw_counts.values()))
total_counts = np.sum(counts)
probabilities = counts / total_counts_v

# Parse qubit states for analysis
q0_q3_states = [state[0] + state[5] for state in states]
q4_states = [state[4] for state in states]

# Parse qubit states for analysis
q0_states = [int(state[0]) for state in states]
q3_states = [int(state[5]) for state in states]
q4_states = [int(state[4]) for state in states]

# Define high-count states threshold
high_count_threshold = 1000
high_count_indices = np.where(counts > high_count_threshold)[0]

# State Distribution Histogram
plt.figure(figsize=(12, 6))
plt. bar(states, counts, color='blue')
plt.title("State Distribution Histogram")
plt.xlabel("States")
plt.ylabel("Counts")
plt.xticks(rotation=90, fontsize=8)
plt.tight_layout()
plt. show()

# Correlations Between Qubits q0 and q3
q0_q3_counts = {}
for state, count in zip(states, counts):
    pair = state[0] + state[5]
    q0_q3_counts[pair] = q0_q3_counts.get(pair, 0) + count

pairs = list(q0_q3_counts.keys())
pair_counts = np.array(list(q0_q3_counts.values()))
plt.figure(figsize=(6, 6))
plt. bar(pairs, pair_counts / total_counts, color='green')
plt.title("Correlations Between Qubits q0 and q3")
plt.xlabel("q0-q3 States")
plt.ylabel("Probability")
plt. show()

# Impact of Delayed-Choice Qubit q4
q4_counts = {"0": 0, "1": 0}
for state, count in zip(states, counts):
    q4 = state[4]
    q4_counts[q4] += count

plt.figure(figsize=(6, 6))
plt. bar(q4_counts.keys(), np.array(list(q4_counts.values())) / total_counts, color='purple')
plt.title("Impact of Delayed-Choice Qubit q4")
plt.xlabel("q4 State")
plt.ylabel("Probability")
plt. show()

# Noise Analysis: High vs. Low Frequency States
high_count_threshold = 100
high_states = [(state, count) for state, count in zip(states, counts) if count > high_count_threshold]
low_states = [(state, count) for state, count in zip(states, counts) if count <= high_count_threshold]

high_states_labels = [state[0] for state in high_states]
high_states_values = np.array([count for _, count in high_states])
low_states_labels = [state[0] for state in low_states]
low_states_values = np.array([count for _, count in low_states])

plt.figure(figsize=(12, 6))
plt. bar(high_states_labels, high_states_values, color='red', label='High Frequency States')
plt. bar(low_states_labels, low_states_values, color='blue', label='Low Frequency States', alpha=0.7)
plt.title("Noise Analysis: High vs. Low Frequency States")
plt.xlabel("States")
plt.ylabel("Counts")
plt.legend()
plt.xticks(rotation=90, fontsize=8)
plt.tight_layout()
plt. show()

# Cumulative Probability Distribution
sorted_probs = np.sort(probabilities)[::-1]
cumulative_probs = np.cumsum(sorted_probs)
plt.figure(figsize=(10, 6))
plt.plot(range(len(cumulative_probs)), cumulative_probs, marker='o', linestyle='-')
plt.title("Cumulative Probability Distribution")
plt.xlabel("State Index (Sorted by Probability)")
plt.ylabel("Cumulative Probability")
plt.grid()
plt. show()

# Heatmap of State Probabilities
plt.figure(figsize=(12, 8))
heatmap_data = np.array(probabilities).reshape(-1, 1)  # Convert to 2D for heatmap
sns.heatmap(heatmap_data, annot=False, cmap="viridis", cbar=True, xticklabels=["Probability"], yticklabels=states)
plt.title("Heatmap of State Probabilities")
plt.ylabel("States")
plt.xlabel("Probability")
plt. show()

# Histogram of High-Count States
high_count_threshold = 1000
high_count_states = [(state, count) for state, count in zip(states, counts) if count > high_count_threshold]
high_labels = [state for state, count in high_count_states]
high_values = [count for state, count in high_count_states]

plt.figure(figsize=(12, 6))
plt. bar(high_labels, high_values, color='orange')
plt.title("Histogram of High-Count States")
plt.xlabel("States")
plt.ylabel("Counts")
plt.xticks(rotation=90, fontsize=8)
plt.tight_layout()
plt. show()

# Probability of Delayed-Choice Correlations
q4_correlations = {"0": {"00": 0, "01": 0, "10": 0, "11": 0}, "1": {"00": 0, "01": 0, "10": 0, "11": 0}}
for state, count in zip(states, counts):
    q4 = state[4]
    q0_q3 = state[0] + state[5]
    q4_correlations[q4][q0_q3] += count

# Normalize correlations
for q4 in q4_correlations:
    total = sum(q4_correlations[q4].values())
    for q0_q3 in q4_correlations[q4]:
        q4_correlations[q4][q0_q3] /= total

fig, axes = plt.subplots(1, 2, figsize=(12, 6), sharey=True)
for idx, (q4, correlations) in enumerate(q4_correlations.items()):
    axes[idx].bar(correlations.keys(), correlations.values(), color=['blue', 'green', 'red', 'purple'])
    axes[idx].set_title(f"Correlations for q4 = {q4}")
    axes[idx].set_xlabel("q0-q3 States")
    axes[idx].set_ylabel("Probability")

plt.tight_layout()
plt. show()

# 3D Bar Plot of State Probabilities
fig = plt.figure(figsize=(10, 7))
ax = fig.add_subplot(111, projection='3d')
x = np.arange(len(states))
y = np.zeros_like(x)
z = np.zeros_like(x)
dx = np.ones_like(x)
dy = np.ones_like(x)
dz = probabilities
ax.bar3d(x, y, z, dx, dy, dz, shade=True)
ax.set_title("3D Bar Plot of State Probabilities")
ax.set_xlabel("State Index")
ax.set_ylabel("Y Index")
ax.set_zlabel("Probability")
plt. show()

# 3D Cumulative Probability Surface
fig = plt.figure(figsize=(12, 8))
ax = fig.add_subplot(111, projection='3d')

# Prepare data for the surface
sorted_probs = np.sort(probabilities)[::-1]
cumulative_probs = np.cumsum(sorted_probs)
x = np.arange(len(sorted_probs))
y = np.arange(2)  
X, Y = np.meshgrid(x, y)  # Create a 2D grid
Z = np.tile(cumulative_probs, (2, 1))  # Make Z 2D by repeating the cumulative_probs

ax.plot_surface(X, Y, Z, cmap='viridis', edgecolor='none')
ax.set_title("3D Cumulative Probability Surface")
ax.set_xlabel("State Index (Sorted by Probability)")
ax.set_ylabel("Y Index")
ax.set_zlabel("Cumulative Probability")
plt. show()

# 3D Delayed-Choice Qubit Correlations
fig = plt.figure(figsize=(10, 7))
ax = fig.add_subplot(111, projection='3d')
ax.scatter(q0_states, q3_states, q4_states, c=probabilities, cmap='plasma', s=100)
ax.set_title("3D Delayed-Choice Qubit Correlations")
ax.set_xlabel("q0 State")
ax.set_ylabel("q3 State")
ax.set_zlabel("q4 State")
plt. show()

# 3D State Transition Path
fig = plt.figure(figsize=(12, 8))
ax = fig.add_subplot(111, projection='3d')

# Extract high-count states and their neighbors
dominant_indices = np.where(probabilities > 0.05)[0]  # Threshold for significant probabilities
x_transitions = []
y_transitions = []
z_transitions = []
weights = []

for idx in dominant_indices:
    state = states[idx]
    prob = probabilities[idx]
    # Generate neighboring transition states by flipping bits
    for flip in range(len(state)):
        neighbor = list(state)
        neighbor[flip] = '1' if neighbor[flip] == '0' else '0'  # Flip the bit
        neighbor_str = ''.join(neighbor)
        if neighbor_str in states:
            neighbor_idx = states.index(neighbor_str)
            x_transitions.append(idx)
            y_transitions.append(neighbor_idx)
            z_transitions.append(0)  
            weights.append(prob * probabilities[neighbor_idx])  # Weighted influence

# Normalize weights for better visualization
weights = np.array(weights)
weights = weights / max(weights)

# Plot state transitions
ax.scatter(x_transitions, y_transitions, z_transitions, c=weights, cmap='coolwarm', s=100 * weights)
for i, (x, y, z, w) in enumerate(zip(x_transitions, y_transitions, z_transitions, weights)):
    ax.plot([x, x], [y, y], [z, w], 'k-', alpha=0.3)  

ax.set_title("3D State Transition Path")
ax.set_xlabel("Dominant State Index")
ax.set_ylabel("Neighbor State Index")
ax.set_zlabel("Transition Weight")
plt. show()

# End