Thermal noise model (quantumlib#4673)

* Remove files for subsequent PRs

* Add coverage for OpId

* Clarify if-case

* snake_case_gamma

* Add insertion noise model.

* Update to use new match methods

* Apply review cleanup

* Introduce Thermal noise model.

* Align with checks

* remove ThermalChannel ref

* Answer review comments

* Document mutli-type behavior.

* qubit_dims -> qubits

* Extract methods, clean docstrings

* Cache kraus_ops_from_rates

* Document edge behaviors.

* Handle other qubits nicely

Co-authored-by: Cirq Bot <[email protected]>

Author: 95-martin-orion

cirq/devices/__init__.py

@@ -52,14 +52,17 @@
     InsertionNoiseModel,
 )
 
+from cirq.devices.thermal_noise_model import (
+    ThermalNoiseModel,
+)
+
 from cirq.devices.noise_utils import (
     OpIdentifier,
     decay_constant_to_xeb_fidelity,
     decay_constant_to_pauli_error,
     pauli_error_to_decay_constant,
     xeb_fidelity_to_decay_constant,
     pauli_error_from_t1,
-    pauli_error_from_depolarization,
     average_error,
     decoherence_pauli_error,
 )

cirq/devices/insertion_noise_model.py

@@ -15,7 +15,7 @@
 import dataclasses
 from typing import TYPE_CHECKING, Dict, List, Optional, Sequence
 
-from cirq import devices, ops
+from cirq import devices
 from cirq.devices import noise_utils
 
 if TYPE_CHECKING:
@@ -32,7 +32,11 @@ class InsertionNoiseModel(devices.NoiseModel):
 
     Args:
         ops_added: a map of gate types (and optionally, qubits they act on) to
-            operations that should be added.
+            operations that should be added. If two gate types provided apply
+            to a target gate, the most specific type will match; if neither
+            type is more specific (e.g. A is a subtype of B, but B defines
+            qubits and A does not) then the first one appering in this dict
+            will match.
         prepend: whether to add the new moment before the current one.
         require_physical_tag: whether to only apply noise to operations tagged
             with PHYSICAL_GATE_TAG.
@@ -63,6 +67,10 @@ def noisy_moment(
                 noise_ops.append(self.ops_added[match_id])
         if not noise_ops:
             return [moment]
+
+        from cirq import circuits
+
+        noise_steps = circuits.Circuit(noise_ops)
         if self.prepend:
-            return [ops.Moment(noise_ops), moment]
-        return [moment, ops.Moment(noise_ops)]
+            return [*noise_steps.moments, moment]
+        return [moment, *noise_steps.moments]

cirq/devices/insertion_noise_model_test.py

@@ -51,6 +51,20 @@ def test_insertion_noise():
     assert model.noisy_moment(moment_3, system_qubits=[q0, q1]) == [moment_3]
 
 
+def test_colliding_noise_qubits():
+    # Check that noise affecting other qubits doesn't cause issues.
+    q0, q1, q2, q3 = cirq.LineQubit.range(4)
+    op_id0 = OpIdentifier(cirq.CZPowGate)
+    model = InsertionNoiseModel({op_id0: cirq.CNOT(q1, q2)}, require_physical_tag=False)
+
+    moment_0 = cirq.Moment(cirq.CZ(q0, q1), cirq.CZ(q2, q3))
+    assert model.noisy_moment(moment_0, system_qubits=[q0, q1, q2, q3]) == [
+        moment_0,
+        cirq.Moment(cirq.CNOT(q1, q2)),
+        cirq.Moment(cirq.CNOT(q1, q2)),
+    ]
+
+
 def test_prepend():
     q0, q1 = cirq.LineQubit.range(2)
     op_id0 = OpIdentifier(cirq.XPowGate, q0)

cirq/devices/noise_utils.py

@@ -13,7 +13,6 @@
 # limitations under the License.
 
 from typing import TYPE_CHECKING, Any, Dict, Tuple, Type, Union
-import warnings
 import numpy as np
 
 from cirq import ops, protocols, value
@@ -47,9 +46,6 @@ def qubits(self) -> Tuple['cirq.Qid', ...]:
     def _predicate(self, *args, **kwargs):
         return self._gate_family._predicate(*args, **kwargs)
 
-    def swapped(self):
-        return OpIdentifier(self.gate_type, *self.qubits[::-1])
-
     def is_proper_subtype_of(self, op_id: 'OpIdentifier'):
         """Returns true if this is contained within op_id, but not equal to it.
 
@@ -175,29 +171,6 @@ def pauli_error_from_t1(t_ns: float, t1_ns: float) -> float:
     return (1 - np.exp(-t_ns / t2)) / 2 + (1 - np.exp(-t_ns / t1_ns)) / 4
 
 
-def pauli_error_from_depolarization(t_ns: float, t1_ns: float, pauli_error: float = 0) -> float:
-    """Calculates the amount of pauli error from depolarization.
-
-    This computes non-T1 error for a specific duration, `t`. If pauli error
-    from T1 decay is more than total pauli error, this returns zero; otherwise,
-    it returns the portion of pauli error not attributable to T1 error.
-
-    Args:
-        t_ns: The duration of the gate in ns.
-        t1_ns: The T1 decay constant in ns.
-        pauli_error: The total pauli error.
-
-    Returns:
-        Calculated Pauli error resulting from depolarization.
-    """
-    t1_pauli_error = pauli_error_from_t1(t_ns, t1_ns)
-    if pauli_error >= t1_pauli_error:
-        return pauli_error - t1_pauli_error
-
-    warnings.warn("Pauli error from T1 decay is greater than total Pauli error", RuntimeWarning)
-    return 0
-
-
 def average_error(decay_constant: float, num_qubits: int = 1) -> float:
     """Calculates the average error from the depolarization decay constant.
 
@@ -213,7 +186,7 @@ def average_error(decay_constant: float, num_qubits: int = 1) -> float:
 
 
 def decoherence_pauli_error(t1_ns: float, tphi_ns: float, gate_time_ns: float) -> float:
-    """The component of Pauli error caused by decoherence.
+    """The component of Pauli error caused by decoherence on a single qubit.
 
     Args:
         t1_ns: T1 time in nanoseconds.

cirq/devices/noise_utils_test.py

@@ -23,7 +23,6 @@
     pauli_error_to_decay_constant,
     xeb_fidelity_to_decay_constant,
     pauli_error_from_t1,
-    pauli_error_from_depolarization,
     average_error,
     decoherence_pauli_error,
 )
@@ -58,7 +57,7 @@ def test_op_id_str():
 def test_op_id_swap():
     q0, q1 = cirq.LineQubit.range(2)
     base_id = OpIdentifier(cirq.CZPowGate, q0, q1)
-    swap_id = base_id.swapped()
+    swap_id = OpIdentifier(base_id.gate_type, *base_id.qubits[::-1])
     assert cirq.CZ(q0, q1) in base_id
     assert cirq.CZ(q0, q1) not in swap_id
     assert cirq.CZ(q1, q0) not in base_id
@@ -125,19 +124,6 @@ def test_pauli_error_from_t1(t, t1_ns, expected_output):
     assert val == expected_output
 
 
-@pytest.mark.parametrize(
-    't,t1_ns,pauli_error,expected_output',
-    [
-        (20, 1e5, 0.01, 0.01 - ((1 - np.exp(-20 / 2e5)) / 2 + (1 - np.exp(-20 / 1e5)) / 4)),
-        # In this case, the formula produces a negative result.
-        (4000, 1e4, 0.01, 0),
-    ],
-)
-def test_pauli_error_from_depolarization(t, t1_ns, pauli_error, expected_output):
-    val = pauli_error_from_depolarization(t, t1_ns, pauli_error)
-    assert val == expected_output
-
-
 @pytest.mark.parametrize(
     'decay_constant,num_qubits,expected_output',
     [