Add simulate method to Circuit

cirq/__init__.py

@@ -39,6 +39,7 @@
     PointOptimizationSummary,
     PointOptimizer,
     TextDiagramDrawer,
+    SimulateCircuitResult,
 )
 
 from cirq.devices import (

cirq/circuits/__init__.py

@@ -40,3 +40,6 @@
     PointOptimizer,
     PointOptimizationSummary,
 )
+from cirq.circuits.simulate_circuit_result import (
+    SimulateCircuitResult,
+)

cirq/circuits/circuit.py

@@ -27,12 +27,13 @@
 from cirq import ops
 from cirq.circuits.insert_strategy import InsertStrategy
 from cirq.circuits.moment import Moment
+from cirq.circuits.simulate_circuit_result import SimulateCircuitResult
 from cirq.circuits.text_diagram_drawer import TextDiagramDrawer
 from cirq.extension import Extensions
-from cirq.ops import QubitId
 
 if TYPE_CHECKING:
     # pylint: disable=unused-import
+    from cirq import study
     from typing import Set
 
 
@@ -425,14 +426,61 @@ def clear_operations_touching(self,
                 self.moments[k] = self.moments[k].without_operations_touching(
                     qubits)
 
-    def qubits(self) -> FrozenSet[QubitId]:
+    def qubits(self) -> FrozenSet[ops.QubitId]:
         """Returns the qubits acted upon by Operations in this circuit."""
         return frozenset(q for m in self.moments for q in m.qubits)
 
+    def simulate(self,
+                 *unnamed_args,
+                 extensions: Extensions = None,
+                 initial_state: Union[int, np.ndarray] = 0,
+                 param_resolver: 'study.ParamResolver' = None,
+                 qubit_order: ops.QubitOrderOrList = ops.QubitOrder.DEFAULT):
+        """Samples measurements and an output vector for this circuit.
+
+        Args:
+            param_resolver: Specifies how to turn parameterized gates into
+                concrete gates.
+            qubit_order: Determines the canonical ordering of the qubits, which
+                defines the order of amplitudes in the state vector (which is
+                relevant when specifying an initial state and interpreting the
+                output state).
+            initial_state: If an int, the state is set to the computational
+                basis state corresponding to this state.
+                Otherwise  if this is a np.ndarray it is the full initial
+                state. In this case it must be the correct size, be normalized
+                (an L2 norm of 1), and be safely castable to a np.complex64.
+            extensions: Extensions that will be applied while trying to
+                decompose the circuit's gates. If None, defaults to a set of
+                extensions capable of decomposing CompositeGate instances and
+                synthesizing KnownMatrixGate instances applied to 1 or 2
+                qubits.
+
+        Returns:
+            A SimulateCircuitResult storing the final state vector as well as
+            keyed measurement results.
+        """
+
+        assert not unnamed_args
+
+        # Local import to avoid circular dependency.
+        from cirq import google
+
+        result = google.XmonSimulator().simulate(
+            circuit=self,
+            param_resolver=param_resolver,
+            qubit_order=qubit_order,
+            initial_state=initial_state,
+            extensions=extensions)
+
+        return SimulateCircuitResult(
+            measurements=result.measurements,
+            final_state=result.final_state)
+
     def to_unitary_matrix(
             self,
             qubit_order: ops.QubitOrderOrList = ops.QubitOrder.DEFAULT,
-            qubits_that_should_be_present: Iterable[QubitId] = (),
+            qubits_that_should_be_present: Iterable[ops.QubitId] = (),
             ignore_terminal_measurements: bool = True,
             ext: Extensions = None) -> np.ndarray:
         """Converts the circuit into a unitary matrix, if possible.
@@ -464,7 +512,7 @@ def to_unitary_matrix(
         qs = ops.QubitOrder.as_qubit_order(qubit_order).order_for(
             self.qubits().union(qubits_that_should_be_present))
         qubit_map = {i: q
-                     for q, i in enumerate(qs)}  # type: Dict[QubitId, int]
+                     for q, i in enumerate(qs)}  # type: Dict[ops.QubitId, int]
         matrix_ops = _flatten_to_known_matrix_ops(self.iter_ops(), ext)
         return _operations_to_unitary_matrix(matrix_ops,
                                              qubit_map,
@@ -589,7 +637,7 @@ def _get_operation_text_diagram_exponent(op: ops.Operation,
 
 def _draw_moment_in_diagram(moment: Moment,
                             ext: Extensions,
-                            qubit_map: Dict[QubitId, int],
+                            qubit_map: Dict[ops.QubitId, int],
                             out_diagram: TextDiagramDrawer,
                             precision: Optional[int]):
     if not moment.operations:
@@ -655,11 +703,11 @@ def _flatten_to_known_matrix_ops(iter_ops: Iterable[ops.Operation],
 
 
 def _operations_to_unitary_matrix(iter_ops: Iterable[ops.Operation],
-                                  qubit_map: Dict[QubitId, int],
+                                  qubit_map: Dict[ops.QubitId, int],
                                   ignore_terminal_measurements: bool,
                                   ext: Extensions) -> np.ndarray:
     total = np.eye(1 << len(qubit_map))
-    measured_qubits = set()  # type: Set[QubitId]
+    measured_qubits = set()  # type: Set[ops.QubitId]
     for op in iter_ops:
         meas_gate = ext.try_cast(op.gate, ops.MeasurementGate)
         if meas_gate is not None:
@@ -688,7 +736,7 @@ def _operations_to_unitary_matrix(iter_ops: Iterable[ops.Operation],
 
 
 def _operation_to_unitary_matrix(op: ops.Operation,
-                                 qubit_map: Dict[QubitId, int],
+                                 qubit_map: Dict[ops.QubitId, int],
                                  ext: Extensions) -> np.ndarray:
     known_matrix_gate = ext.try_cast(op.gate, ops.KnownMatrixGate)
     if known_matrix_gate is None:

cirq/circuits/circuit_test.py

@@ -1078,3 +1078,61 @@ def test_iter_ops():
         ops.H(a)]
 
     assert list(c.iter_ops()) == expected
+
+
+def test_simulate():
+    a = ops.NamedQubit('a')
+    b = ops.NamedQubit('b')
+
+    c = Circuit.from_ops(ops.X(a), ops.CNOT(a, b))
+
+    # Classical operations.
+    assert c.simulate().approx_eq(
+        cirq.SimulateCircuitResult(
+            measurements={},
+            final_state=np.array([0, 0, 0, 1])),
+        atol=1e-6)
+
+    # Custom initial computational basis state.
+    assert c.simulate(initial_state=1).approx_eq(
+        cirq.SimulateCircuitResult(
+            measurements={},
+            final_state=np.array([0, 0, 1, 0])),
+        atol=1e-6)
+
+    # Custom initial vector.
+    assert c.simulate(initial_state=np.array([0, 0.8, 0, -0.6],
+                                             dtype=np.float32)).approx_eq(
+            cirq.SimulateCircuitResult(
+                measurements={},
+                final_state=np.array([0, -0.6, 0.8, 0])),
+            atol=1e-6)
+
+    # Empty Circuit
+    assert Circuit().simulate().approx_eq(
+        cirq.SimulateCircuitResult({}, np.array([1])))
+
+    # Superposing operation.
+    assert Circuit.from_ops(ops.H(a)).simulate().approx_eq(
+        cirq.SimulateCircuitResult(
+            measurements={},
+            final_state=np.array([np.sqrt(0.5), np.sqrt(0.5)])),
+        atol=1e-6)
+
+    # Collapsing measurement samples.
+    c2 = Circuit.from_ops(ops.H(a),
+                          ops.H(b),
+                          ops.CZ(a, b),
+                          ops.MeasurementGate('a').on(a))
+    for _ in range(10):
+        r2 = c2.simulate()
+        if r2.measurements['a'][0]:
+            assert r2.approx_eq(cirq.SimulateCircuitResult(
+                    measurements={'a': np.array([True])},
+                    final_state=np.array([0, 0, np.sqrt(0.5), -np.sqrt(0.5)])),
+                atol=1e-6)
+        else:
+            assert r2.approx_eq(cirq.SimulateCircuitResult(
+                    measurements={'a': np.array([False])},
+                    final_state=np.array([np.sqrt(0.5), np.sqrt(0.5), 0, 0])),
+                atol=1e-6)

cirq/circuits/simulate_circuit_result.py

@@ -0,0 +1,72 @@
+from typing import Dict, Any, Iterable
+
+import numpy as np
+from cirq import linalg
+
+
+class SimulateCircuitResult:
+    """Measurements and final state vector from a circuit simulation.
+
+    Attributes:
+        measurements: A dictionary from measurement gate key to measurement
+            results. Each measurement result value is a numpy ndarray of actual
+            boolean measurement results (ordered by the qubits acted on by the
+            measurement gate.)
+        final_state: The state vector output by the circuit. The final state
+            of the simulated system.
+    """
+
+    def __init__(self,
+                 measurements: Dict[str, np.ndarray],
+                 final_state: np.ndarray) -> None:
+        """
+        Args:
+            measurements: A dictionary from measurement gate key to measurement
+                results. Each measurement result value is a numpy ndarray of
+                actual boolean measurement results (ordered by the qubits acted
+                on by the measurement gate.)
+            final_state: The state vector output by the circuit. The final
+                state of the simulated system.
+        """
+        self.measurements = measurements
+        self.final_state = final_state
+
+    def __repr__(self):
+        return ('SimulateCircuitResult(measurements={!r}, '
+                'final_state={!r})').format(self.measurements,
+                                            self.final_state)
+
+    def __str__(self):
+        return 'measurements: {}\nfinal_state: {}'.format(
+            _keyed_iterable_bitstrings(self.measurements),
+            self.final_state)
+
+    def approx_eq(self,
+                  other: 'SimulateCircuitResult',
+                  atol: float=0,
+                  ignore_global_phase=True) -> bool:
+        if len(self.measurements) != len(other.measurements):
+            return False
+
+        for k, v in self.measurements.items():
+            other_v = other.measurements.get(k)
+            if (other_v is None or
+                    len(other_v) != len(v) or
+                    np.any(v != other_v)):
+                return False
+
+        cmp_vector = (linalg.allclose_up_to_global_phase if ignore_global_phase
+                      else np.allclose)
+        return cmp_vector(self.final_state, other.final_state, atol=atol)
+
+
+def _iterable_bitstring(vals: Iterable[Any]) -> str:
+    return ''.join('1' if v else '0' for v in vals)
+
+
+def _keyed_iterable_bitstrings(keyed_vals: Dict[str, Iterable[Any]]) -> str:
+    if not keyed_vals:
+        return '(none)'
+    results = sorted(
+        (key, _iterable_bitstring(val)) for key, val in keyed_vals.items())
+    return ' '.join(['{}={}'.format(key, val) for key, val in results])

cirq/circuits/simulate_circuit_result_test.py

@@ -0,0 +1,101 @@
+# Copyright 2018 The Cirq Developers
+#
+# Licensed under the Apache License, Version 2.0 (the "License");
+# you may not use this file except in compliance with the License.
+# You may obtain a copy of the License at
+#
+#     https://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+import itertools
+import numpy as np
+
+import cirq
+
+
+def test_init():
+    measurements = {'a': np.array([True])}
+    final_state = np.array([0, 1])
+    c = cirq.SimulateCircuitResult(measurements=measurements,
+                                   final_state=final_state)
+    assert c.measurements is measurements
+    assert c.final_state is final_state
+
+
+def test_str():
+    empty = cirq.SimulateCircuitResult({}, np.array([1]))
+    assert str(empty) == 'measurements: (none)\nfinal_state: [1]'
+
+    multi = cirq.SimulateCircuitResult(
+        {'a': np.array([True, True]), 'b': np.array([False])},
+        np.array([0, 1, 0, 0]))
+    assert str(multi) == 'measurements: a=11 b=0\nfinal_state: [0 1 0 0]'
+
+
+@cirq.testing.only_test_in_python3
+def test_repr():
+    multi = cirq.SimulateCircuitResult(
+        {},
+        np.array([0, 1, 0, 0]))
+    expected = ("SimulateCircuitResult(measurements={}, "
+                "final_state=array([0, 1, 0, 0]))")
+    assert repr(multi) == expected
+
+
+def test_approx_eq():
+    empty = cirq.SimulateCircuitResult({}, np.array([1]))
+
+    assert empty.approx_eq(empty, atol=1e-2)
+    assert empty.approx_eq(empty, atol=1e-6)
+    assert empty.approx_eq(empty, atol=1e-2, ignore_global_phase=False)
+    assert empty.approx_eq(empty, atol=1e-6, ignore_global_phase=False)
+
+    empty_neg = cirq.SimulateCircuitResult({}, np.array([-1]))
+    assert empty.approx_eq(empty_neg, atol=1e-2)
+    assert empty.approx_eq(empty_neg, atol=1e-6)
+    assert not empty.approx_eq(empty_neg, atol=1e-2, ignore_global_phase=False)
+    assert not empty.approx_eq(empty_neg, atol=1e-6, ignore_global_phase=False)
+
+    empty_near = cirq.SimulateCircuitResult({}, np.array([0.999]))
+    assert empty.approx_eq(empty_near, atol=1e-2)
+    assert not empty.approx_eq(empty_near, atol=1e-6)
+    assert empty.approx_eq(empty_near, atol=1e-2, ignore_global_phase=False)
+    assert not empty.approx_eq(empty_near,
+                               atol=1e-6,
+                               ignore_global_phase=False)
+
+    just_on = cirq.SimulateCircuitResult({},
+                                         np.array([0, 1]))
+    near_on = cirq.SimulateCircuitResult({},
+                                         np.array([0, 0.999]))
+    just_off = cirq.SimulateCircuitResult({},
+                                         np.array([1, 0]))
+    collapse_on = cirq.SimulateCircuitResult({'a': np.array([True])},
+                                             np.array([0, 1]))
+    collapse_off = cirq.SimulateCircuitResult({'a': np.array([False])},
+                                              np.array([1, 0]))
+    other_collapse_on = cirq.SimulateCircuitResult({'b': np.array([True])},
+                                                   np.array([0, 1]))
+    multi_collapse_on = cirq.SimulateCircuitResult(
+        {'a': np.array([True, True])},
+        np.array([0, 1]))
+
+    approx_eq_groups = [
+        [just_on, near_on],
+        [just_off],
+        [collapse_on],
+        [collapse_off],
+        [other_collapse_on],
+        [multi_collapse_on],
+    ]
+    for g1, g2 in itertools.product(approx_eq_groups, repeat=2):
+        for e1, e2 in itertools.product(g1, g2):
+            for g in [False, True]:
+                assert (g1 is g2) == e1.approx_eq(e2,
+                                                  atol=1e-2,
+                                                  ignore_global_phase=g)