From d4f7fb4c7e82bd748a1d1ae78415262210306ae6 Mon Sep 17 00:00:00 2001 From: Bram Evert Date: Fri, 11 Sep 2026 09:38:05 +0000 Subject: [PATCH 1/3] Support affine gate-parameter expressions A gate angle may now be any affine expression a*theta + b in a single memory reference (RX(theta[0]/2 + pi), RZ(-2*phi[1]), RX(theta + theta)), which is what quilc emits when it compiles a parametric program. The scale and offset travel with the ParametricGate and are per-member data of a gate batch, so the batching key is unchanged and such gates still share one vmap with plain references. Products, powers and functions of references are rejected with an error naming the parameter. Co-Authored-By: Claude Fable 5.1 Claude-Session: https://claude.ai/code/session_01JfZKmgqinMhR4BS87y4G4F --- docs/source/simulation_architecture.rst | 7 +- pyquil/simulation/_resolver.py | 111 +++++++++++++++---- pyquil/simulation/_simulator.py | 14 ++- test/unit/test_affine_parameters.py | 136 ++++++++++++++++++++++++ test/unit/test_density_matrix.py | 7 +- 5 files changed, 247 insertions(+), 28 deletions(-) create mode 100644 test/unit/test_affine_parameters.py diff --git a/docs/source/simulation_architecture.rst b/docs/source/simulation_architecture.rst index 47390df87..4b9b2bdee 100644 --- a/docs/source/simulation_architecture.rst +++ b/docs/source/simulation_architecture.rst @@ -216,7 +216,12 @@ Expansion does several things at once: *not* resolved to a number. It is wrapped in a ``ParametricGate`` that, given :math:`\theta`, constructs the gate matrix. This keeps gate construction inside the traced/differentiated graph, which is what makes ``jax.grad`` with - respect to gate angles work. + respect to gate angles work. An angle may be any *affine* expression + :math:`a\,\theta_i + b` in a single memory reference -- ``RX(theta[0]/2 + pi)``, + ``RZ(-2*phi[1])`` -- which is what quilc emits when it compiles a parametric program; + the scale and offset travel with the gate, so such gates still batch with plain + ``RX(theta[1])``. Products or functions of references (``theta[0]*theta[1]``, + ``SIN(theta[0])``) are rejected with an error naming the parameter. * **DEFCIRCUIT and cycle expansion.** ``DEFCIRCUIT`` bodies are expanded with formal-argument substitution. When a circuit invocation matches a diff --git a/pyquil/simulation/_resolver.py b/pyquil/simulation/_resolver.py index afe538d24..145cd657a 100644 --- a/pyquil/simulation/_resolver.py +++ b/pyquil/simulation/_resolver.py @@ -54,7 +54,7 @@ NoiseModelLike, ) from pyquil.quil import Program -from pyquil.quilatom import MemoryReference, Qubit, _contained_mrefs, substitute +from pyquil.quilatom import Add, Div, MemoryReference, Mul, Qubit, Sub, _contained_mrefs, substitute from pyquil.quilbase import ( AbstractInstruction, ArithmeticBinaryOp, @@ -85,6 +85,50 @@ ParameterRef: TypeAlias = tuple[str, int] +#: ``(reference, scale, offset)``: the value ``scale * reference + offset``, or a constant when +#: ``reference`` is ``None``. +_AffineForm: TypeAlias = tuple[MemoryReference | None, float, float] + + +def _affine_form(expression: Any) -> _AffineForm | None: + """Write a gate parameter as ``scale * theta + offset`` in one memory reference. + + Returns ``None`` when the expression is not of that form: a product or quotient of two + references, a power, a function such as ``SIN``, or two *different* references. + + :param expression: A number, a :class:`~pyquil.quilatom.MemoryReference`, or an + arithmetic expression over them. + """ + if isinstance(expression, MemoryReference): + return expression, 1.0, 0.0 + if not _contained_mrefs(expression): + return None, 0.0, float(np.real(expression)) + if isinstance(expression, (Add, Sub)): + left, right = _affine_form(expression.op1), _affine_form(expression.op2) + if left is None or right is None: + return None + (ref_l, a_l, b_l), (ref_r, a_r, b_r) = left, right + if ref_l is not None and ref_r is not None and ref_l != ref_r: + return None + sign = 1.0 if isinstance(expression, Add) else -1.0 + return ref_l if ref_l is not None else ref_r, a_l + sign * a_r, b_l + sign * b_r + if isinstance(expression, (Mul, Div)): + left, right = _affine_form(expression.op1), _affine_form(expression.op2) + if left is None or right is None: + return None + (ref_l, a_l, b_l), (ref_r, a_r, b_r) = left, right + if isinstance(expression, Div): + if ref_r is not None or b_r == 0: + return None + return ref_l, a_l / b_r, b_l / b_r + if ref_l is not None and ref_r is not None: + return None + if ref_l is None: + return ref_r, b_l * a_r, b_l * b_r + return ref_l, a_l * b_r, b_l * b_r + return None + + @dataclass(frozen=True, slots=True) class ParametricGate: """A parametric gate whose matrix depends on runtime parameters. @@ -93,22 +137,36 @@ class ParametricGate: constructor and parameter layout are exposed so that gates of the same kind can be built together in one vectorised operation. + A gate argument is either a literal number or an affine function ``scale * theta + offset`` + of one slot of the parameter vector -- which covers a bare memory reference (``RX(theta[0])``) + as well as the arithmetic quilc emits when it compiles parametric programs + (``RX(theta[0]/2 + pi)``). + :param gate_fn: The quax gate constructor (e.g. ``qx.gates.RX``), or a parametric ``DEFGATE`` callable. :param param_indices: For each gate argument, its slot in the flat parameter vector, or ``-1`` when the argument is a literal number. Gates that read the same memory reference share a slot; see :func:`expand_program`. :param concrete_values: For each gate argument, its literal value (``nan`` for a slot). + :param scales: For each gate argument, the factor multiplying the slot value (``1`` for a + bare reference; unused for a literal). + :param offsets: For each gate argument, the constant added to the scaled slot value (``0`` + for a bare reference; unused for a literal). """ gate_fn: Callable[..., qx.Operator] param_indices: tuple[int, ...] concrete_values: tuple[float, ...] + scales: tuple[float, ...] + offsets: tuple[float, ...] def __call__(self, params: Array) -> qx.Unitary: """Build the gate for one parameter vector.""" resolved: list[Any] = [ - params[pi] if pi >= 0 else cv for pi, cv in zip(self.param_indices, self.concrete_values, strict=True) + params[pi] * scale + offset if pi >= 0 else cv + for pi, cv, scale, offset in zip( + self.param_indices, self.concrete_values, self.scales, self.offsets, strict=True + ) ] result = self.gate_fn(*resolved) if not isinstance(result, qx.Unitary): @@ -337,37 +395,46 @@ def _resolve_gate(inst: Gate) -> tuple[ExpandedOp, tuple[int, ...]]: param_indices: list[int] = [] concrete_values: list[float] = [] + scales: list[float] = [] + offsets: list[float] = [] for p in inst.params: - mrefs = _contained_mrefs(p) # type: ignore[arg-type] - if not mrefs: - # A concrete number: a compile-time constant for this gate. - param_indices.append(-1) - concrete_values.append(float(np.real(p))) - elif not isinstance(p, MemoryReference): - # An arithmetic expression over one or more memory regions, e.g. - # ``RX(theta[0] / 2) 0``. Each ParametricGate argument maps to a single - # slot of the flat parameter vector, which is what lets the simulator - # batch same-shaped gates under one ``jax.vmap``; an arbitrary - # expression would have to become part of that batching key. + form = _affine_form(p) + if form is None: + # Each ParametricGate argument is an affine function of a single slot of the + # parameter vector; that is what lets the simulator batch same-shaped gates + # under one ``jax.vmap`` with per-gate scale and offset arrays. Anything + # else -- a product of two references, SIN(theta), theta^2 -- would need + # its own traced graph. + mrefs = _contained_mrefs(p) # type: ignore[arg-type] raise ValueError( - f"Gate parameter {p} in {inst.out()!r} is an expression over memory " - f"region(s) {sorted(m.name for m in mrefs)}, which is not supported. " - "Pass the parameter directly (e.g. RX(theta[0]) with the division folded " - "into the value you bind), or substitute concrete values into the program " - "before simulating." + f"Gate parameter {p} in {inst.out()!r} is not an affine expression " + f"(a * theta + b) in a single memory reference; it involves " + f"{sorted(str(m) for m in mrefs)}. Only such expressions are supported: " + "rewrite the program, or substitute concrete values before simulating." ) - elif p.name in measure_regs: + ref, scale, offset = form + if ref is None: + # A concrete number: a compile-time constant for this gate. + param_indices.append(-1) + concrete_values.append(offset) + scales.append(1.0) + offsets.append(0.0) + elif ref.name in measure_regs: # Classically-conditioned angle: the value is only known mid-circuit. raise ValueError( f"Gate parameter {p} in {inst.out()!r} reads memory region " - f"'{p.name}', which is written by a MEASURE in this program. " + f"'{ref.name}', which is written by a MEASURE in this program. " "Feed-forward (classically-conditioned) parameters are not supported." ) else: - param_indices.append(slots.setdefault((p.name, p.offset), len(slots))) + param_indices.append(slots.setdefault((ref.name, ref.offset), len(slots))) concrete_values.append(float("nan")) + scales.append(scale) + offsets.append(offset) - return ParametricGate(gate_def, tuple(param_indices), tuple(concrete_values)), qubits + return ParametricGate( + gate_def, tuple(param_indices), tuple(concrete_values), tuple(scales), tuple(offsets) + ), qubits # Fixed gate → resolve to Unitary now. unitary = get_instruction_unitary(inst, custom_gates=custom_gates) diff --git a/pyquil/simulation/_simulator.py b/pyquil/simulation/_simulator.py index 00b4bbe67..4db774398 100644 --- a/pyquil/simulation/_simulator.py +++ b/pyquil/simulation/_simulator.py @@ -626,6 +626,11 @@ class _GateBatch: positions: list[int] = field(default_factory=list) #: Parameter-vector index for each free argument, one list per member. param_indices: list[list[int]] = field(default_factory=list) + #: Affine coefficients of each free argument (``scale * params[index] + offset``), one list + #: per member. They are per-member data rather than part of the batch key, so + #: ``RX(theta[0] / 2)`` and ``RX(theta[1])`` share one vmap. + param_scales: list[list[float]] = field(default_factory=list) + param_offsets: list[list[float]] = field(default_factory=list) def builder(self) -> Callable[[Array], Array]: """Return ``params -> (n_members, width, width)`` embedded gate matrices.""" @@ -635,6 +640,8 @@ def builder(self) -> Callable[[Array], Array]: target_dims, group_positions = self.target_dims, self.group_positions width, as_superop = self.width, self.as_superop param_indices = jnp.asarray(self.param_indices) # (n_members, n_free) + scales = jnp.asarray(self.param_scales, dtype=float) + offsets = jnp.asarray(self.param_offsets, dtype=float) def single(free_values: Array) -> Array: args: list[Any] = [None] * n_args @@ -648,7 +655,7 @@ def single(free_values: Array) -> Array: return _embed_op_to_group(gate, target_dims, group_positions, width, as_superop=as_superop) batched = jax.vmap(single) - return lambda params: batched(params[param_indices]) + return lambda params: batched(params[param_indices] * scales + offsets) def _make_group_fold(group_start: list[int], n_ops: int, width: int) -> Callable[[Array], Array]: @@ -763,7 +770,10 @@ def _build_vectorized_operator_constructor( ) batches[key] = batch batch.positions.append(pos) - batch.param_indices.append([pi for pi in op.param_indices if pi >= 0]) + free = [j for j, pi in enumerate(op.param_indices) if pi >= 0] + batch.param_indices.append([op.param_indices[j] for j in free]) + batch.param_scales.append([op.scales[j] for j in free]) + batch.param_offsets.append([op.offsets[j] for j in free]) else: # Constant operations are embedded once, eagerly. In superoperator mode this # also covers the non-unitary ops a noise model contributes (channel SuperOps, diff --git a/test/unit/test_affine_parameters.py b/test/unit/test_affine_parameters.py new file mode 100644 index 000000000..259cfbd9d --- /dev/null +++ b/test/unit/test_affine_parameters.py @@ -0,0 +1,136 @@ +############################################################################## +# Copyright 2026 Rigetti Computing +# +# 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 +# +# http://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. +############################################################################## +"""Affine gate-parameter expressions (``a * theta + b``) in the quax-based simulators.""" + +import jax +import jax.numpy as jnp +import numpy as np +import pytest + +from pyquil.gates import RX, RY, RZ, H +from pyquil.quil import Program +from pyquil.quilatom import MemoryReference, quil_sin +from pyquil.quilbase import Declare +from pyquil.simulation._resolver import _affine_form, expand_program +from pyquil.simulation._simulator import DensityMatrixSimulator, PureStateVectorSimulator + +THETA = MemoryReference("theta", 0) +PHI = MemoryReference("theta", 1) + + +def _program(*gates): + program = Program() + program += Declare("theta", "REAL", 2) + for gate in gates: + program += gate + return program + + +def _state(program, **memory): + sim = PureStateVectorSimulator(program) + return np.asarray(sim.compute(sim.linearize(memory)).matrix).reshape(-1) + + +class TestAffineForm: + @pytest.mark.parametrize( + ("expression", "expected"), + [ + (THETA, (THETA, 1.0, 0.0)), + (THETA / 2 + np.pi, (THETA, 0.5, np.pi)), + (2 * PHI, (PHI, 2.0, 0.0)), + (-THETA, (THETA, -1.0, 0.0)), + (np.pi - PHI, (PHI, -1.0, np.pi)), + (THETA + THETA, (THETA, 2.0, 0.0)), + ((3 * THETA - 1) / 4, (THETA, 0.75, -0.25)), + (0.7, (None, 0.0, 0.7)), + ], + ) + def test_recognises_affine_expressions(self, expression, expected): + ref, scale, offset = _affine_form(expression) + assert ref == expected[0] + assert scale == pytest.approx(expected[1]) + assert offset == pytest.approx(expected[2]) + + @pytest.mark.parametrize( + "expression", + [THETA * PHI, THETA * THETA, quil_sin(THETA), THETA**2, 2 / THETA, THETA + PHI], + ids=["product", "square", "sin", "power", "reciprocal", "two_references"], + ) + def test_rejects_non_affine_expressions(self, expression): + assert _affine_form(expression) is None + + +class TestSimulation: + @pytest.mark.parametrize( + ("gate", "literal"), + [ + (RX(THETA / 2 + np.pi, 0), lambda t, p: RX(t / 2 + np.pi, 0)), + (RZ(2 * PHI, 0), lambda t, p: RZ(2 * p, 0)), + (RX(-THETA, 0), lambda t, p: RX(-t, 0)), + (RY(np.pi - PHI, 0), lambda t, p: RY(np.pi - p, 0)), + (RX(THETA + THETA, 0), lambda t, p: RX(2 * t, 0)), + ], + ids=["half_plus_pi", "double", "negated", "pi_minus", "same_reference_twice"], + ) + def test_matches_the_literal_gate(self, gate, literal): + theta, phi = 0.37, 1.21 + got = _state(_program(H(0), gate), theta=[theta, phi]) + expected = _state(Program(H(0), literal(theta, phi))) + np.testing.assert_allclose(got, expected, atol=1e-12) + + def test_parsed_quilc_style_program(self): + """The arithmetic quilc emits for a compiled parametric program simulates directly.""" + program = Program("DECLARE theta REAL[1]\nRX((theta[0]/2)+pi) 0\nRZ(-2*theta[0]) 0") + theta = 0.6 + got = _state(program, theta=[theta]) + expected = _state(Program(RX(theta / 2 + np.pi, 0), RZ(-2 * theta, 0))) + np.testing.assert_allclose(got, expected, atol=1e-12) + + def test_shares_a_slot_and_batches_with_a_plain_reference(self): + """``RX(theta[0]/2)`` and ``RX(theta[1])`` use one slot each and one gate batch.""" + program = _program(RX(THETA / 2, 0), RX(PHI, 1), RX(THETA, 2)) + sim = PureStateVectorSimulator(program) + assert sim.parameters == (("theta", 0), ("theta", 1)) + theta, phi = 0.5, 0.9 + got = _state(program, theta=[theta, phi]) + expected = _state(Program(RX(theta / 2, 0), RX(phi, 1), RX(theta, 2))) + np.testing.assert_allclose(got, expected, atol=1e-12) + + def test_grad_and_jit(self): + """d/dtheta P(|1>) for RX(theta/2 + pi) is -sin(theta/2)/4 (the pi flips the population).""" + sim = DensityMatrixSimulator(_program(RX(THETA / 2 + np.pi, 0))) + + def excited(params): + return jnp.real(sim.compute(params).matrix[1, 1]) + + theta = 0.8 + params = jnp.array([theta]) + np.testing.assert_allclose(jax.jit(excited)(params), np.cos(theta / 4) ** 2, atol=1e-12) + grad = jax.grad(excited)(params) + np.testing.assert_allclose(grad, [-np.sin(theta / 2) / 4], atol=1e-12) + + def test_expand_program_records_scale_and_offset(self): + ops, _, parameters = expand_program(_program(RX(3 * THETA - 1, 0))) + (gate,) = ops + assert parameters == (("theta", 0),) + assert gate.param_indices == (0,) + assert gate.scales == (3.0,) + assert gate.offsets == (-1.0,) + + @pytest.mark.parametrize("expression", [THETA * PHI, quil_sin(THETA), THETA**2], ids=["product", "sin", "power"]) + def test_non_affine_parameter_reports_clearly(self, expression): + with pytest.raises(ValueError, match="not an affine expression"): + PureStateVectorSimulator(_program(RX(expression, 0))) diff --git a/test/unit/test_density_matrix.py b/test/unit/test_density_matrix.py index 9c63c446a..d1771c522 100644 --- a/test/unit/test_density_matrix.py +++ b/test/unit/test_density_matrix.py @@ -684,9 +684,10 @@ def test_unsupported_modifier_reports_clearly(self, program, modifier): with pytest.raises(ValueError, match=f"modifiers are not supported.*{modifier}"): _dm(program) - def test_expression_valued_parameter_reports_clearly(self): - program = Program(Declare("theta", "REAL", 1), RX(MemoryReference("theta", 0) / 2, 0)) - with pytest.raises(ValueError, match="expression over memory"): + def test_non_affine_parameter_reports_clearly(self): + theta = MemoryReference("theta", 0) + program = Program(Declare("theta", "REAL", 1), RX(theta * theta, 0)) + with pytest.raises(ValueError, match="not an affine expression"): DensityMatrixSimulator(program, qubits=[0]) def test_feed_forward_parameter_reports_clearly(self): From 13826d5d6c1e8ac757776c9a5d8f84fad255348f Mon Sep 17 00:00:00 2001 From: Bram Evert Date: Fri, 11 Sep 2026 10:07:36 +0000 Subject: [PATCH 2/3] Support every Quil arithmetic expression as a gate argument Replace the affine-only parser with a compiler from pyQuil expression trees to JAX functions: + - * / ^, SIN, COS, SQRT, EXP, CIS, and real or complex literals. A ParameterExpression carries the slots it reads and a shape key (references replaced by %0, %1, ...); gate batches are keyed by shape, so SIN(theta[0]) and SIN(theta[1]) still share one vmap and differ only in the slots gathered per member. Complex-valued arguments (CIS, a complex literal) are accepted for DEFGATE gates and rejected with a clear error for the built-in gates, whose angles are real. Parameter-free expressions such as SIN(pi/4) fold to literals. Co-Authored-By: Claude Fable 5.1 Claude-Session: https://claude.ai/code/session_01JfZKmgqinMhR4BS87y4G4F --- docs/source/simulation_architecture.rst | 13 +- pyquil/simulation/_resolver.py | 233 ++++++++++++++---------- pyquil/simulation/_simulator.py | 72 ++++---- test/unit/test_affine_parameters.py | 136 -------------- test/unit/test_density_matrix.py | 8 +- test/unit/test_parameter_expressions.py | 233 ++++++++++++++++++++++++ test/unit/test_resolver.py | 2 +- 7 files changed, 421 insertions(+), 276 deletions(-) delete mode 100644 test/unit/test_affine_parameters.py create mode 100644 test/unit/test_parameter_expressions.py diff --git a/docs/source/simulation_architecture.rst b/docs/source/simulation_architecture.rst index 4b9b2bdee..5b6d933ea 100644 --- a/docs/source/simulation_architecture.rst +++ b/docs/source/simulation_architecture.rst @@ -216,12 +216,13 @@ Expansion does several things at once: *not* resolved to a number. It is wrapped in a ``ParametricGate`` that, given :math:`\theta`, constructs the gate matrix. This keeps gate construction inside the traced/differentiated graph, which is what makes ``jax.grad`` with - respect to gate angles work. An angle may be any *affine* expression - :math:`a\,\theta_i + b` in a single memory reference -- ``RX(theta[0]/2 + pi)``, - ``RZ(-2*phi[1])`` -- which is what quilc emits when it compiles a parametric program; - the scale and offset travel with the gate, so such gates still batch with plain - ``RX(theta[1])``. Products or functions of references (``theta[0]*theta[1]``, - ``SIN(theta[0])``) are rejected with an error naming the parameter. + respect to gate angles work. An angle may be any Quil arithmetic expression over memory + references -- ``+ - * / ^``, the functions ``SIN``, ``COS``, ``SQRT``, ``EXP`` and ``CIS``, + and real or complex literals -- so ``RX(theta[0]/2 + pi)``, the form quilc emits when it + compiles a parametric program, and ``RZ(2*SIN(phi[1]))`` both simulate directly. Gates whose + expressions have the same *shape* (``SIN(theta[0])`` and ``SIN(theta[1])``) are still built in + one vectorised batch. A complex-valued argument (``CIS``, a complex literal) is accepted only + by a ``DEFGATE`` gate, since the built-in gates take real angles. * **DEFCIRCUIT and cycle expansion.** ``DEFCIRCUIT`` bodies are expanded with formal-argument substitution. When a circuit invocation matches a diff --git a/pyquil/simulation/_resolver.py b/pyquil/simulation/_resolver.py index 145cd657a..ebd6a241e 100644 --- a/pyquil/simulation/_resolver.py +++ b/pyquil/simulation/_resolver.py @@ -54,7 +54,20 @@ NoiseModelLike, ) from pyquil.quil import Program -from pyquil.quilatom import Add, Div, MemoryReference, Mul, Qubit, Sub, _contained_mrefs, substitute +from pyquil.quilatom import ( + Add, + BinaryExp, + Div, + Function, + MemoryReference, + Mul, + Parameter, + Pow, + Qubit, + Sub, + _contained_mrefs, + substitute, +) from pyquil.quilbase import ( AbstractInstruction, ArithmeticBinaryOp, @@ -85,48 +98,104 @@ ParameterRef: TypeAlias = tuple[str, int] -#: ``(reference, scale, offset)``: the value ``scale * reference + offset``, or a constant when -#: ``reference`` is ``None``. -_AffineForm: TypeAlias = tuple[MemoryReference | None, float, float] +#: The Quil arithmetic functions, as JAX functions. ``CIS(x)`` is ``exp(i x)``. +_FUNCTIONS: dict[str, Callable[[Array], Array]] = { + "SIN": jnp.sin, + "COS": jnp.cos, + "SQRT": jnp.sqrt, + "EXP": jnp.exp, + "CIS": lambda x: jnp.exp(1j * x), +} +_BINARY_OPERATORS: dict[type[BinaryExp], Callable[[Array, Array], Array]] = { + Add: jnp.add, + Sub: jnp.subtract, + Mul: jnp.multiply, + Div: jnp.divide, + Pow: jnp.power, +} + + +@dataclass(frozen=True, slots=True) +class ParameterExpression: + """A gate argument given as a Quil arithmetic expression over memory references. + + Any expression Quil allows is supported: ``+ - * / ^``, the functions ``SIN``, ``COS``, + ``SQRT``, ``EXP`` and ``CIS``, and real or complex literals -- for example + ``RX(theta[0]/2 + pi)``, ``RZ(2*SIN(phi[1]))`` or ``CPHASE(CIS(theta[0])) 0 1`` for a + ``DEFGATE`` taking a complex parameter. Calling an instance with the flat parameter + vector evaluates it with JAX, so it can be jitted and differentiated through. + + :param slot_indices: Slots of the parameter vector the expression reads, in order of first + appearance. + :param key: The expression with its references replaced by ``%0``, ``%1``, ... in that + order. Two expressions with the same key have the same shape and constants and differ + only in which slots they read, which is what lets the simulator evaluate them together + in one vectorised operation. + :param is_complex: Whether the value may be complex (the expression contains ``CIS`` or a + complex literal). Everything else is evaluated in real arithmetic; note that Quil + would evaluate ``SQRT`` of a negative number or a fractional power of one as complex, + which real arithmetic reports as ``nan``. + """ + + slot_indices: tuple[int, ...] + key: str + is_complex: bool + _fn: Callable[[Array], Array] + + def evaluate(self, values: Array) -> Array: + """Evaluate the expression given the values of its slots, in :attr:`slot_indices` order.""" + return self._fn(values) + + def __call__(self, params: Array) -> Array: + """Evaluate the expression for one parameter vector.""" + return self.evaluate(params[jnp.asarray(self.slot_indices)]) -def _affine_form(expression: Any) -> _AffineForm | None: - """Write a gate parameter as ``scale * theta + offset`` in one memory reference. +def _literal_value(expression: Any) -> float | complex: + """Evaluate a parameter-free expression (a number, ``pi/2``, ``SIN(pi/4)``, ``1i``) to a scalar.""" + value = complex(np.asarray(expression, dtype=complex).item()) if not isinstance(expression, complex) else expression + return value.real if value.imag == 0 else value - Returns ``None`` when the expression is not of that form: a product or quotient of two - references, a power, a function such as ``SIN``, or two *different* references. - :param expression: A number, a :class:`~pyquil.quilatom.MemoryReference`, or an - arithmetic expression over them. +def _compile_expression(expression: Any, slot_of: Callable[[MemoryReference], int]) -> ParameterExpression: + """Compile a Quil expression over memory references into a :class:`ParameterExpression`. + + :param expression: The gate parameter; must contain at least one memory reference. + :param slot_of: Maps a memory reference to its slot in the parameter vector. + :raises ValueError: If the expression contains an unbound ``DEFGATE`` parameter. """ - if isinstance(expression, MemoryReference): - return expression, 1.0, 0.0 - if not _contained_mrefs(expression): - return None, 0.0, float(np.real(expression)) - if isinstance(expression, (Add, Sub)): - left, right = _affine_form(expression.op1), _affine_form(expression.op2) - if left is None or right is None: - return None - (ref_l, a_l, b_l), (ref_r, a_r, b_r) = left, right - if ref_l is not None and ref_r is not None and ref_l != ref_r: - return None - sign = 1.0 if isinstance(expression, Add) else -1.0 - return ref_l if ref_l is not None else ref_r, a_l + sign * a_r, b_l + sign * b_r - if isinstance(expression, (Mul, Div)): - left, right = _affine_form(expression.op1), _affine_form(expression.op2) - if left is None or right is None: - return None - (ref_l, a_l, b_l), (ref_r, a_r, b_r) = left, right - if isinstance(expression, Div): - if ref_r is not None or b_r == 0: - return None - return ref_l, a_l / b_r, b_l / b_r - if ref_l is not None and ref_r is not None: - return None - if ref_l is None: - return ref_r, b_l * a_r, b_l * b_r - return ref_l, a_l * b_r, b_l * b_r - return None + references: list[MemoryReference] = [] + is_complex = False + + def build(node: Any) -> tuple[Callable[[Array], Array], str]: + nonlocal is_complex + if isinstance(node, MemoryReference): + if node not in references: + references.append(node) + index = references.index(node) + return (lambda values, index=index: values[index]), f"%{index}" + if isinstance(node, Parameter): + raise ValueError(f"Unbound DEFGATE parameter {node} in gate argument {expression}.") + if isinstance(node, BinaryExp): + left, left_key = build(node.op1) + right, right_key = build(node.op2) + operator = _BINARY_OPERATORS[type(node)] + return ( + lambda values: operator(left(values), right(values)) + ), f"({left_key}{node.operator.strip()}{right_key})" + if isinstance(node, Function): + inner, inner_key = build(node.expression) + if node.name not in _FUNCTIONS: + raise ValueError(f"Unknown Quil function {node.name!r} in gate argument {expression}.") + function = _FUNCTIONS[node.name] + is_complex = is_complex or node.name == "CIS" + return (lambda values: function(inner(values))), f"{node.name}({inner_key})" + literal = _literal_value(node) + is_complex = is_complex or isinstance(literal, complex) + return (lambda values, literal=literal: jnp.asarray(literal)), repr(literal) + + fn, key = build(expression) + return ParameterExpression(tuple(slot_of(ref) for ref in references), key, is_complex, fn) @dataclass(frozen=True, slots=True) @@ -134,40 +203,22 @@ class ParametricGate: """A parametric gate whose matrix depends on runtime parameters. Calling an instance with the flat parameter vector returns the gate's ``qx.Unitary``. The - constructor and parameter layout are exposed so that gates of the same kind can be built + constructor and argument layout are exposed so that gates of the same kind can be built together in one vectorised operation. - A gate argument is either a literal number or an affine function ``scale * theta + offset`` - of one slot of the parameter vector -- which covers a bare memory reference (``RX(theta[0])``) - as well as the arithmetic quilc emits when it compiles parametric programs - (``RX(theta[0]/2 + pi)``). - :param gate_fn: The quax gate constructor (e.g. ``qx.gates.RX``), or a parametric ``DEFGATE`` callable. - :param param_indices: For each gate argument, its slot in the flat parameter vector, or - ``-1`` when the argument is a literal number. Gates that read the same memory - reference share a slot; see :func:`expand_program`. - :param concrete_values: For each gate argument, its literal value (``nan`` for a slot). - :param scales: For each gate argument, the factor multiplying the slot value (``1`` for a - bare reference; unused for a literal). - :param offsets: For each gate argument, the constant added to the scaled slot value (``0`` - for a bare reference; unused for a literal). + :param arguments: One entry per gate argument: a literal number, or a + :class:`ParameterExpression` reading the parameter vector. Gates that read the same + memory reference share a slot; see :func:`expand_program`. """ gate_fn: Callable[..., qx.Operator] - param_indices: tuple[int, ...] - concrete_values: tuple[float, ...] - scales: tuple[float, ...] - offsets: tuple[float, ...] + arguments: tuple[float | complex | ParameterExpression, ...] def __call__(self, params: Array) -> qx.Unitary: """Build the gate for one parameter vector.""" - resolved: list[Any] = [ - params[pi] * scale + offset if pi >= 0 else cv - for pi, cv, scale, offset in zip( - self.param_indices, self.concrete_values, self.scales, self.offsets, strict=True - ) - ] + resolved: list[Any] = [arg(params) if isinstance(arg, ParameterExpression) else arg for arg in self.arguments] result = self.gate_fn(*resolved) if not isinstance(result, qx.Unitary): result = qx.Unitary.from_matrix(result.matrix, result.dims) @@ -393,48 +444,34 @@ def _resolve_gate(inst: Gate) -> tuple[ExpandedOp, tuple[int, ...]]: if isinstance(gate_def, qx.Unitary): raise ValueError(f"Gate '{gate_name}' is not parametric but {inst.out()!r} passes parameters.") - param_indices: list[int] = [] - concrete_values: list[float] = [] - scales: list[float] = [] - offsets: list[float] = [] + arguments: list[float | complex | ParameterExpression] = [] for p in inst.params: - form = _affine_form(p) - if form is None: - # Each ParametricGate argument is an affine function of a single slot of the - # parameter vector; that is what lets the simulator batch same-shaped gates - # under one ``jax.vmap`` with per-gate scale and offset arrays. Anything - # else -- a product of two references, SIN(theta), theta^2 -- would need - # its own traced graph. - mrefs = _contained_mrefs(p) # type: ignore[arg-type] - raise ValueError( - f"Gate parameter {p} in {inst.out()!r} is not an affine expression " - f"(a * theta + b) in a single memory reference; it involves " - f"{sorted(str(m) for m in mrefs)}. Only such expressions are supported: " - "rewrite the program, or substitute concrete values before simulating." - ) - ref, scale, offset = form - if ref is None: - # A concrete number: a compile-time constant for this gate. - param_indices.append(-1) - concrete_values.append(offset) - scales.append(1.0) - offsets.append(0.0) - elif ref.name in measure_regs: + if not _contained_mrefs(p): # type: ignore[arg-type] + # A literal: a compile-time constant for this gate. + arguments.append(_literal_value(p)) + continue + feed_forward = [m for m in _contained_mrefs(p) if m.name in measure_regs] # type: ignore[arg-type] + if feed_forward: # Classically-conditioned angle: the value is only known mid-circuit. raise ValueError( f"Gate parameter {p} in {inst.out()!r} reads memory region " - f"'{ref.name}', which is written by a MEASURE in this program. " + f"'{feed_forward[0].name}', which is written by a MEASURE in this program. " "Feed-forward (classically-conditioned) parameters are not supported." ) - else: - param_indices.append(slots.setdefault((ref.name, ref.offset), len(slots))) - concrete_values.append(float("nan")) - scales.append(scale) - offsets.append(offset) - - return ParametricGate( - gate_def, tuple(param_indices), tuple(concrete_values), tuple(scales), tuple(offsets) - ), qubits + expression = _compile_expression(p, lambda ref: slots.setdefault((ref.name, ref.offset), len(slots))) + if ( + expression.is_complex + and gate_name in qx.gates.QUANTUM_GATES + and gate_name not in (custom_gates or {}) + ): + raise ValueError( + f"Gate parameter {p} in {inst.out()!r} is complex-valued (it contains CIS or a complex " + f"literal), but the built-in gate {gate_name} takes real angles. Complex arguments are " + "only supported for DEFGATE gates." + ) + arguments.append(expression) + + return ParametricGate(gate_def, tuple(arguments)), qubits # Fixed gate → resolve to Unitary now. unitary = get_instruction_unitary(inst, custom_gates=custom_gates) diff --git a/pyquil/simulation/_simulator.py b/pyquil/simulation/_simulator.py index 4db774398..640be76c3 100644 --- a/pyquil/simulation/_simulator.py +++ b/pyquil/simulation/_simulator.py @@ -75,6 +75,7 @@ from pyquil.simulation._circuit import Circuit, CircuitOp, Group, MergePlan from pyquil.simulation._resolver import ( ExpandedOp, + ParameterExpression, ParameterRef, ParametricGate, Resolution, @@ -598,10 +599,10 @@ def _embed_op_to_group( @dataclass class _GateBatch: - """A set of gates sharing one constructor, concrete layout, and embedding. + """A set of gates sharing one constructor, argument layout, and embedding. Members differ only in which entries of the parameter vector feed their - free arguments, so all of them are built with a single ``jax.vmap``. This + expression-valued arguments, so all of them are built with a single ``jax.vmap``. This keeps the traced graph proportional to the number of distinct gate *kinds* rather than the number of gates. @@ -612,8 +613,12 @@ class _GateBatch: gate_fn: Callable[..., qx.Operator] n_args: int - #: ``(slot, value)`` for each compile-time-constant argument. - concrete_args: tuple[tuple[int, float], ...] + #: ``(position, value)`` for each literal argument. + literal_args: tuple[tuple[int, float | complex], ...] + #: ``(position, expression)`` for each expression-valued argument. Every member's + #: expressions have the same keys, so the first member's serve as the template: they are + #: evaluated on each member's own slot values. + expression_args: tuple[tuple[int, ParameterExpression], ...] #: Per-qudit dimensions of the merge group each member embeds into. target_dims: tuple[int, ...] #: Positions within the group occupied by the gate's qudits. @@ -624,38 +629,38 @@ class _GateBatch: as_superop: bool #: Sorted-array positions this batch fills, one per member. positions: list[int] = field(default_factory=list) - #: Parameter-vector index for each free argument, one list per member. - param_indices: list[list[int]] = field(default_factory=list) - #: Affine coefficients of each free argument (``scale * params[index] + offset``), one list - #: per member. They are per-member data rather than part of the batch key, so - #: ``RX(theta[0] / 2)`` and ``RX(theta[1])`` share one vmap. - param_scales: list[list[float]] = field(default_factory=list) - param_offsets: list[list[float]] = field(default_factory=list) + #: Parameter-vector slots read by each member: the slot indices of its expression + #: arguments, concatenated in argument order. + slot_indices: list[list[int]] = field(default_factory=list) def builder(self) -> Callable[[Array], Array]: """Return ``params -> (n_members, width, width)`` embedded gate matrices.""" - concrete = {slot for slot, _ in self.concrete_args} - free_slots = [j for j in range(self.n_args) if j not in concrete] - gate_fn, n_args, concrete_args = self.gate_fn, self.n_args, self.concrete_args + gate_fn, n_args, literal_args, expression_args = ( + self.gate_fn, + self.n_args, + self.literal_args, + self.expression_args, + ) target_dims, group_positions = self.target_dims, self.group_positions width, as_superop = self.width, self.as_superop - param_indices = jnp.asarray(self.param_indices) # (n_members, n_free) - scales = jnp.asarray(self.param_scales, dtype=float) - offsets = jnp.asarray(self.param_offsets, dtype=float) + slot_indices = jnp.asarray(self.slot_indices, dtype=jnp.int32) # (n_members, n_slots) - def single(free_values: Array) -> Array: + def single(values: Array) -> Array: args: list[Any] = [None] * n_args - for slot, val in concrete_args: - args[slot] = val - for k, slot in enumerate(free_slots): - args[slot] = free_values[k] + for position, value in literal_args: + args[position] = value + offset = 0 + for position, expression in expression_args: + count = len(expression.slot_indices) + args[position] = expression.evaluate(values[offset : offset + count]) + offset += count gate = gate_fn(*args) if not isinstance(gate, qx.Unitary): gate = qx.Unitary.from_matrix(gate.matrix, gate.dims) return _embed_op_to_group(gate, target_dims, group_positions, width, as_superop=as_superop) batched = jax.vmap(single) - return lambda params: batched(params[param_indices] * scales + offsets) + return lambda params: batched(params[slot_indices]) def _make_group_fold(group_start: list[int], n_ops: int, width: int) -> Callable[[Array], Array]: @@ -755,14 +760,22 @@ def _build_vectorized_operator_constructor( # Key by embedding *type* (op dims + group dims + positions), not # physical qubits: embeddings that trace to the same graph share a vmap. embed_key = (tuple(dims[q] for q in op_sub), target_dims, group_positions) - concrete_args = tuple((j, op.concrete_values[j]) for j, pi in enumerate(op.param_indices) if pi < 0) - key = (id(op.gate_fn), concrete_args, embed_key) + literal_args = tuple( + (j, arg) for j, arg in enumerate(op.arguments) if not isinstance(arg, ParameterExpression) + ) + expression_args = tuple( + (j, arg) for j, arg in enumerate(op.arguments) if isinstance(arg, ParameterExpression) + ) + # Expressions enter the key by *shape* only, so ``SIN(theta[0])`` and ``SIN(theta[1])`` + # share a batch and differ in the slots they read. + key = (id(op.gate_fn), literal_args, tuple((j, arg.key) for j, arg in expression_args), embed_key) batch = batches.get(key) if batch is None: batch = _GateBatch( gate_fn=op.gate_fn, - n_args=len(op.param_indices), - concrete_args=concrete_args, + n_args=len(op.arguments), + literal_args=literal_args, + expression_args=expression_args, target_dims=target_dims, group_positions=group_positions, width=width, @@ -770,10 +783,7 @@ def _build_vectorized_operator_constructor( ) batches[key] = batch batch.positions.append(pos) - free = [j for j, pi in enumerate(op.param_indices) if pi >= 0] - batch.param_indices.append([op.param_indices[j] for j in free]) - batch.param_scales.append([op.scales[j] for j in free]) - batch.param_offsets.append([op.offsets[j] for j in free]) + batch.slot_indices.append([slot for _, arg in expression_args for slot in arg.slot_indices]) else: # Constant operations are embedded once, eagerly. In superoperator mode this # also covers the non-unitary ops a noise model contributes (channel SuperOps, diff --git a/test/unit/test_affine_parameters.py b/test/unit/test_affine_parameters.py deleted file mode 100644 index 259cfbd9d..000000000 --- a/test/unit/test_affine_parameters.py +++ /dev/null @@ -1,136 +0,0 @@ -############################################################################## -# Copyright 2026 Rigetti Computing -# -# 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 -# -# http://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. -############################################################################## -"""Affine gate-parameter expressions (``a * theta + b``) in the quax-based simulators.""" - -import jax -import jax.numpy as jnp -import numpy as np -import pytest - -from pyquil.gates import RX, RY, RZ, H -from pyquil.quil import Program -from pyquil.quilatom import MemoryReference, quil_sin -from pyquil.quilbase import Declare -from pyquil.simulation._resolver import _affine_form, expand_program -from pyquil.simulation._simulator import DensityMatrixSimulator, PureStateVectorSimulator - -THETA = MemoryReference("theta", 0) -PHI = MemoryReference("theta", 1) - - -def _program(*gates): - program = Program() - program += Declare("theta", "REAL", 2) - for gate in gates: - program += gate - return program - - -def _state(program, **memory): - sim = PureStateVectorSimulator(program) - return np.asarray(sim.compute(sim.linearize(memory)).matrix).reshape(-1) - - -class TestAffineForm: - @pytest.mark.parametrize( - ("expression", "expected"), - [ - (THETA, (THETA, 1.0, 0.0)), - (THETA / 2 + np.pi, (THETA, 0.5, np.pi)), - (2 * PHI, (PHI, 2.0, 0.0)), - (-THETA, (THETA, -1.0, 0.0)), - (np.pi - PHI, (PHI, -1.0, np.pi)), - (THETA + THETA, (THETA, 2.0, 0.0)), - ((3 * THETA - 1) / 4, (THETA, 0.75, -0.25)), - (0.7, (None, 0.0, 0.7)), - ], - ) - def test_recognises_affine_expressions(self, expression, expected): - ref, scale, offset = _affine_form(expression) - assert ref == expected[0] - assert scale == pytest.approx(expected[1]) - assert offset == pytest.approx(expected[2]) - - @pytest.mark.parametrize( - "expression", - [THETA * PHI, THETA * THETA, quil_sin(THETA), THETA**2, 2 / THETA, THETA + PHI], - ids=["product", "square", "sin", "power", "reciprocal", "two_references"], - ) - def test_rejects_non_affine_expressions(self, expression): - assert _affine_form(expression) is None - - -class TestSimulation: - @pytest.mark.parametrize( - ("gate", "literal"), - [ - (RX(THETA / 2 + np.pi, 0), lambda t, p: RX(t / 2 + np.pi, 0)), - (RZ(2 * PHI, 0), lambda t, p: RZ(2 * p, 0)), - (RX(-THETA, 0), lambda t, p: RX(-t, 0)), - (RY(np.pi - PHI, 0), lambda t, p: RY(np.pi - p, 0)), - (RX(THETA + THETA, 0), lambda t, p: RX(2 * t, 0)), - ], - ids=["half_plus_pi", "double", "negated", "pi_minus", "same_reference_twice"], - ) - def test_matches_the_literal_gate(self, gate, literal): - theta, phi = 0.37, 1.21 - got = _state(_program(H(0), gate), theta=[theta, phi]) - expected = _state(Program(H(0), literal(theta, phi))) - np.testing.assert_allclose(got, expected, atol=1e-12) - - def test_parsed_quilc_style_program(self): - """The arithmetic quilc emits for a compiled parametric program simulates directly.""" - program = Program("DECLARE theta REAL[1]\nRX((theta[0]/2)+pi) 0\nRZ(-2*theta[0]) 0") - theta = 0.6 - got = _state(program, theta=[theta]) - expected = _state(Program(RX(theta / 2 + np.pi, 0), RZ(-2 * theta, 0))) - np.testing.assert_allclose(got, expected, atol=1e-12) - - def test_shares_a_slot_and_batches_with_a_plain_reference(self): - """``RX(theta[0]/2)`` and ``RX(theta[1])`` use one slot each and one gate batch.""" - program = _program(RX(THETA / 2, 0), RX(PHI, 1), RX(THETA, 2)) - sim = PureStateVectorSimulator(program) - assert sim.parameters == (("theta", 0), ("theta", 1)) - theta, phi = 0.5, 0.9 - got = _state(program, theta=[theta, phi]) - expected = _state(Program(RX(theta / 2, 0), RX(phi, 1), RX(theta, 2))) - np.testing.assert_allclose(got, expected, atol=1e-12) - - def test_grad_and_jit(self): - """d/dtheta P(|1>) for RX(theta/2 + pi) is -sin(theta/2)/4 (the pi flips the population).""" - sim = DensityMatrixSimulator(_program(RX(THETA / 2 + np.pi, 0))) - - def excited(params): - return jnp.real(sim.compute(params).matrix[1, 1]) - - theta = 0.8 - params = jnp.array([theta]) - np.testing.assert_allclose(jax.jit(excited)(params), np.cos(theta / 4) ** 2, atol=1e-12) - grad = jax.grad(excited)(params) - np.testing.assert_allclose(grad, [-np.sin(theta / 2) / 4], atol=1e-12) - - def test_expand_program_records_scale_and_offset(self): - ops, _, parameters = expand_program(_program(RX(3 * THETA - 1, 0))) - (gate,) = ops - assert parameters == (("theta", 0),) - assert gate.param_indices == (0,) - assert gate.scales == (3.0,) - assert gate.offsets == (-1.0,) - - @pytest.mark.parametrize("expression", [THETA * PHI, quil_sin(THETA), THETA**2], ids=["product", "sin", "power"]) - def test_non_affine_parameter_reports_clearly(self, expression): - with pytest.raises(ValueError, match="not an affine expression"): - PureStateVectorSimulator(_program(RX(expression, 0))) diff --git a/test/unit/test_density_matrix.py b/test/unit/test_density_matrix.py index d1771c522..f412fbd11 100644 --- a/test/unit/test_density_matrix.py +++ b/test/unit/test_density_matrix.py @@ -58,7 +58,7 @@ ) from pyquil.noise._noise_model import NoiseModel from pyquil.quil import Program -from pyquil.quilatom import MemoryReference +from pyquil.quilatom import MemoryReference, quil_cis from pyquil.quilbase import Declare, Gate, ResetQubit from pyquil.simulation._reference import ReferenceDensitySimulator, ReferenceWavefunctionSimulator from pyquil.simulation._simulator import DensityMatrixSimulator @@ -684,10 +684,10 @@ def test_unsupported_modifier_reports_clearly(self, program, modifier): with pytest.raises(ValueError, match=f"modifiers are not supported.*{modifier}"): _dm(program) - def test_non_affine_parameter_reports_clearly(self): + def test_complex_valued_parameter_for_builtin_gate_reports_clearly(self): theta = MemoryReference("theta", 0) - program = Program(Declare("theta", "REAL", 1), RX(theta * theta, 0)) - with pytest.raises(ValueError, match="not an affine expression"): + program = Program(Declare("theta", "REAL", 1), RX(quil_cis(theta), 0)) + with pytest.raises(ValueError, match="complex-valued"): DensityMatrixSimulator(program, qubits=[0]) def test_feed_forward_parameter_reports_clearly(self): diff --git a/test/unit/test_parameter_expressions.py b/test/unit/test_parameter_expressions.py new file mode 100644 index 000000000..2c933b9b0 --- /dev/null +++ b/test/unit/test_parameter_expressions.py @@ -0,0 +1,233 @@ +############################################################################## +# Copyright 2026 Rigetti Computing +# +# 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 +# +# http://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. +############################################################################## +"""Quil arithmetic expressions as gate arguments in the quax-based simulators.""" + +import jax +import jax.numpy as jnp +import numpy as np +import pytest + +from pyquil.gates import PHASE, RX, RY, RZ, H +from pyquil.quil import Program +from pyquil.quilatom import MemoryReference, Parameter, quil_cis, quil_cos, quil_exp, quil_sin, quil_sqrt +from pyquil.quilbase import Declare, DefGate, Gate +from pyquil.simulation._resolver import ParameterExpression, _compile_expression, expand_program +from pyquil.simulation._simulator import DensityMatrixSimulator, PureStateVectorSimulator + +THETA = MemoryReference("theta", 0) +PHI = MemoryReference("theta", 1) + + +def _program(*gates): + program = Program() + program += Declare("theta", "REAL", 2) + for gate in gates: + program += gate + return program + + +def _state(program, **memory): + sim = PureStateVectorSimulator(program) + return np.asarray(sim.compute(sim.linearize(memory)).matrix).reshape(-1) + + +def _compile(expression): + return _compile_expression(expression, lambda ref: ref.offset) + + +class TestCompileExpression: + @pytest.mark.parametrize( + ("expression", "key", "slots"), + [ + (THETA, "%0", (0,)), + (THETA / 2 + np.pi, "((%0/2.0)+3.141592653589793)", (0,)), + (2 * PHI, "(2.0*%0)", (1,)), + (-THETA, "(-1.0*%0)", (0,)), + (THETA + THETA, "(%0+%0)", (0,)), + (THETA * PHI, "(%0*%1)", (0, 1)), + (PHI - THETA, "(%0-%1)", (1, 0)), + (quil_sin(THETA), "SIN(%0)", (0,)), + (quil_cos(2 * THETA), "COS((2.0*%0))", (0,)), + (quil_sqrt(THETA), "SQRT(%0)", (0,)), + (quil_exp(-THETA), "EXP((-1.0*%0))", (0,)), + (quil_cis(THETA), "CIS(%0)", (0,)), + (THETA**2, "(%0^2.0)", (0,)), + ], + ) + def test_key_and_slots(self, expression, key, slots): + compiled = _compile(expression) + assert compiled.key == key + assert compiled.slot_indices == slots + + def test_same_shape_different_slots_share_a_key(self): + assert _compile(quil_sin(THETA)).key == _compile(quil_sin(PHI)).key + assert _compile(THETA / 2).key == _compile(PHI / 2).key + assert _compile(THETA / 2).key != _compile(THETA / 3).key + + @pytest.mark.parametrize( + ("expression", "is_complex"), + [(THETA, False), (quil_sin(THETA), False), (THETA**2, False), (quil_cis(THETA), True), (1j * THETA, True)], + ) + def test_complex_flag(self, expression, is_complex): + assert _compile(expression).is_complex is is_complex + + @pytest.mark.parametrize( + ("expression", "theta", "expected"), + [ + (THETA / 2 + np.pi, 0.4, 0.2 + np.pi), + (quil_sin(THETA) * quil_cos(THETA), 0.7, np.sin(0.7) * np.cos(0.7)), + (quil_sqrt(THETA), 0.81, 0.9), + (quil_exp(-THETA), 1.5, np.exp(-1.5)), + (THETA**2 - 2**THETA, 1.5, 1.5**2 - 2**1.5), + (quil_cis(THETA), 0.3, np.exp(0.3j)), + ((1 + 2j) * THETA, 0.5, 0.5 + 1j), + (quil_sqrt(2) * THETA, 1.0, np.sqrt(2)), + ], + ) + def test_evaluates_like_quil(self, expression, theta, expected): + compiled = _compile(expression) + np.testing.assert_allclose(compiled(jnp.array([theta, 0.0])), expected, atol=1e-12) + + def test_unbound_defgate_parameter_reports_clearly(self): + with pytest.raises(ValueError, match="Unbound DEFGATE parameter"): + _compile(THETA + Parameter("p")) + + +class TestSimulation: + @pytest.mark.parametrize( + ("gate", "literal"), + [ + (RX(THETA / 2 + np.pi, 0), lambda t, p: RX(t / 2 + np.pi, 0)), + (RZ(2 * PHI, 0), lambda t, p: RZ(2 * p, 0)), + (RX(-THETA, 0), lambda t, p: RX(-t, 0)), + (RY(np.pi - PHI, 0), lambda t, p: RY(np.pi - p, 0)), + (RX(THETA + THETA, 0), lambda t, p: RX(2 * t, 0)), + (RX(THETA * PHI, 0), lambda t, p: RX(t * p, 0)), + (RX(quil_sin(THETA), 0), lambda t, p: RX(np.sin(t), 0)), + (RY(np.pi * quil_cos(PHI), 0), lambda t, p: RY(np.pi * np.cos(p), 0)), + (RZ(quil_sqrt(THETA), 0), lambda t, p: RZ(np.sqrt(t), 0)), + (RX(quil_exp(-THETA), 0), lambda t, p: RX(np.exp(-t), 0)), + (RX(THETA**2, 0), lambda t, p: RX(t**2, 0)), + (RX(quil_sqrt(2) * THETA, 0), lambda t, p: RX(np.sqrt(2) * t, 0)), + ], + ids=[ + "half_plus_pi", + "double", + "negated", + "pi_minus", + "same_reference_twice", + "product_of_references", + "sin", + "pi_cos", + "sqrt", + "exp", + "power", + "sqrt2_constant", + ], + ) + def test_matches_the_literal_gate(self, gate, literal): + theta, phi = 0.37, 1.21 + got = _state(_program(H(0), gate), theta=[theta, phi]) + expected = _state(Program(H(0), literal(theta, phi))) + np.testing.assert_allclose(got, expected, atol=1e-12) + + def test_parsed_quilc_style_program(self): + """The arithmetic quilc emits for a compiled parametric program simulates directly.""" + program = Program( + "DECLARE theta REAL[1]\nRX((theta[0]/2)+pi) 0\nRZ(-2*theta[0]) 0\nRY(SIN(theta[0])*COS(pi/3)) 0" + ) + theta = 0.6 + got = _state(program, theta=[theta]) + expected = _state( + Program(RX(theta / 2 + np.pi, 0), RZ(-2 * theta, 0), RY(np.sin(theta) * np.cos(np.pi / 3), 0)) + ) + np.testing.assert_allclose(got, expected, atol=1e-12) + + def test_cis_and_complex_literals_in_a_defgate_argument(self): + """A DEFGATE with a complex parameter accepts CIS and complex arithmetic.""" + z = Parameter("z") + defgate = DefGate("CPH", [[1, 0], [0, z]], [z]) + theta = 0.8 + + def run(argument): + program = _program() + program += defgate + program += H(0) + program += Gate("CPH", [argument], [0]) + return _state(program, theta=[theta, 0.0]) + + expected = _state(Program(H(0), PHASE(theta, 0))) + np.testing.assert_allclose(run(quil_cis(THETA)), expected, atol=1e-12) + # 1i * CIS(theta) == CIS(theta + pi/2) + np.testing.assert_allclose( + run(1j * quil_cis(THETA)), _state(Program(H(0), PHASE(theta + np.pi / 2, 0))), atol=1e-12 + ) + # A complex literal times a real reference, folded through a real result: (1i)*(-1i) == 1. + np.testing.assert_allclose(run(quil_cis(1j * (-1j) * THETA)), expected, atol=1e-12) + + def test_complex_argument_to_a_builtin_gate_reports_clearly(self): + with pytest.raises(ValueError, match="complex-valued"): + PureStateVectorSimulator(_program(RX(quil_cis(THETA), 0))) + with pytest.raises(ValueError, match="complex-valued"): + PureStateVectorSimulator(_program(RX(1j * THETA, 0))) + + def test_feed_forward_parameter_still_rejected(self): + from pyquil.gates import MEASURE + + program = Program() + program += Declare("ro", "BIT", 1) + program += MEASURE(0, ("ro", 0)) + program += RX(2 * MemoryReference("ro", 0), 0) + with pytest.raises(ValueError, match="written by a MEASURE"): + DensityMatrixSimulator(program) + + def test_same_shape_expressions_share_one_slot_each_and_batch(self): + program = _program(RX(quil_sin(THETA), 0), RX(quil_sin(PHI), 1), RX(THETA, 2)) + sim = PureStateVectorSimulator(program) + assert sim.parameters == (("theta", 0), ("theta", 1)) + theta, phi = 0.5, 0.9 + got = _state(program, theta=[theta, phi]) + expected = _state(Program(RX(np.sin(theta), 0), RX(np.sin(phi), 1), RX(theta, 2))) + np.testing.assert_allclose(got, expected, atol=1e-12) + + def test_grad_and_jit(self): + """d/dtheta P(|1>) for RX(SIN(theta)) is sin(sin theta) cos(theta) / 2.""" + sim = DensityMatrixSimulator(_program(RX(quil_sin(THETA), 0))) + + def excited(params): + return jnp.real(sim.compute(params).matrix[1, 1]) + + theta = 0.8 + params = jnp.array([theta]) + np.testing.assert_allclose(jax.jit(excited)(params), np.sin(np.sin(theta) / 2) ** 2, atol=1e-12) + grad = jax.grad(excited)(params) + np.testing.assert_allclose(grad, [np.sin(np.sin(theta)) * np.cos(theta) / 2], atol=1e-12) + + def test_expand_program_records_the_expression(self): + ops, _, parameters = expand_program(_program(RX(3 * THETA - 1, 0))) + (gate,) = ops + assert parameters == (("theta", 0),) + (argument,) = gate.arguments + assert isinstance(argument, ParameterExpression) + assert argument.slot_indices == (0,) + assert argument.key == "((3.0*%0)-1.0)" + + def test_literal_expressions_are_folded(self): + """A parameter-free expression such as SIN(pi/4) is a literal, not a parameter.""" + program = Program("RX(SIN(pi/4)) 0") + sim = PureStateVectorSimulator(program) + assert sim.parameters == () + np.testing.assert_allclose(_state(program), _state(Program(RX(np.sin(np.pi / 4), 0))), atol=1e-12) diff --git a/test/unit/test_resolver.py b/test/unit/test_resolver.py index 496e29e86..5dd911963 100644 --- a/test/unit/test_resolver.py +++ b/test/unit/test_resolver.py @@ -430,7 +430,7 @@ def test_repeated_reference_shares_a_slot(self): program = Program(Declare("theta", "REAL", 1), RX(theta0, 0), RZ(theta0, 1)) ops, _, parameters = expand_program(program) assert parameters == (("theta", 0),) - assert [op.param_indices for op in ops] == [(0,), (0,)] + assert [op.arguments[0].slot_indices for op in ops] == [(0,), (0,)] def test_distinct_references_get_distinct_slots_in_first_use_order(self): program = Program( From 908432d49ff34660f9cde9004a96613a4b17bd0a Mon Sep 17 00:00:00 2001 From: Bram Evert Date: Mon, 14 Sep 2026 10:06:30 +0000 Subject: [PATCH 3/3] Address the review of the expression compiler Eric's review of #1869: * ``ParameterExpression._fn`` was private for no reason. The compiled closure is now the public ``evaluate`` field, documented; the one-line wrapper method is gone. * Name the two entry points' arguments for what they are: ``evaluate`` takes the narrowed vector of just the slots the expression reads, ``__call__`` takes the circuit-wide parameter vector and gathers from it. ``_GateBatch.single`` takes ``slot_values``. * Lift the ``build`` closure to a module-level recursive ``_build_expression(node, references)``: the reference accumulator is an explicit argument and ``is_complex`` propagates up by return value, so the ``nonlocal`` is gone. Compiler errors now name the node, and ``expand_program`` adds the instruction that contains it. Also, from a read-through of the branch: * A complex *literal* beside an expression argument escaped the complex-argument guard, so ``FSIM(theta[0], 1i) 0 1`` was accepted and failed much later with quax's ``TypeError: evolve is not implemented for Operator``. Both argument paths are checked now, against a single ``is_builtin`` decided where the gate is resolved rather than recomputed from ``QUANTUM_GATES``. * ``ParameterExpression``'s docstring offered ``CPHASE(CIS(theta[0]))`` as the DEFGATE example, but CPHASE is built-in, so that exact program raises. Use a DEFGATE name. * Document in the architecture guide, as the docstring already does, that SQRT and ^ are evaluated in real arithmetic and give nan where Quil would give a complex result. Co-Authored-By: Claude Opus 5 (1M context) Claude-Session: https://claude.ai/code/session_01M6WgBnWYqp3DecvdG8pH2b --- docs/source/simulation_architecture.rst | 4 +- pyquil/simulation/_resolver.py | 146 ++++++++++++++---------- pyquil/simulation/_simulator.py | 7 +- test/unit/test_parameter_expressions.py | 11 +- 4 files changed, 103 insertions(+), 65 deletions(-) diff --git a/docs/source/simulation_architecture.rst b/docs/source/simulation_architecture.rst index 5b6d933ea..4bc1d2b3d 100644 --- a/docs/source/simulation_architecture.rst +++ b/docs/source/simulation_architecture.rst @@ -222,7 +222,9 @@ Expansion does several things at once: compiles a parametric program, and ``RZ(2*SIN(phi[1]))`` both simulate directly. Gates whose expressions have the same *shape* (``SIN(theta[0])`` and ``SIN(theta[1])``) are still built in one vectorised batch. A complex-valued argument (``CIS``, a complex literal) is accepted only - by a ``DEFGATE`` gate, since the built-in gates take real angles. + by a ``DEFGATE`` gate, since the built-in gates take real angles; everything else is evaluated + in real arithmetic, so ``SQRT`` of a negative number or a fractional power of one gives ``nan`` + where Quil would give a complex result. * **DEFCIRCUIT and cycle expansion.** ``DEFCIRCUIT`` bodies are expanded with formal-argument substitution. When a circuit invocation matches a diff --git a/pyquil/simulation/_resolver.py b/pyquil/simulation/_resolver.py index ebd6a241e..187c56bb1 100644 --- a/pyquil/simulation/_resolver.py +++ b/pyquil/simulation/_resolver.py @@ -121,9 +121,13 @@ class ParameterExpression: Any expression Quil allows is supported: ``+ - * / ^``, the functions ``SIN``, ``COS``, ``SQRT``, ``EXP`` and ``CIS``, and real or complex literals -- for example - ``RX(theta[0]/2 + pi)``, ``RZ(2*SIN(phi[1]))`` or ``CPHASE(CIS(theta[0])) 0 1`` for a - ``DEFGATE`` taking a complex parameter. Calling an instance with the flat parameter - vector evaluates it with JAX, so it can be jitted and differentiated through. + ``RX(theta[0]/2 + pi)``, ``RZ(2*SIN(phi[1]))``, or ``CPH(CIS(theta[0])) 0`` for a + ``DEFGATE CPH(%z)`` taking a complex parameter. Evaluation goes through JAX, so an + expression can be jitted and differentiated through. + + There are two ways in. :attr:`evaluate` takes the *narrowed* vector of only the values this + expression reads, which is what the simulator gathers once per vectorised batch; calling the + instance takes the whole circuit-wide parameter vector and gathers from it first. :param slot_indices: Slots of the parameter vector the expression reads, in order of first appearance. @@ -135,19 +139,21 @@ class ParameterExpression: complex literal). Everything else is evaluated in real arithmetic; note that Quil would evaluate ``SQRT`` of a negative number or a fractional power of one as complex, which real arithmetic reports as ``nan``. + :param evaluate: Evaluates the expression from the values of :attr:`slot_indices`, in that + order. """ slot_indices: tuple[int, ...] key: str is_complex: bool - _fn: Callable[[Array], Array] - - def evaluate(self, values: Array) -> Array: - """Evaluate the expression given the values of its slots, in :attr:`slot_indices` order.""" - return self._fn(values) + evaluate: Callable[[Array], Array] def __call__(self, params: Array) -> Array: - """Evaluate the expression for one parameter vector.""" + """Evaluate the expression from the whole circuit-wide parameter vector. + + Gathers :attr:`slot_indices` out of ``params`` and hands the narrowed vector to + :attr:`evaluate`. + """ return self.evaluate(params[jnp.asarray(self.slot_indices)]) @@ -157,45 +163,58 @@ def _literal_value(expression: Any) -> float | complex: return value.real if value.imag == 0 else value +def _build_expression(node: Any, references: list[MemoryReference]) -> tuple[Callable[[Array], Array], str, bool]: + """Compile one node of a Quil expression tree. + + :param node: The node to compile; recursion handles its operands. + :param references: The memory references seen so far, in first-appearance order. New ones are + appended, and a reference's position here is the index into the narrowed value vector that + the compiled closure reads. + :returns: The closure, the shape key, and whether the value may be complex. + :raises ValueError: If the node is an unbound ``DEFGATE`` parameter or an unknown function. + The message names only the offending node; :func:`expand_program` adds the instruction. + """ + if isinstance(node, MemoryReference): + if node not in references: + references.append(node) + index = references.index(node) + return (lambda slot_values, index=index: slot_values[index]), f"%{index}", False + if isinstance(node, Parameter): + raise ValueError(f"Unbound DEFGATE parameter {node}.") + if isinstance(node, BinaryExp): + left, left_key, left_complex = _build_expression(node.op1, references) + right, right_key, right_complex = _build_expression(node.op2, references) + operator = _BINARY_OPERATORS[type(node)] + return ( + (lambda slot_values: operator(left(slot_values), right(slot_values))), + f"({left_key}{node.operator.strip()}{right_key})", + left_complex or right_complex, + ) + if isinstance(node, Function): + if node.name not in _FUNCTIONS: + raise ValueError(f"Unknown Quil function {node.name!r}.") + function = _FUNCTIONS[node.name] + inner, inner_key, inner_complex = _build_expression(node.expression, references) + return ( + (lambda slot_values: function(inner(slot_values))), + f"{node.name}({inner_key})", + inner_complex or node.name == "CIS", + ) + literal = _literal_value(node) + return (lambda slot_values, literal=literal: jnp.asarray(literal)), repr(literal), isinstance(literal, complex) + + def _compile_expression(expression: Any, slot_of: Callable[[MemoryReference], int]) -> ParameterExpression: """Compile a Quil expression over memory references into a :class:`ParameterExpression`. :param expression: The gate parameter; must contain at least one memory reference. :param slot_of: Maps a memory reference to its slot in the parameter vector. - :raises ValueError: If the expression contains an unbound ``DEFGATE`` parameter. + :raises ValueError: If the expression contains an unbound ``DEFGATE`` parameter or an unknown + function. """ references: list[MemoryReference] = [] - is_complex = False - - def build(node: Any) -> tuple[Callable[[Array], Array], str]: - nonlocal is_complex - if isinstance(node, MemoryReference): - if node not in references: - references.append(node) - index = references.index(node) - return (lambda values, index=index: values[index]), f"%{index}" - if isinstance(node, Parameter): - raise ValueError(f"Unbound DEFGATE parameter {node} in gate argument {expression}.") - if isinstance(node, BinaryExp): - left, left_key = build(node.op1) - right, right_key = build(node.op2) - operator = _BINARY_OPERATORS[type(node)] - return ( - lambda values: operator(left(values), right(values)) - ), f"({left_key}{node.operator.strip()}{right_key})" - if isinstance(node, Function): - inner, inner_key = build(node.expression) - if node.name not in _FUNCTIONS: - raise ValueError(f"Unknown Quil function {node.name!r} in gate argument {expression}.") - function = _FUNCTIONS[node.name] - is_complex = is_complex or node.name == "CIS" - return (lambda values: function(inner(values))), f"{node.name}({inner_key})" - literal = _literal_value(node) - is_complex = is_complex or isinstance(literal, complex) - return (lambda values, literal=literal: jnp.asarray(literal)), repr(literal) - - fn, key = build(expression) - return ParameterExpression(tuple(slot_of(ref) for ref in references), key, is_complex, fn) + evaluate, key, is_complex = _build_expression(expression, references) + return ParameterExpression(tuple(slot_of(ref) for ref in references), key, is_complex, evaluate) @dataclass(frozen=True, slots=True) @@ -436,9 +455,9 @@ def _resolve_gate(inst: Gate) -> tuple[ExpandedOp, tuple[int, ...]]: if any(_contained_mrefs(p) for p in inst.params): # type: ignore[arg-type] gate_name = inst.name if custom_gates is not None and gate_name in custom_gates: - gate_def = custom_gates[gate_name] + gate_def, is_builtin = custom_gates[gate_name], False elif gate_name in qx.gates.QUANTUM_GATES: - gate_def = qx.gates.QUANTUM_GATES[gate_name] + gate_def, is_builtin = qx.gates.QUANTUM_GATES[gate_name], True else: raise KeyError(f"Unknown gate '{gate_name}'.") if isinstance(gate_def, qx.Unitary): @@ -446,30 +465,37 @@ def _resolve_gate(inst: Gate) -> tuple[ExpandedOp, tuple[int, ...]]: arguments: list[float | complex | ParameterExpression] = [] for p in inst.params: + argument: float | complex | ParameterExpression if not _contained_mrefs(p): # type: ignore[arg-type] # A literal: a compile-time constant for this gate. - arguments.append(_literal_value(p)) - continue - feed_forward = [m for m in _contained_mrefs(p) if m.name in measure_regs] # type: ignore[arg-type] - if feed_forward: - # Classically-conditioned angle: the value is only known mid-circuit. - raise ValueError( - f"Gate parameter {p} in {inst.out()!r} reads memory region " - f"'{feed_forward[0].name}', which is written by a MEASURE in this program. " - "Feed-forward (classically-conditioned) parameters are not supported." - ) - expression = _compile_expression(p, lambda ref: slots.setdefault((ref.name, ref.offset), len(slots))) - if ( - expression.is_complex - and gate_name in qx.gates.QUANTUM_GATES - and gate_name not in (custom_gates or {}) - ): + argument = _literal_value(p) + is_complex = isinstance(argument, complex) + else: + feed_forward = [m for m in _contained_mrefs(p) if m.name in measure_regs] # type: ignore[arg-type] + if feed_forward: + # Classically-conditioned angle: the value is only known mid-circuit. + raise ValueError( + f"Gate parameter {p} in {inst.out()!r} reads memory region " + f"'{feed_forward[0].name}', which is written by a MEASURE in this program. " + "Feed-forward (classically-conditioned) parameters are not supported." + ) + try: + argument = _compile_expression( + p, lambda ref: slots.setdefault((ref.name, ref.offset), len(slots)) + ) + except ValueError as error: + # The compiler sees one expression; name the instruction it came from. + raise ValueError(f"Gate parameter {p} in {inst.out()!r}: {error}") from error + is_complex = argument.is_complex + # Checked for literals too: quax's built-in constructors take real angles, and a + # complex one otherwise surfaces much later as an opaque error from deep in quax. + if is_complex and is_builtin: raise ValueError( f"Gate parameter {p} in {inst.out()!r} is complex-valued (it contains CIS or a complex " f"literal), but the built-in gate {gate_name} takes real angles. Complex arguments are " "only supported for DEFGATE gates." ) - arguments.append(expression) + arguments.append(argument) return ParametricGate(gate_def, tuple(arguments)), qubits diff --git a/pyquil/simulation/_simulator.py b/pyquil/simulation/_simulator.py index 640be76c3..7fbf3a2f8 100644 --- a/pyquil/simulation/_simulator.py +++ b/pyquil/simulation/_simulator.py @@ -617,7 +617,7 @@ class _GateBatch: literal_args: tuple[tuple[int, float | complex], ...] #: ``(position, expression)`` for each expression-valued argument. Every member's #: expressions have the same keys, so the first member's serve as the template: they are - #: evaluated on each member's own slot values. + #: evaluated on the narrowed vector of slot values gathered for each member. expression_args: tuple[tuple[int, ParameterExpression], ...] #: Per-qudit dimensions of the merge group each member embeds into. target_dims: tuple[int, ...] @@ -645,14 +645,15 @@ def builder(self) -> Callable[[Array], Array]: width, as_superop = self.width, self.as_superop slot_indices = jnp.asarray(self.slot_indices, dtype=jnp.int32) # (n_members, n_slots) - def single(values: Array) -> Array: + def single(slot_values: Array) -> Array: + """Build one member's embedded matrix from the slot values it reads, in order.""" args: list[Any] = [None] * n_args for position, value in literal_args: args[position] = value offset = 0 for position, expression in expression_args: count = len(expression.slot_indices) - args[position] = expression.evaluate(values[offset : offset + count]) + args[position] = expression.evaluate(slot_values[offset : offset + count]) offset += count gate = gate_fn(*args) if not isinstance(gate, qx.Unitary): diff --git a/test/unit/test_parameter_expressions.py b/test/unit/test_parameter_expressions.py index 2c933b9b0..7e481f4fd 100644 --- a/test/unit/test_parameter_expressions.py +++ b/test/unit/test_parameter_expressions.py @@ -20,7 +20,7 @@ import numpy as np import pytest -from pyquil.gates import PHASE, RX, RY, RZ, H +from pyquil.gates import FSIM, PHASE, RX, RY, RZ, H from pyquil.quil import Program from pyquil.quilatom import MemoryReference, Parameter, quil_cis, quil_cos, quil_exp, quil_sin, quil_sqrt from pyquil.quilbase import Declare, DefGate, Gate @@ -184,6 +184,15 @@ def test_complex_argument_to_a_builtin_gate_reports_clearly(self): with pytest.raises(ValueError, match="complex-valued"): PureStateVectorSimulator(_program(RX(1j * THETA, 0))) + def test_complex_literal_beside_an_expression_argument_reports_clearly(self): + """A complex *literal* is rejected too, not only a complex expression.""" + with pytest.raises(ValueError, match="complex-valued"): + PureStateVectorSimulator(_program(FSIM(THETA, 1j, 0, 1))) + + def test_unbound_defgate_parameter_names_the_instruction(self): + with pytest.raises(ValueError, match=r"RX\(.*\) 0'?: Unbound DEFGATE parameter"): + PureStateVectorSimulator(_program(RX(THETA + Parameter("p"), 0))) + def test_feed_forward_parameter_still_rejected(self): from pyquil.gates import MEASURE