Cirq Quantum Programming in Python — Deep Dive

Cirq’s Architecture

Cirq is built around immutable data structures. Circuits, moments, and operations are all frozen objects — modifications create new instances. This functional style makes circuit transformations predictable and composable.

Circuit Internals

import cirq

# GridQubits match Google's Sycamore chip layout
q0 = cirq.GridQubit(0, 0)
q1 = cirq.GridQubit(0, 1)
q2 = cirq.GridQubit(1, 0)

circuit = cirq.Circuit([
    cirq.Moment([cirq.H(q0), cirq.H(q2)]),        # Parallel H gates
    cirq.Moment([cirq.CNOT(q0, q1)]),               # Entangle 0-1
    cirq.Moment([cirq.CNOT(q2, q1)]),               # Entangle 2-1
    cirq.Moment([cirq.measure(q0, q1, q2, key='result')])
])

print(circuit)

The text representation shows a timeline diagram:

(0, 0): ───H───@───────M('result')───

(0, 1): ───────X───X───M─────────────

(1, 0): ───H───────@───M─────────────

Insert Strategies in Detail

circuit = cirq.Circuit()
circuit.append(cirq.H(q0), strategy=cirq.InsertStrategy.NEW)
circuit.append(cirq.H(q1), strategy=cirq.InsertStrategy.EARLIEST)
# H(q1) gets packed into the same moment as H(q0) since they don't conflict

circuit.append(cirq.CNOT(q0, q1), strategy=cirq.InsertStrategy.NEW_THEN_INLINE)

InsertStrategy.EARLIEST is crucial for circuit optimization — it minimizes depth by parallelizing non-conflicting operations.

Custom Gate Definition

Cirq allows defining custom gates with full control over their matrix representation:

import numpy as np

class RotationGate(cirq.Gate):
    def __init__(self, theta):
        super().__init__()
        self._theta = theta

    def _num_qubits_(self):
        return 1

    def _unitary_(self):
        c = np.cos(self._theta / 2)
        s = np.sin(self._theta / 2)
        return np.array([[c, -s], [s, c]])

    def _circuit_diagram_info_(self, args):
        return f'R({self._theta:.2f})'

    def _decompose_(self, qubits):
        # Decompose into native gates for compilation
        return [cirq.ry(self._theta).on(qubits[0])]

gate = RotationGate(np.pi / 4)
circuit = cirq.Circuit(gate.on(q0))

The _decompose_ method is important — it tells the compiler how to break your custom gate into primitives the hardware understands.

Simulation Deep Dive

Statevector Simulation

simulator = cirq.Simulator()

# Get the full statevector at each step
result = simulator.simulate(circuit)
print(result.final_state_vector)

# Step through the circuit moment by moment
for step in simulator.simulate_moment_steps(circuit):
    print(f"State after moment: {step.state_vector()}")

Density Matrix Simulation

For mixed states and noise modeling:

dm_simulator = cirq.DensityMatrixSimulator()
result = dm_simulator.simulate(circuit)
print(result.final_density_matrix)

Sampling (Mimicking Hardware)

result = simulator.run(circuit, repetitions=10000)
counts = result.histogram(key='result')
# {0: 2503, 3: 2497, 5: 2500, 6: 2500}  (example for GHZ-like state)

Noise Modeling

Cirq provides detailed noise channels for realistic simulation:

# Depolarizing noise after every gate
noise_model = cirq.ConstantQubitNoiseModel(
    qubit_noise_gate=cirq.depolarize(p=0.01)
)

noisy_simulator = cirq.DensityMatrixSimulator(noise=noise_model)
result = noisy_simulator.run(circuit, repetitions=10000)

Custom Noise Channels

class T1Decay(cirq.Gate):
    """Amplitude damping channel modeling T1 decay."""
    def __init__(self, gamma):
        self._gamma = gamma

    def _num_qubits_(self):
        return 1

    def _kraus_(self):
        g = self._gamma
        return [
            np.array([[1, 0], [0, np.sqrt(1 - g)]]),
            np.array([[0, np.sqrt(g)], [0, 0]])
        ]

# Apply after long idle periods
circuit_with_noise = cirq.Circuit([
    cirq.Moment([cirq.H(q0)]),
    cirq.Moment([T1Decay(0.02).on(q0)]),  # Idle decay
    cirq.Moment([cirq.measure(q0, key='m')])
])

Realistic Device Noise

For Google hardware, use cirq_google.engine.virtual_engine_factory to create a simulator with calibration-based noise models that mirror actual device behavior.

Hardware-Native Compilation

Google’s Sycamore processor natively supports the √iSWAP gate. All two-qubit operations must be decomposed into this:

import cirq_google

# Compile for Sycamore's native gate set
sycamore_gateset = cirq_google.SycamoreGateset()

# Validate that a circuit uses only native gates
try:
    sycamore_gateset.validate_circuit(circuit)
except ValueError as e:
    print(f"Non-native gates found: {e}")

# Use the built-in optimizer
optimized = cirq.optimize_for_target_gateset(
    circuit,
    gateset=cirq_google.SycamoreGateset()
)

Qubit Routing for Grid Topology

device = cirq_google.Sycamore
device_graph = device.metadata.qubit_set

# Map logical qubits to physical qubits
router = cirq.RouteCQC(device.metadata.nx_graph)
routed_circuit, mapping = router.route_circuit(circuit)

The router inserts SWAP gates to move qubit states between non-adjacent positions. Minimizing SWAP count is critical — each SWAP decomposes into three CX gates, each introducing noise.

Advanced Patterns

Quantum Error Detection Circuit

def bit_flip_detection():
    """Encodes one logical qubit into three physical qubits."""
    data = cirq.LineQubit(0)
    anc1, anc2 = cirq.LineQubit(1), cirq.LineQubit(2)

    circuit = cirq.Circuit([
        # Encode
        cirq.Moment([cirq.CNOT(data, anc1)]),
        cirq.Moment([cirq.CNOT(data, anc2)]),
        # Syndrome measurement (detect single bit flip)
        cirq.Moment([cirq.CNOT(data, cirq.LineQubit(3)),
                      cirq.CNOT(anc1, cirq.LineQubit(3))]),
        cirq.Moment([cirq.CNOT(anc1, cirq.LineQubit(4)),
                      cirq.CNOT(anc2, cirq.LineQubit(4))]),
        cirq.Moment([cirq.measure(cirq.LineQubit(3), cirq.LineQubit(4),
                                   key='syndrome')])
    ])
    return circuit

Variational Quantum Eigensolver (VQE) Structure

import sympy

def hardware_efficient_ansatz(qubits, depth):
    """Ansatz suitable for near-term hardware."""
    params = []
    moments = []

    for d in range(depth):
        # Single-qubit rotation layer
        layer_params = []
        single_ops = []
        for i, q in enumerate(qubits):
            theta = sympy.Symbol(f'θ_{d}_{i}')
            phi = sympy.Symbol(f'φ_{d}_{i}')
            layer_params.extend([theta, phi])
            single_ops.append(cirq.ry(theta).on(q))
        moments.append(cirq.Moment(single_ops))

        rz_ops = [cirq.rz(sympy.Symbol(f'φ_{d}_{i}')).on(q)
                   for i, q in enumerate(qubits)]
        moments.append(cirq.Moment(rz_ops))

        # Entangling layer — nearest-neighbor CZ
        entangle_ops = [cirq.CZ(qubits[i], qubits[i+1])
                        for i in range(0, len(qubits) - 1, 2)]
        moments.append(cirq.Moment(entangle_ops))

        params.extend(layer_params)

    return cirq.Circuit(moments), params

Parameter Resolution for Optimization

import scipy.optimize

def cost_function(param_values, circuit, params, simulator, target):
    resolver = cirq.ParamResolver(dict(zip(params, param_values)))
    result = simulator.simulate(circuit, param_resolver=resolver)
    # Compute expectation value against target observable
    state = result.final_state_vector
    return np.real(state.conj() @ target @ state)

# Classical optimization loop
initial_params = np.random.uniform(0, 2*np.pi, len(params))
result = scipy.optimize.minimize(
    cost_function,
    initial_params,
    args=(circuit, params, simulator, hamiltonian_matrix),
    method='COBYLA'
)

Performance and Scaling

Simulation Limits

QubitsMemory (statevector)Time (single simulation)
20~16 MB< 1 second
25~512 MBseconds
30~16 GBminutes
35~512 GBimpractical locally

For larger circuits, use Google’s cloud simulators (qsim) which leverage GPU acceleration and can handle 32-40 qubits efficiently.

qsim Integration

import qsimcirq

qsim_simulator = qsimcirq.QSimSimulator()
result = qsim_simulator.simulate(circuit)

qsim uses optimized C++ and GPU kernels, achieving 10-100x speedups over Cirq’s pure-Python simulator for circuits above 20 qubits.

Google Quantum AI Integration

import cirq_google

# Connect to Google's quantum computing service
engine = cirq_google.Engine(project_id='your-project')

# List available processors
processors = engine.list_processors()

# Submit a job
job = engine.run(
    program=circuit,
    processor_ids=['rainbow'],
    repetitions=10000
)
results = job.results()

Access requires a Google Cloud project with the Quantum Engine API enabled. Research access is available through Google’s quantum computing partnerships program.

Debugging Toolkit

# Circuit statistics
print(f"Depth: {len(circuit)}")
print(f"Gate count: {sum(1 for _ in circuit.all_operations())}")
print(f"Two-qubit gates: {sum(1 for op in circuit.all_operations() if len(op.qubits) == 2)}")

# Unitary of the entire circuit (small circuits only)
unitary = cirq.unitary(circuit)
print(f"Unitary shape: {unitary.shape}")

# Check if a circuit is valid for a device
device.validate_circuit(circuit)

One thing to remember: Cirq’s strength is its honest representation of quantum hardware constraints — moments model real gate timing, grid qubits model physical topology, and native gate compilation reflects actual device capabilities. This transparency is what makes it the framework of choice for pushing the boundaries of near-term quantum hardware.

pythonquantum-computingcirqgoogle

See Also