Use this Skill for quantum computing experiments with Qiskit: qubit circuits, gates, measurement, VQE/QAOA variational algorithms, noise models, and statevector vs shot simulation.
Install with Codex or Claude Copy this prompt, paste it into Codex, Claude, or another assistant, and let it review the skill page and install it for you.
A direct command skips the review prompt. Inspect the source before running it.
Use this Skill for quantum computing experiments with Qiskit: qubit circuits, gates, measurement, VQE/QAOA variational algorithms, noise models, and statevector vs shot simulation.
TL;DR — Build and simulate quantum circuits with Qiskit. Use statevector
or shot-based simulation, add realistic noise models, implement VQE for ground
state energy estimation, and QAOA for combinatorial optimization (MaxCut).
When to Use
Use this Skill when you need to:
Design and test quantum gate circuits (H, CNOT, Ry, Rz, T, S, SWAP)
Simulate quantum circuits on classical hardware (Aer statevector / QASM)
"""
Build a Bell state circuit, simulate with statevector and shot-based simulators,
and plot the measurement histogram.
"""import numpy as np
import matplotlib.pyplot as plt
from qiskit import QuantumCircuit, transpile
from qiskit_aer import AerSimulator
from qiskit.visualization import plot_histogram
defcreate_bell_circuit() -> QuantumCircuit:
"""
Create a 2-qubit Bell state circuit |Φ+⟩ = (|00⟩ + |11⟩)/√2.
Circuit:
q0: ──H──●──
q1: ─────X──
"""
qc = QuantumCircuit(2, 2)
qc.h(0) # Hadamard on qubit 0
qc.cx(0, 1) # CNOT: control=0, target=1
qc.measure([0, 1], [0, 1])
return qc
defsimulate_statevector(qc_no_measure: QuantumCircuit) -> np.ndarray:
"""
Get exact statevector amplitudes before measurement.
Args:
qc_no_measure: Circuit without measurement instructions.
Returns:
Complex statevector array of length 2^n.
"""
sim = AerSimulator(method="statevector")
qc_sv = qc_no_measure.copy()
qc_sv.save_statevector()
job = sim.run(transpile(qc_sv, sim))
result = job.result()
sv = np.array(result.get_statevector())
return sv
defsimulate_shots(qc: QuantumCircuit, n_shots: int = 4096) -> dict:
"""
Run shot-based simulation and return measurement counts.
Args:
qc: Circuit with measurement gates.
n_shots: Number of repetitions.
Returns:
Dictionary of bitstring: count.
"""
sim = AerSimulator(method="qasm")
job = sim.run(transpile(qc, sim), shots=n_shots)
return job.result().get_counts()
defplot_counts(counts: dict, title: str = "Bell State", output: str = "bell.png") -> None:
"""Plot measurement histogram."""
fig, ax = plt.subplots(figsize=(5, 4))
states = sorted(counts.keys())
vals = [counts[s] / sum(counts.values()) for s in states]
ax.bar(states, vals, color="#4C72B0")
ax.set_xlabel("Measurement outcome")
ax.set_ylabel("Probability")
ax.set_title(title)
fig.tight_layout()
fig.savefig(output, dpi=150)
print(f"Histogram saved to {output}")
if __name__ == "__main__":
# Statevector
qc_bare = QuantumCircuit(2)
qc_bare.h(0)
qc_bare.cx(0, 1)
sv = simulate_statevector(qc_bare)
print("Statevector amplitudes:", sv)
print(f"|00⟩ probability: {abs(sv[0])**2:.4f}")
print(f"|11⟩ probability: {abs(sv[3])**2:.4f}")
# Shot-based
qc = create_bell_circuit()
counts = simulate_shots(qc, n_shots=8192)
print("Shot counts:", counts)
plot_counts(counts, title="Bell State |Φ+⟩", output="bell_histogram.png")
Step 2 — VQE for H₂ Ground State Energy
"""
VQE: Variational Quantum Eigensolver for estimating the ground state energy
of a simple 2-qubit Hamiltonian representing H₂ at fixed bond length.
H = -1.0523 * II + 0.3979 * ZI - 0.3979 * IZ - 0.0112 * ZZ + 0.1809 * XX
(Jordan-Wigner encoding of H₂ minimal basis Hamiltonian)
"""import numpy as np
from qiskit.circuit.library import EfficientSU2
from qiskit.quantum_info import SparsePauliOp
from qiskit_aer.primitives import Estimator
from qiskit_algorithms import VQE
from qiskit_algorithms.optimizers import COBYLA, SPSA
import matplotlib.pyplot as plt
defbuild_h2_hamiltonian() -> SparsePauliOp:
"""
Construct the H₂ Hamiltonian in the minimal (STO-3G) basis at R=0.735 Å.
Returns:
SparsePauliOp representing the qubit Hamiltonian (energy in Hartree).
"""# Pauli coefficients from Jordan-Wigner transformation# Strings are in reverse qubit order (Qiskit convention)
H = SparsePauliOp.from_list([
("II", -1.0523732),
("IZ", 0.3979374),
("ZI", -0.3979374),
("ZZ", -0.0112801),
("XX", 0.1809312),
])
return H
defrun_vqe(
hamiltonian: SparsePauliOp,
reps: int = 2,
max_iter: int = 200,
seed: int = 42,
) -> dict:
"""
Run VQE with EfficientSU2 ansatz to find ground state energy.
Args:
hamiltonian: Target Hamiltonian as SparsePauliOp.
reps: Number of repetition layers in EfficientSU2.
max_iter: Maximum optimizer iterations.
seed: Random seed for reproducibility.
Returns:
Dictionary with keys: energy, optimal_params, num_evals.
"""
n_qubits = hamiltonian.num_qubits
# Ansatz: hardware-efficient EfficientSU2
ansatz = EfficientSU2(n_qubits, reps=reps, entanglement="linear")
print(f"Ansatz parameters: {ansatz.num_parameters}")
# Estimator primitive (Aer)
estimator = Estimator(approximation=True)
# COBYLA optimizer (gradient-free)
optimizer = COBYLA(maxiter=max_iter)
np.random.seed(seed)
initial_params = np.random.uniform(-np.pi, np.pi, ansatz.num_parameters)
# Track convergence
energies = []
defcallback(nfev, params, energy, stddev):
energies.append(energy)
iflen(energies) % 50 == 0:
print(f" Iter {len(energies):4d}: E = {energy:.6f} Ha")
vqe = VQE(estimator=estimator, ansatz=ansatz, optimizer=optimizer,
callback=callback, initial_point=initial_params)
result = vqe.compute_minimum_eigenvalue(hamiltonian)
print(f"\nVQE ground state energy: {result.eigenvalue.real:.6f} Hartree")
print(f" FCI (exact) reference: -1.137270 Hartree")
print(f" Error: {abs(result.eigenvalue.real - (-1.137270)):.2e} Ha")
# Plot convergence
fig, ax = plt.subplots(figsize=(7, 4))
ax.plot(energies, color="#2196F3", linewidth=1.2)
ax.axhline(-1.137270, color="red", linestyle="--", linewidth=0.8, label="FCI energy")
ax.set_xlabel("Optimizer iteration")
ax.set_ylabel("Energy (Hartree)")
ax.set_title("VQE Convergence for H₂")
ax.legend()
fig.tight_layout()
fig.savefig("vqe_convergence.png", dpi=150)
return {
"energy": result.eigenvalue.real,
"optimal_params": result.optimal_parameters,
"num_evals": result.cost_function_evals,
"convergence": energies,
}
if __name__ == "__main__":
H = build_h2_hamiltonian()
result = run_vqe(H, reps=2, max_iter=300)
print(f"Converged in {result['num_evals']} function evaluations")
Step 3 — QAOA for MaxCut on a 4-Node Graph
"""
QAOA: Quantum Approximate Optimization Algorithm for the MaxCut problem.
MaxCut: partition graph nodes into two sets S and V\S to maximize the number
of edges between them.
Cost Hamiltonian: H_C = Σ_{(i,j)∈E} (I - Z_i Z_j) / 2
"""import numpy as np
import matplotlib.pyplot as plt
import networkx as nx
from qiskit import QuantumCircuit
from qiskit.quantum_info import SparsePauliOp
from qiskit_aer.primitives import Sampler
from qiskit_algorithms import QAOA
from qiskit_algorithms.optimizers import COBYLA
defbuild_maxcut_hamiltonian(graph: nx.Graph) -> SparsePauliOp:
"""
Build the MaxCut cost Hamiltonian from a NetworkX graph.
H_C = Σ_{(i,j)∈E} (I - Z_i Z_j) / 2
Args:
graph: NetworkX graph with integer node labels 0..n-1.
Returns:
SparsePauliOp for the MaxCut Hamiltonian.
"""
n = graph.number_of_nodes()
pauli_list = []
for u, v in graph.edges():
# Z_u Z_v term: Pauli string of length n with Z at positions u and v
zz_str = ["I"] * n
zz_str[u] = "Z"
zz_str[v] = "Z"# Qiskit uses reversed qubit ordering
pauli_list.append(("".join(reversed(zz_str)), -0.5))
pauli_list.append(("I" * n, 0.5))
return SparsePauliOp.from_list(pauli_list)
defrun_qaoa_maxcut(
graph: nx.Graph,
p: int = 2,
max_iter: int = 200,
) -> dict:
"""
Run QAOA to solve MaxCut on the given graph.
Args:
graph: NetworkX graph.
p: QAOA depth (number of cost+mixer layers).
max_iter: Maximum optimizer iterations.
Returns:
Dictionary with keys: best_cut, best_bitstring, counts, energy.
"""
hamiltonian = build_maxcut_hamiltonian(graph)
n = graph.number_of_nodes()
sampler = Sampler()
optimizer = COBYLA(maxiter=max_iter)
qaoa = QAOA(sampler=sampler, optimizer=optimizer, reps=p)
result = qaoa.compute_minimum_eigenvalue(hamiltonian)
# Get best measurement outcome
counts = result.best_measurement
best_bitstring = result.best_measurement["bitstring"]
partition = [int(b) for b in best_bitstring]
# Count cut edges
cut_value = sum(
1for u, v in graph.edges() if partition[u] != partition[v]
)
max_possible = graph.number_of_edges()
print(f"QAOA (p={p}) MaxCut result:")
print(f" Best bitstring: {best_bitstring}")
print(f" Cut edges: {cut_value} / {max_possible}")
print(f" Approximation ratio: {cut_value / max_possible:.3f}")
return {
"best_bitstring": best_bitstring,
"cut_value": cut_value,
"energy": result.eigenvalue.real,
}
defplot_maxcut_result(graph: nx.Graph, bitstring: str, output: str = "maxcut.png") -> None:
"""Visualize the MaxCut partition on the graph."""
partition = [int(b) for b in bitstring]
colors = ["#E91E63"if p == 0else"#2196F3"for p in partition]
fig, ax = plt.subplots(figsize=(5, 4))
pos = nx.spring_layout(graph, seed=42)
nx.draw(graph, pos, ax=ax, node_color=colors, with_labels=True,
node_size=600, font_color="white", font_weight="bold")
cut_edges = [(u, v) for u, v in graph.edges() if partition[u] != partition[v]]
nx.draw_networkx_edges(graph, pos, edgelist=cut_edges, ax=ax,
edge_color="green", width=2.5, style="dashed")
ax.set_title(f"MaxCut: {len(cut_edges)} edges cut\n(cut edges in dashed green)")
fig.tight_layout()
fig.savefig(output, dpi=150)
print(f"MaxCut visualization saved to {output}")
if __name__ == "__main__":
# 4-node cycle graph (max cut = 4)
G = nx.cycle_graph(4)
G.add_edge(0, 2) # 5-edge graph with max cut = 4
result = run_qaoa_maxcut(G, p=2)
plot_maxcut_result(G, result["best_bitstring"], "maxcut_result.png")
Advanced Usage
Noise Model Simulation
"""
Add a realistic noise model to a circuit and compare noisy vs ideal results.
"""import numpy as np
from qiskit import QuantumCircuit, transpile
from qiskit_aer import AerSimulator
from qiskit_aer.noise import NoiseModel, depolarizing_error, thermal_relaxation_error
defbuild_noise_model(
p_depol_1q: float = 0.001,
p_depol_2q: float = 0.01,
T1: float = 50e3, # ns
T2: float = 70e3, # ns
gate_time_1q: float = 50, # ns
gate_time_2q: float = 300, # ns) -> NoiseModel:
"""
Build a noise model with depolarizing + thermal relaxation errors.
Args:
p_depol_1q: Depolarizing probability for single-qubit gates.
p_depol_2q: Depolarizing probability for two-qubit gates.
T1: Longitudinal relaxation time (ns).
T2: Transverse relaxation time (ns).
gate_time_1q: Duration of single-qubit gates (ns).
gate_time_2q: Duration of two-qubit gates (ns).
Returns:
Qiskit NoiseModel.
"""
noise_model = NoiseModel()
# Thermal relaxation (T1, T2 decoherence)
error_1q = thermal_relaxation_error(T1, T2, gate_time_1q)
error_2q = thermal_relaxation_error(T1, T2, gate_time_2q).expand(
thermal_relaxation_error(T1, T2, gate_time_2q))
# Depolarizing noise
dep_1q = depolarizing_error(p_depol_1q, 1)
dep_2q = depolarizing_error(p_depol_2q, 2)
# Combine: apply both thermal + depolarizing
noise_model.add_all_qubit_quantum_error(error_1q.compose(dep_1q), ["h", "ry", "rz"])
noise_model.add_all_qubit_quantum_error(error_2q.compose(dep_2q), ["cx"])
return noise_model
defcompare_ideal_vs_noisy(n_qubits: int = 4, n_shots: int = 8192) -> None:
"""Build a GHZ state circuit and compare ideal and noisy measurements."""
qc = QuantumCircuit(n_qubits, n_qubits)
qc.h(0)
for i inrange(n_qubits - 1):
qc.cx(i, i + 1)
qc.measure_all()
# Ideal simulation
ideal_sim = AerSimulator()
ideal_counts = ideal_sim.run(transpile(qc, ideal_sim), shots=n_shots).result().get_counts()
# Noisy simulation
noise_model = build_noise_model()
noisy_sim = AerSimulator(noise_model=noise_model)
noisy_counts = noisy_sim.run(transpile(qc, noisy_sim), shots=n_shots).result().get_counts()
# Compare fidelity
all_states = set(ideal_counts) | set(noisy_counts)
ideal_probs = {s: ideal_counts.get(s, 0) / n_shots for s in all_states}
noisy_probs = {s: noisy_counts.get(s, 0) / n_shots for s in all_states}
fidelity = sum(np.sqrt(ideal_probs[s] * noisy_probs[s]) for s in all_states) ** 2print(f"GHZ state ({n_qubits} qubits) fidelity (noisy vs ideal): {fidelity:.4f}")
Circuit Transpilation and Optimization
"""
Transpile a circuit to a specific basis gate set with optimization.
"""from qiskit import QuantumCircuit, transpile
from qiskit.transpiler import PassManager
from qiskit_aer import AerSimulator
deftranspile_and_compare(qc: QuantumCircuit) -> None:
"""Show circuit depth before and after transpilation."""
backend = AerSimulator()
print(f"Original circuit depth: {qc.depth()}")
print(f"Original gate count: {qc.count_ops()}")
for opt_level in [0, 1, 2, 3]:
qc_t = transpile(qc, backend=backend, optimization_level=opt_level)
print(f"Opt level {opt_level}: depth={qc_t.depth()}, "f"cx gates={qc_t.count_ops().get('cx', 0)}")
Troubleshooting
Error
Cause
Fix
ProviderError: No Backend found
Old IBM provider API
Use qiskit-ibm-runtime for real hardware
ValueError: The instruction save_statevector is not supported
Wrong simulator method
Use AerSimulator(method="statevector")
VQE not converging
Poor initial parameters or ansatz
Try SPSA optimizer; increase reps; use parameter shift gradient