# Core Qiskit
pip install qiskit

# With visualization tools
pip install qiskit[visualization]

# Quantum chemistry extension
pip install qiskit-nature qiskit-nature-pyscf

# Machine learning extension
pip install qiskit-machine-learning

# Optimization extension
pip install qiskit-optimization

# Full installation
pip install 'qiskit[all]' qiskit-nature qiskit-machine-learning qiskit-optimization

# Core imports
from qiskit import QuantumCircuit, QuantumRegister, ClassicalRegister
from qiskit import transpile, assemble
from qiskit.providers.aer import AerSimulator
from qiskit.visualization import plot_histogram, plot_bloch_multivector

# Quantum algorithms
from qiskit.algorithms import VQE, QAOA, Grover, Shor
from qiskit.algorithms.optimizers import SLSQP, COBYLA, SPSA

# Quantum info
from qiskit.quantum_info import Statevector, DensityMatrix, Operator
from qiskit.quantum_info import entropy, entanglement_of_formation

# Circuit library
from qiskit.circuit.library import QFT, RealAmplitudes, EfficientSU2

from qiskit import QuantumCircuit
from qiskit.providers.aer import AerSimulator

# Create circuit
qc = QuantumCircuit(2, 2)

# Add gates
qc.h(0)              # Hadamard on qubit 0
qc.cx(0, 1)          # CNOT from 0 to 1

# Measure
qc.measure([0, 1], [0, 1])

# Simulate
simulator = AerSimulator()
job = simulator.run(qc, shots=1000)
result = job.result()
counts = result.get_counts()

print(f"Results: {counts}")

from qiskit import QuantumCircuit
from qiskit.algorithms import VQE
from qiskit.algorithms.optimizers import SLSQP
from qiskit.circuit.library import RealAmplitudes
from qiskit.primitives import Estimator
from qiskit.quantum_info import SparsePauliOp

# Define Hamiltonian
hamiltonian = SparsePauliOp(['ZZ', 'IZ', 'ZI'], coeffs=[1.0, -0.5, -0.5])

# Create ansatz
ansatz = RealAmplitudes(num_qubits=2, reps=1)

# Setup VQE
optimizer = SLSQP(maxiter=100)
estimator = Estimator()
vqe = VQE(estimator, ansatz, optimizer)

# Run
result = vqe.compute_minimum_eigenvalue(hamiltonian)
print(f"Ground state energy: {result.eigenvalue:.6f}")

from qiskit import QuantumCircuit
from qiskit.providers.aer import AerSimulator
from qiskit.primitives import Estimator

# ❌ BAD: No measurement
qc = QuantumCircuit(2, 2)
qc.h(0)
qc.cx(0, 1)
# Forgot qc.measure()!

# ✅ GOOD: Always measure when needed
qc = QuantumCircuit(2, 2)
qc.h(0)
qc.cx(0, 1)
qc.measure([0, 1], [0, 1])

# ❌ BAD: Using deprecated execute()
from qiskit import execute
result = execute(qc, backend, shots=1024).result()

# ✅ GOOD: Use new run() method
simulator = AerSimulator()
job = simulator.run(qc, shots=1024)
result = job.result()

# ❌ BAD: Assuming perfect measurement
counts = result.get_counts()
# Assuming exactly 50/50 split!
assert counts['00'] == 512

# ✅ GOOD: Handle statistical variation
counts = result.get_counts()
ratio = counts.get('00', 0) / sum(counts.values())
print(f"Measured |00⟩ with probability {ratio:.3f}")

# ❌ BAD: Not checking circuit properties
qc = QuantumCircuit(20)  # Many qubits!
# Adding many gates...
# Trying to simulate without checking depth/size!

# ✅ GOOD: Check circuit properties
qc = QuantumCircuit(20)
# ... add gates ...
print(f"Circuit depth: {qc.depth()}")
print(f"Gate count: {len(qc.data)}")
print(f"Qubits: {qc.num_qubits}")

# ❌ BAD: Ignoring transpilation
job = backend.run(qc)  # May fail on real hardware!

# ✅ GOOD: Transpile for backend
from qiskit import transpile
transpiled_qc = transpile(qc, backend=backend, optimization_level=3)
job = backend.run(transpiled_qc)

from qiskit import QuantumCircuit, QuantumRegister, ClassicalRegister
import numpy as np

# Method 1: Simple initialization
qc = QuantumCircuit(3, 3)  # 3 qubits, 3 classical bits

# Method 2: Using registers
qr = QuantumRegister(3, 'q')
cr = ClassicalRegister(3, 'c')
qc = QuantumCircuit(qr, cr)

# Method 3: Multiple registers
qr1 = QuantumRegister(2, 'data')
qr2 = QuantumRegister(1, 'ancilla')
cr = ClassicalRegister(2, 'meas')
qc = QuantumCircuit(qr1, qr2, cr)

print(f"Number of qubits: {qc.num_qubits}")
print(f"Number of classical bits: {qc.num_clbits}")
print(f"Circuit depth: {qc.depth()}")

from qiskit import QuantumCircuit
import numpy as np

qc = QuantumCircuit(1)

# Pauli gates
qc.x(0)   # Pauli X (NOT gate)
qc.y(0)   # Pauli Y
qc.z(0)   # Pauli Z

# Hadamard gate
qc.h(0)   # Creates superposition

# Phase gates
qc.s(0)   # S gate (π/2 phase)
qc.t(0)   # T gate (π/4 phase)
qc.sdg(0) # S dagger
qc.tdg(0) # T dagger

# Rotation gates
qc.rx(np.pi/4, 0)  # Rotation around X
qc.ry(np.pi/4, 0)  # Rotation around Y
qc.rz(np.pi/4, 0)  # Rotation around Z

# General rotation
qc.u(np.pi/4, np.pi/2, np.pi, 0)  # U gate

# Identity (wait)
qc.id(0)

print(f"Gate count: {len(qc.data)}")

from qiskit import QuantumCircuit
import numpy as np

qc = QuantumCircuit(2)

# CNOT (Controlled-NOT)
qc.cx(0, 1)  # Control: 0, Target: 1

# Other controlled gates
qc.cy(0, 1)  # Controlled-Y
qc.cz(0, 1)  # Controlled-Z
qc.ch(0, 1)  # Controlled-Hadamard

# SWAP gate
qc.swap(0, 1)

# Controlled phase
qc.cp(np.pi/4, 0, 1)

# Controlled-U
qc.cu(np.pi/4, np.pi/2, np.pi, 0, 0, 1)

# Toffoli (CCX) - needs 3 qubits
qc_3 = QuantumCircuit(3)
qc_3.ccx(0, 1, 2)  # Controls: 0,1, Target: 2

print(qc.draw())

from qiskit import QuantumCircuit
from qiskit.providers.aer import AerSimulator
from qiskit.quantum_info import Statevector

# Bell state (maximally entangled)
def create_bell_state():
    qc = QuantumCircuit(2)
    qc.h(0)
    qc.cx(0, 1)
    return qc

bell = create_bell_state()
state = Statevector.from_instruction(bell)
print(f"Bell state: {state}")

# GHZ state (3-qubit entanglement)
def create_ghz_state(n):
    qc = QuantumCircuit(n)
    qc.h(0)
    for i in range(n-1):
        qc.cx(i, i+1)
    return qc

ghz = create_ghz_state(3)
state_ghz = Statevector.from_instruction(ghz)
print(f"GHZ state: {state_ghz}")

# W state (another type of 3-qubit entanglement)
def create_w_state():
    qc = QuantumCircuit(3)
    qc.ry(1.9106, 0)
    qc.ch(0, 1)
    qc.x(0)
    qc.cy(0, 1)
    qc.ccx(0, 1, 2)
    qc.x(0)
    return qc

w = create_w_state()
print(f"W state circuit depth: {w.depth()}")

from qiskit import QuantumCircuit
from qiskit.circuit import Parameter, ParameterVector
import numpy as np

# Single parameter
theta = Parameter('θ')
qc = QuantumCircuit(1)
qc.ry(theta, 0)

# Bind parameter
bound_qc = qc.bind_parameters({theta: np.pi/4})
print(f"Unbound: {qc}")
print(f"Bound: {bound_qc}")

# Multiple parameters
params = ParameterVector('θ', 4)
qc_param = QuantumCircuit(2)
qc_param.ry(params[0], 0)
qc_param.ry(params[1], 1)
qc_param.cx(0, 1)
qc_param.ry(params[2], 0)
qc_param.ry(params[3], 1)

# Bind all parameters
values = [np.pi/4, np.pi/3, np.pi/2, np.pi/6]
bound = qc_param.bind_parameters(dict(zip(params, values)))

print(f"Number of parameters: {qc_param.num_parameters}")

from qiskit import QuantumCircuit
from qiskit.algorithms import VQE
from qiskit.algorithms.optimizers import SLSQP, COBYLA
from qiskit.circuit.library import RealAmplitudes, EfficientSU2
from qiskit.primitives import Estimator
from qiskit.quantum_info import SparsePauliOp
import numpy as np

# Define Hamiltonian (e.g., H2 molecule)
# H = -1.05 * ZZ + 0.39 * XX - 0.39 * YY - 0.01 * ZI
hamiltonian = SparsePauliOp(
    ['ZZ', 'XX', 'YY', 'ZI', 'IZ'],
    coeffs=[-1.05, 0.39, -0.39, -0.01, -0.01]
)

# Create ansatz (variational form)
ansatz = RealAmplitudes(num_qubits=2, reps=2)

# Alternative: EfficientSU2
# ansatz = EfficientSU2(num_qubits=2, reps=2)

# Setup optimizer
optimizer = SLSQP(maxiter=100)

# Create estimator primitive
estimator = Estimator()

# Setup and run VQE
vqe = VQE(estimator, ansatz, optimizer)
result = vqe.compute_minimum_eigenvalue(hamiltonian)

print(f"Ground state energy: {result.eigenvalue:.6f}")
print(f"Optimal parameters: {result.optimal_parameters}")
print(f"Optimizer evaluations: {result.cost_function_evals}")

# Get optimal circuit
optimal_circuit = ansatz.bind_parameters(result.optimal_point)
print(f"\nOptimal circuit depth: {optimal_circuit.depth()}")

from qiskit import QuantumCircuit
from qiskit.algorithms import QAOA
from qiskit.algorithms.optimizers import COBYLA
from qiskit.primitives import Sampler
from qiskit.quantum_info import SparsePauliOp
import numpy as np

# Max-Cut problem on a triangle graph
# Cost Hamiltonian: H = (1-Z0*Z1)/2 + (1-Z1*Z2)/2 + (1-Z0*Z2)/2
cost_hamiltonian = SparsePauliOp(
    ['ZZ', 'ZI', 'IZ'],
    coeffs=[0.5, -0.5, -0.5]
)

# Setup QAOA
optimizer = COBYLA(maxiter=100)
sampler = Sampler()
qaoa = QAOA(sampler, optimizer, reps=2)

# Run QAOA
result = qaoa.compute_minimum_eigenvalue(cost_hamiltonian)

print(f"Optimal objective: {result.eigenvalue:.6f}")
print(f"Optimal parameters: {result.optimal_parameters}")

# Most probable solution
from qiskit.result import QuasiDistribution
prob_dist = result.eigenstate
if isinstance(prob_dist, QuasiDistribution):
    most_likely = max(prob_dist, key=prob_dist.get)
    print(f"Most likely solution: {bin(most_likely)[2:].zfill(2)}")

from qiskit import QuantumCircuit
from qiskit.algorithms import Grover, AmplificationProblem
from qiskit.circuit.library import GroverOperator
from qiskit.primitives import Sampler
import numpy as np

def grover_search(marked_states, n_qubits):
    """
    Grover's algorithm to find marked states.
    
    Args:
        marked_states: List of marked states (e.g., ['101', '110'])
        n_qubits: Number of qubits
    """
    # Create oracle that marks the target states
    oracle = QuantumCircuit(n_qubits)
    
    # Mark state |101⟩ by flipping phase
    if '101' in marked_states:
        oracle.x([0, 2])  # Flip to make 101 -> 111
        oracle.h(2)
        oracle.ccx(0, 1, 2)
        oracle.h(2)
        oracle.x([0, 2])  # Flip back
    
    # Grover diffusion operator
    diffusion = QuantumCircuit(n_qubits)
    diffusion.h(range(n_qubits))
    diffusion.x(range(n_qubits))
    diffusion.h(n_qubits - 1)
    diffusion.mcx(list(range(n_qubits - 1)), n_qubits - 1)
    diffusion.h(n_qubits - 1)
    diffusion.x(range(n_qubits))
    diffusion.h(range(n_qubits))
    
    # Complete Grover iteration
    grover_op = oracle.compose(diffusion)
    
    # Setup and run
    problem = AmplificationProblem(oracle, is_good_state=marked_states)
    grover = Grover(sampler=Sampler())
    result = grover.amplify(problem)
    
    return result

# Example: Search for |101⟩ in 3-qubit space
result = grover_search(['101'], 3)
print(f"Top measurement: {result.top_measurement}")
print(f"Oracle evaluations: {result.oracle_evaluation}")

# Simpler example: 2-qubit search
qc = QuantumCircuit(2, 2)
qc.h([0, 1])  # Superposition

# Oracle marks |11⟩
qc.cz(0, 1)

# Diffusion
qc.h([0, 1])
qc.x([0, 1])
qc.cz(0, 1)
qc.x([0, 1])
qc.h([0, 1])

qc.measure([0, 1], [0, 1])

from qiskit.providers.aer import AerSimulator
simulator = AerSimulator()
job = simulator.run(qc, shots=1000)
counts = job.result().get_counts()
print(f"Grover results: {counts}")

from qiskit import QuantumCircuit
from qiskit.circuit.library import QFT
from qiskit.providers.aer import AerSimulator
import numpy as np

def quantum_phase_estimation(unitary, n_counting_qubits):
    """
    Estimate the phase of an eigenstate of a unitary operator.
    
    Args:
        unitary: Unitary operator as a QuantumCircuit
        n_counting_qubits: Precision qubits
    """
    n_system_qubits = unitary.num_qubits
    
    qc = QuantumCircuit(n_counting_qubits + n_system_qubits, n_counting_qubits)
    
    # Initialize eigenstate (simplified: use |1⟩)
    qc.x(n_counting_qubits)
    
    # Hadamard on counting qubits
    for i in range(n_counting_qubits):
        qc.h(i)
    
    # Controlled-U^(2^i) operations
    repetitions = 1
    for counting_qubit in range(n_counting_qubits):
        for _ in range(repetitions):
            qc.append(unitary.control(), [counting_qubit] + 
                     list(range(n_counting_qubits, n_counting_qubits + n_system_qubits)))
        repetitions *= 2
    
    # Inverse QFT
    qc.append(QFT(n_counting_qubits, inverse=True), range(n_counting_qubits))
    
    # Measure counting qubits
    qc.measure(range(n_counting_qubits), range(n_counting_qubits))
    
    return qc

# Example: Phase estimation for T gate (phase π/4)
t_gate = QuantumCircuit(1)
t_gate.t(0)

qpe_circuit = quantum_phase_estimation(t_gate, n_counting_qubits=3)

simulator = AerSimulator()
job = simulator.run(qpe_circuit, shots=1000)
counts = job.result().get_counts()
print(f"Phase estimation results: {counts}")

# Most common result gives phase
most_common = max(counts, key=counts.get)
phase_estimate = int(most_common, 2) / (2**3)
print(f"Estimated phase: {phase_estimate * 2 * np.pi:.4f} rad")
print(f"True phase: {np.pi/4:.4f} rad")

# Requires: pip install qiskit-nature qiskit-nature-pyscf

from qiskit_nature.units import DistanceUnit
from qiskit_nature.second_q.drivers import PySCFDriver
from qiskit_nature.second_q.mappers import JordanWignerMapper, ParityMapper
from qiskit_nature.second_q.circuit.library import HartreeFock, UCCSD
from qiskit.algorithms import VQE
from qiskit.algorithms.optimizers import SLSQP
from qiskit.primitives import Estimator
import numpy as np

# Define molecule (H2)
driver = PySCFDriver(
    atom='H 0 0 0; H 0 0 0.735',
    basis='sto3g',
    charge=0,
    spin=0,
    unit=DistanceUnit.ANGSTROM
)

# Get problem
problem = driver.run()
hamiltonian = problem.hamiltonian.second_q_op()

# Map to qubits
mapper = JordanWignerMapper()
qubit_op = mapper.map(hamiltonian)

print(f"Number of qubits required: {qubit_op.num_qubits}")

# Initial state (Hartree-Fock)
num_particles = problem.num_particles
num_spatial_orbitals = problem.num_spatial_orbitals

init_state = HartreeFock(
    num_spatial_orbitals=num_spatial_orbitals,
    num_particles=num_particles,
    qubit_mapper=mapper
)

# Ansatz (UCCSD)
ansatz = UCCSD(
    num_spatial_orbitals=num_spatial_orbitals,
    num_particles=num_particles,
    qubit_mapper=mapper,
    initial_state=init_state
)

# VQE calculation
optimizer = SLSQP(maxiter=100)
estimator = Estimator()

vqe = VQE(estimator, ansatz, optimizer)
result = vqe.compute_minimum_eigenvalue(qubit_op)

print(f"\nVQE Ground State Energy: {result.eigenvalue:.6f} Ha")
print(f"Reference HF energy: {problem.reference_energy:.6f} Ha")

from qiskit_nature.second_q.algorithms import GroundStateEigensolver, ExcitedStatesEigensolver
from qiskit_nature.second_q.algorithms import QEOM
from qiskit_nature.second_q.mappers import JordanWignerMapper
from qiskit_nature.second_q.drivers import PySCFDriver
from qiskit.algorithms import VQE
from qiskit.algorithms.optimizers import SLSQP
from qiskit.primitives import Estimator

# Setup molecule
driver = PySCFDriver(atom='H 0 0 0; H 0 0 0.735', basis='sto3g')
problem = driver.run()

# Mapper
mapper = JordanWignerMapper()

# Ground state solver
optimizer = SLSQP(maxiter=100)
estimator = Estimator()

# Simple ansatz for demonstration
from qiskit.circuit.library import RealAmplitudes
ansatz = RealAmplitudes(num_qubits=4, reps=2)

vqe = VQE(estimator, ansatz, optimizer)
ground_solver = GroundStateEigensolver(mapper, vqe)

# Calculate ground state
ground_result = ground_solver.solve(problem)
print(f"Ground state energy: {ground_result.total_energies[0]:.6f} Ha")

# Excited states with QEOM
qeom = QEOM(ground_solver, 'sd')  # singles and doubles
excited_solver = ExcitedStatesEigensolver(mapper, qeom)

# Calculate excited states
excited_result = excited_solver.solve(problem)
print(f"\nExcited state energies:")
for i, energy in enumerate(excited_result.total_energies):
    print(f"State {i}: {energy:.6f} Ha")

from qiskit_nature.second_q.drivers import PySCFDriver
from qiskit_nature.second_q.mappers import JordanWignerMapper
from qiskit_nature.second_q.algorithms import GroundStateEigensolver
from qiskit.algorithms import VQE
from qiskit.algorithms.optimizers import SLSQP
from qiskit.primitives import Estimator
from qiskit.circuit.library import RealAmplitudes
import numpy as np

def calculate_pes(distances, molecule_template, basis='sto3g'):
    """
    Calculate potential energy surface for a diatomic molecule.
    
    Args:
        distances: Array of interatomic distances
        molecule_template: Molecule string with {} for distance
        basis: Basis set
    """
    energies = []
    
    for distance in distances:
        # Setup molecule at this distance
        atom_string = molecule_template.format(distance)
        driver = PySCFDriver(atom=atom_string, basis=basis)
        problem = driver.run()
        
        # Map to qubits
        mapper = JordanWignerMapper()
        
        # Simple VQE calculation
        ansatz = RealAmplitudes(num_qubits=4, reps=1)
        optimizer = SLSQP(maxiter=50)
        estimator = Estimator()
        
        vqe = VQE(estimator, ansatz, optimizer)
        solver = GroundStateEigensolver(mapper, vqe)
        
        # Calculate energy
        result = solver.solve(problem)
        energies.append(result.total_energies[0])
        
        print(f"Distance {distance:.2f} Å: Energy = {energies[-1]:.6f} Ha")
    
    return np.array(energies)

# Calculate PES for H2
distances = np.linspace(0.5, 2.5, 5)  # Angstroms
molecule_template = 'H 0 0 0; H 0 0 {}'

energies = calculate_pes(distances, molecule_template)

# Find equilibrium bond length
min_idx = np.argmin(energies)
print(f"\nEquilibrium bond length: {distances[min_idx]:.2f} Å")
print(f"Minimum energy: {energies[min_idx]:.6f} Ha")

from qiskit import QuantumCircuit
from qiskit.quantum_info import SparsePauliOp, Statevector
from qiskit.synthesis import SuzukiTrotter
import numpy as np

def trotter_simulation(hamiltonian, time, n_steps):
    """
    Simulate time evolution under Hamiltonian using Trotter decomposition.
    
    Args:
        hamiltonian: SparsePauliOp Hamiltonian
        time: Total evolution time
        n_steps: Number of Trotter steps
    """
    dt = time / n_steps
    
    # Create Trotter circuit
    n_qubits = hamiltonian.num_qubits
    qc = QuantumCircuit(n_qubits)
    
    # Initial state (e.g., |0⟩)
    # Apply Trotter steps
    for _ in range(n_steps):
        # For each Pauli term in Hamiltonian
        for pauli, coeff in zip(hamiltonian.paulis, hamiltonian.coeffs):
            # Apply exp(-i * coeff * dt * pauli)
            angle = 2 * float(coeff.real) * dt
            
            # Simple implementation for single Pauli terms
            pauli_str = str(pauli)
            
            # Apply appropriate rotation
            if 'Z' in pauli_str and pauli_str.count('I') == n_qubits - 1:
                idx = pauli_str.index('Z')
                qc.rz(angle, n_qubits - 1 - idx)
            elif 'X' in pauli_str and pauli_str.count('I') == n_qubits - 1:
                idx = pauli_str.index('X')
                qc.rx(angle, n_qubits - 1 - idx)
    
    return qc

# Example: Ising model Hamiltonian
hamiltonian = SparsePauliOp(['ZZ', 'XI', 'IX'], coeffs=[1.0, -0.5, -0.5])

# Simulate for time t=1.0 with 10 Trotter steps
qc = trotter_simulation(hamiltonian, time=1.0, n_steps=10)

print(f"Trotter circuit depth: {qc.depth()}")
print(f"Gate count: {len(qc.data)}")

# Get final state
state = Statevector.from_instruction(qc)
print(f"Final state vector: {state[:4]}")  # First 4 amplitudes

from qiskit import QuantumCircuit
from qiskit.algorithms import VQE
from qiskit.algorithms.optimizers import SLSQP
from qiskit.circuit.library import RealAmplitudes
from qiskit.primitives import Estimator
from qiskit.quantum_info import SparsePauliOp
import numpy as np

def simulate_dynamics(hamiltonian, initial_state_circuit, times, ansatz_reps=2):
    """
    Simulate quantum dynamics using VQE at different times.
    
    Args:
        hamiltonian: Time-evolution Hamiltonian
        initial_state_circuit: Circuit preparing initial state
        times: Array of times to simulate
        ansatz_reps: Repetitions in ansatz
    """
    n_qubits = hamiltonian.num_qubits
    results = []
    
    for t in times:
        # Time-evolved Hamiltonian: H(t) = exp(-iHt) ρ(0) exp(iHt)
        # Approximate with VQE
        
        ansatz = RealAmplitudes(num_qubits=n_qubits, reps=ansatz_reps)
        
        # Combine initial state + ansatz
        full_circuit = initial_state_circuit.copy()
        full_circuit.compose(ansatz, inplace=True)
        
        optimizer = SLSQP(maxiter=100)
        estimator = Estimator()
        
        vqe = VQE(estimator, full_circuit, optimizer)
        result = vqe.compute_minimum_eigenvalue(hamiltonian)
        
        results.append({
            'time': t,
            'energy': result.eigenvalue,
            'parameters': result.optimal_parameters
        })
        
        print(f"Time {t:.2f}: Energy = {result.eigenvalue:.6f}")
    
    return results

# Example: Simulate Heisenberg model
hamiltonian = SparsePauliOp(
    ['XX', 'YY', 'ZZ'],
    coeffs=[0.5, 0.5, 0.5]
)

# Initial state: |01⟩
initial_state = QuantumCircuit(2)
initial_state.x(1)

# Simulate at different times
times = np.linspace(0, 2, 5)
results = simulate_dynamics(hamiltonian, initial_state, times)

# Requires: pip install qiskit-machine-learning

from qiskit_machine_learning.kernels import FidelityQuantumKernel
from qiskit.circuit.library import ZZFeatureMap
from qiskit.primitives import Sampler
from sklearn.svm import SVC
import numpy as np

# Generate synthetic data
np.random.seed(42)
X_train = np.random.rand(20, 2) * 2 * np.pi
y_train = (X_train[:, 0] + X_train[:, 1] > np.pi).astype(int)

X_test = np.random.rand(10, 2) * 2 * np.pi
y_test = (X_test[:, 0] + X_test[:, 1] > np.pi).astype(int)

# Create feature map
feature_map = ZZFeatureMap(feature_dimension=2, reps=2)

# Create quantum kernel
sampler = Sampler()
quantum_kernel = FidelityQuantumKernel(feature_map=feature_map, fidelity=sampler)

# Compute kernel matrix
kernel_matrix_train = quantum_kernel.evaluate(x_vec=X_train)
kernel_matrix_test = quantum_kernel.evaluate(x_vec=X_test, y_vec=X_train)

print(f"Training kernel matrix shape: {kernel_matrix_train.shape}")
print(f"Test kernel matrix shape: {kernel_matrix_test.shape}")

# Use with SVM
svc = SVC(kernel='precomputed')
svc.fit(kernel_matrix_train, y_train)

# Predict
predictions = svc.predict(kernel_matrix_test)
accuracy = np.mean(predictions == y_test)

print(f"\nQuantum Kernel SVM Accuracy: {accuracy:.2%}")

from qiskit_machine_learning.algorithms import VQC
from qiskit.circuit.library import RealAmplitudes, ZZFeatureMap
from qiskit.algorithms.optimizers import COBYLA
from qiskit.primitives import Sampler
import numpy as np

# Generate training data
np.random.seed(42)
X_train = np.random.rand(40, 2) * 2 - 1  # [-1, 1]
y_train = (X_train[:, 0]**2 + X_train[:, 1]**2 < 0.5).astype(int)

X_test = np.random.rand(20, 2) * 2 - 1
y_test = (X_test[:, 0]**2 + X_test[:, 1]**2 < 0.5).astype(int)

# Feature map
feature_map = ZZFeatureMap(2)

# Ansatz
ansatz = RealAmplitudes(2, reps=2)

# Create VQC
optimizer = COBYLA(maxiter=100)
sampler = Sampler()

vqc = VQC(
    sampler=sampler,
    feature_map=feature_map,
    ansatz=ansatz,
    optimizer=optimizer,
)

# Train
vqc.fit(X_train, y_train)

# Predict
train_score = vqc.score(X_train, y_train)
test_score = vqc.score(X_test, y_test)

print(f"Training accuracy: {train_score:.2%}")
print(f"Test accuracy: {test_score:.2%}")

# Predictions
predictions = vqc.predict(X_test)
print(f"Sample predictions: {predictions[:5]}")
print(f"True labels: {y_test[:5]}")

from qiskit_machine_learning.neural_networks import SamplerQNN
from qiskit_machine_learning.algorithms import NeuralNetworkClassifier
from qiskit.circuit import QuantumCircuit, Parameter
from qiskit.circuit.library import RealAmplitudes
from qiskit.algorithms.optimizers import COBYLA
from qiskit.primitives import Sampler
import numpy as np

# Create parameterized quantum circuit
def create_qnn_circuit(num_qubits, num_features):
    qc = QuantumCircuit(num_qubits)
    
    # Feature map parameters
    feature_params = [Parameter(f'x{i}') for i in range(num_features)]
    
    # Encode features
    for i, param in enumerate(feature_params):
        if i < num_qubits:
            qc.ry(param, i)
    
    # Variational part
    ansatz = RealAmplitudes(num_qubits, reps=1)
    qc.compose(ansatz, inplace=True)
    
    return qc, feature_params, ansatz.parameters

num_qubits = 2
num_features = 2

qc, feature_params, weight_params = create_qnn_circuit(num_qubits, num_features)

# Create QNN
sampler = Sampler()
qnn = SamplerQNN(
    circuit=qc,
    input_params=feature_params,
    weight_params=weight_params,
    sampler=sampler,
)

# Create classifier
optimizer = COBYLA(maxiter=50)
classifier = NeuralNetworkClassifier(
    neural_network=qnn,
    optimizer=optimizer,
)

# Training data
X_train = np.random.rand(30, 2)
y_train = (X_train[:, 0] > X_train[:, 1]).astype(int)

# Train
classifier.fit(X_train, y_train)

# Test
X_test = np.random.rand(10, 2)
y_test = (X_test[:, 0] > X_test[:, 1]).astype(int)

score = classifier.score(X_test, y_test)
print(f"QNN Classifier accuracy: {score:.2%}")

from qiskit.providers.aer import AerSimulator
from qiskit.providers.aer.noise import NoiseModel
from qiskit.providers.aer.noise import depolarizing_error, thermal_relaxation_error
from qiskit import QuantumCircuit
import numpy as np

# Create noise model
noise_model = NoiseModel()

# Depolarizing error on single-qubit gates
error_1q = depolarizing_error(0.001, 1)  # 0.1% error
noise_model.add_all_qubit_quantum_error(error_1q, ['u1', 'u2', 'u3', 'rx', 'ry', 'rz'])

# Depolarizing error on two-qubit gates
error_2q = depolarizing_error(0.01, 2)  # 1% error
noise_model.add_all_qubit_quantum_error(error_2q, ['cx', 'cz', 'swap'])

# Thermal relaxation error
# T1 = 50 microseconds, T2 = 70 microseconds, gate_time = 0.1 microseconds
t1 = 50e3  # ns
t2 = 70e3  # ns
gate_time = 100  # ns

error_thermal = thermal_relaxation_error(t1, t2, gate_time)
noise_model.add_all_qubit_quantum_error(error_thermal, ['id'])

# Measurement error
from qiskit.providers.aer.noise import ReadoutError
prob_meas0_prep1 = 0.05  # Prob of measuring 0 when prepared in 1
prob_meas1_prep0 = 0.10  # Prob of measuring 1 when prepared in 0

readout_error = ReadoutError([[1 - prob_meas1_prep0, prob_meas1_prep0],
                               [prob_meas0_prep1, 1 - prob_meas0_prep1]])
noise_model.add_all_qubit_readout_error(readout_error)

print(f"Noise model: {noise_model}")

# Run circuit with noise
qc = QuantumCircuit(2, 2)
qc.h(0)
qc.cx(0, 1)
qc.measure([0, 1], [0, 1])

# Noiseless simulation
simulator_ideal = AerSimulator()
job_ideal = simulator_ideal.run(qc, shots=1000)
counts_ideal = job_ideal.result().get_counts()

# Noisy simulation
simulator_noisy = AerSimulator(noise_model=noise_model)
job_noisy = simulator_noisy.run(qc, shots=1000)
counts_noisy = job_noisy.result().get_counts()

print(f"\nIdeal results: {counts_ideal}")
print(f"Noisy results: {counts_noisy}")

from qiskit import QuantumCircuit, transpile
from qiskit.providers.aer import AerSimulator
from qiskit.providers.aer.noise import NoiseModel, depolarizing_error
import numpy as np
from scipy.optimize import curve_fit

def run_circuit_with_noise(circuit, noise_level, shots=1000):
    """Run circuit with scaled noise."""
    # Create noise model
    noise_model = NoiseModel()
    error = depolarizing_error(noise_level, 1)
    noise_model.add_all_qubit_quantum_error(error, ['u1', 'u2', 'u3'])
    
    # Simulate
    simulator = AerSimulator(noise_model=noise_model)
    job = simulator.run(circuit, shots=shots)
    result = job.result()
    
    # Extract expectation value (simplified)
    counts = result.get_counts()
    expectation = sum(counts.get(bitstring, 0) * (-1)**bitstring.count('1') 
                     for bitstring in counts) / shots
    
    return expectation

def zero_noise_extrapolation(circuit, noise_levels, shots=1000):
    """
    Perform zero-noise extrapolation.
    
    Args:
        circuit: Quantum circuit
        noise_levels: List of noise scaling factors
        shots: Number of shots per simulation
    """
    expectations = []
    
    for noise in noise_levels:
        exp_val = run_circuit_with_noise(circuit, noise, shots)
        expectations.append(exp_val)
        print(f"Noise {noise:.4f}: Expectation = {exp_val:.4f}")
    
    # Fit exponential model: E(λ) = A + B * exp(-C * λ)
    def exp_model(x, a, b, c):
        return a + b * np.exp(-c * x)
    
    # Fit
    popt, _ = curve_fit(exp_model, noise_levels, expectations)
    
    # Extrapolate to zero noise
    zero_noise_value = exp_model(0, *popt)
    
    print(f"\nZero-noise extrapolated value: {zero_noise_value:.4f}")
    
    return zero_noise_value, expectations

# Example circuit
qc = QuantumCircuit(2, 2)
qc.h(0)
qc.cx(0, 1)
qc.measure([0, 1], [0, 1])

# Run ZNE
noise_levels = np.array([0.001, 0.005, 0.01, 0.02])
zne_value, measured_values = zero_noise_extrapolation(qc, noise_levels)

from qiskit import QuantumCircuit
from qiskit.providers.aer import AerSimulator
from qiskit.providers.aer.noise import NoiseModel, ReadoutError
from qiskit.result import marginal_counts
from qiskit.ignis.mitigation.measurement import (
    complete_meas_cal,
    CompleteMeasFitter
)

# Create readout noise
prob_meas0_prep1 = 0.1
prob_meas1_prep0 = 0.05

readout_error = ReadoutError([[1 - prob_meas1_prep0, prob_meas1_prep0],
                               [prob_meas0_prep1, 1 - prob_meas0_prep1]])

noise_model = NoiseModel()
noise_model.add_all_qubit_readout_error(readout_error)

# Generate calibration circuits
n_qubits = 2
meas_calibs, state_labels = complete_meas_cal(qr=list(range(n_qubits)))

# Run calibration
simulator = AerSimulator(noise_model=noise_model)
cal_results = []

for circuit in meas_calibs:
    job = simulator.run(circuit, shots=1000)
    cal_results.append(job.result())

# Fit the calibration
meas_fitter = CompleteMeasFitter(cal_results, state_labels)

# Get mitigation matrix
print("Calibration matrix:")
print(meas_fitter.cal_matrix)

# Run actual circuit
qc = QuantumCircuit(2, 2)
qc.h(0)
qc.cx(0, 1)
qc.measure([0, 1], [0, 1])

job = simulator.run(qc, shots=1000)
result = job.result()
noisy_counts = result.get_counts()

# Apply mitigation
mitigated_counts = meas_fitter.filter.apply(noisy_counts)

print(f"\nNoisy counts: {noisy_counts}")
print(f"Mitigated counts: {mitigated_counts}")