# Source code for pennylane.templates.subroutines.qsvt

# Copyright 2018-2023 Xanadu Quantum Technologies Inc.

# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at

# Unless required by applicable law or agreed to in writing, software
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
"""
Contains the QSVT template and qsvt wrapper function.
"""
# pylint: disable=too-many-arguments
import copy
import pennylane as qml
from pennylane.queuing import QueuingManager
from pennylane.ops import BlockEncode, PCPhase
from pennylane.operation import AnyWires, Operation
from pennylane import numpy as np

[docs]def qsvt(A, angles, wires, convention=None):
r"""Implements the
quantum singular value transformation <https://arxiv.org/abs/1806.01838>__ (QSVT) circuit.

.. note ::

:class:~.BlockEncode and :class:~.PCPhase used in this implementation of QSVT
are matrix-based operators and well-suited for simulators.
To implement QSVT with user-defined circuits for the block encoding and
projector-controlled phase shifts, use the :class:~.QSVT template.

Given a matrix :math:A, and a list of angles :math:\vec{\phi}, this function applies a
circuit for the quantum singular value transformation using :class:~.BlockEncode and
:class:~.PCPhase.

When the number of angles is even (:math:d is odd), the QSVT circuit is defined as:

.. math::

U_{QSVT} = \tilde{\Pi}_{\phi_1}U\left[\prod^{(d-1)/2}_{k=1}\Pi_{\phi_{2k}}U^\dagger
\tilde{\Pi}_{\phi_{2k+1}}U\right]\Pi_{\phi_{d+1}},

and when the number of angles is odd (:math:d is even):

.. math::

U_{QSVT} = \left[\prod^{d/2}_{k=1}\Pi_{\phi_{2k-1}}U^\dagger\tilde{\Pi}_{\phi_{2k}}U\right]
\Pi_{\phi_{d+1}}.

Here, :math:U denotes a block encoding of :math:A via :class:~.BlockEncode and
:math:\Pi_\phi denotes a projector-controlled phase shift with angle :math:\phi
via :class:~.PCPhase.

This circuit applies a polynomial transformation (:math:Poly^{SV}) to the singular values of
the block encoded matrix:

.. math::

\begin{align}
U_{QSVT}(A, \phi) &=
\begin{bmatrix}
Poly^{SV}(A) & \cdot \\
\cdot & \cdot
\end{bmatrix}.
\end{align}

The polynomial transformation is determined by a combination of the block encoding and choice of angles,
:math:\vec{\phi}. The convention used by :class:~.BlockEncode is commonly refered to as the
reflection convention or :math:R convention. Another equivalent convention for the block encoding is
the :math:Wx or rotation convention.

Depending on the choice of convention for blockencoding, the same phase angles will produce different
polynomial transformations. We provide the functionality to swap between blockencoding conventions and
to transform the phase angles accordingly using the :code:convention keyword argument.

.. seealso::

:class:~.QSVT and A Grand Unification of Quantum Algorithms <https://arxiv.org/pdf/2105.02859.pdf>_.

Args:
A (tensor_like): the general :math:(n \times m) matrix to be encoded
angles (tensor_like): a list of angles by which to shift to obtain the desired polynomial
wires (Iterable[int, str], Wires): the wires the template acts on
convention (string): can be set to "Wx" to convert quantum signal processing angles in the
Wx convention to QSVT angles.

**Example**

To implement QSVT in a circuit, we can use the following method:

>>> dev = qml.device("default.qubit", wires=2)
>>> A = np.array([[0.1, 0.2], [0.3, 0.4]])
>>> angles = np.array([0.1, 0.2, 0.3])
>>> @qml.qnode(dev)
... def example_circuit(A):
...     qml.qsvt(A, angles, wires=[0, 1])
...     return qml.expval(qml.PauliZ(wires=0))

The resulting circuit implements QSVT.

>>> print(qml.draw(example_circuit)(A))
0: ─╭QSVT─┤  <Z>
1: ─╰QSVT─┤

To see the implementation details, we can expand the circuit:

>>> q_script = qml.tape.QuantumScript(ops=[qml.qsvt(A, angles, wires=[0, 1])])
>>> print(q_script.expand().draw(decimals=2))
0: ─╭∏_ϕ(0.30)─╭BlockEncode(M0)─╭∏_ϕ(0.20)─╭BlockEncode(M0)†─╭∏_ϕ(0.10)─┤
1: ─╰∏_ϕ(0.30)─╰BlockEncode(M0)─╰∏_ϕ(0.20)─╰BlockEncode(M0)†─╰∏_ϕ(0.10)─┤
"""
if qml.math.shape(A) == () or qml.math.shape(A) == (1,):
A = qml.math.reshape(A, [1, 1])

c, r = qml.math.shape(A)

with qml.QueuingManager.stop_recording():
UA = BlockEncode(A, wires=wires)
projectors = []

if convention == "Wx":
angles = _qsp_to_qsvt(angles)
global_phase = (len(angles) - 1) % 4

if global_phase:
with qml.QueuingManager.stop_recording():
global_phase_op = qml.exp(
qml.Identity(wires=wires), 0.5j * np.pi * (4 - global_phase)
)

for idx, phi in enumerate(angles):
dim = c if idx % 2 else r
with qml.QueuingManager.stop_recording():
projectors.append(PCPhase(phi, dim=dim, wires=wires))

projectors = projectors[::-1]  # reverse order to match equation

if convention == "Wx":
return qml.prod(global_phase_op, QSVT(UA, projectors))
return QSVT(UA, projectors)

[docs]class QSVT(Operation):
r"""QSVT(UA,projectors)
Implements the
quantum singular value transformation <https://arxiv.org/abs/1806.01838>__ (QSVT) circuit.

.. note ::

This template allows users to define hardware-compatible block encoding and
projector-controlled phase shift circuits. For a QSVT implementation that is
tailored for simulators see :func:~.qsvt .

Given an :class:~.Operator :math:U, which block encodes the matrix :math:A, and a list of
projector-controlled phase shift operations :math:\vec{\Pi}_\phi, this template applies a
circuit for the quantum singular value transformation as follows.

When the number of projector-controlled phase shifts is even (:math:d is odd), the QSVT
circuit is defined as:

.. math::

U_{QSVT} = \tilde{\Pi}_{\phi_1}U\left[\prod^{(d-1)/2}_{k=1}\Pi_{\phi_{2k}}U^\dagger
\tilde{\Pi}_{\phi_{2k+1}}U\right]\Pi_{\phi_{d+1}}.

And when the number of projector-controlled phase shifts is odd (:math:d is even):

.. math::

U_{QSVT} = \left[\prod^{d/2}_{k=1}\Pi_{\phi_{2k-1}}U^\dagger\tilde{\Pi}_{\phi_{2k}}U\right]
\Pi_{\phi_{d+1}}.

This circuit applies a polynomial transformation (:math:Poly^{SV}) to the singular values of
the block encoded matrix:

.. math::

\begin{align}
U_{QSVT}(A, \vec{\phi}) &=
\begin{bmatrix}
Poly^{SV}(A) & \cdot \\
\cdot & \cdot
\end{bmatrix}.
\end{align}

.. seealso::

:func:~.qsvt and A Grand Unification of Quantum Algorithms <https://arxiv.org/pdf/2105.02859.pdf>_.

Args:
UA (Operator): the block encoding circuit, specified as an :class:~.Operator,
like :class:~.BlockEncode
projectors (Sequence[Operator]): a list of projector-controlled phase
shifts that implement the desired polynomial

Raises:
ValueError: if the input block encoding is not an operator

**Example**

To implement QSVT in a circuit, we can use the following method:

>>> dev = qml.device("default.qubit", wires=2)
>>> A = np.array([[0.1]])
>>> block_encode = qml.BlockEncode(A, wires=[0, 1])
>>> shifts = [qml.PCPhase(i + 0.1, dim=1, wires=[0, 1]) for i in range(3)]
>>> @qml.qnode(dev)
>>> def example_circuit():
...    qml.QSVT(block_encode, shifts)
...    return qml.expval(qml.PauliZ(wires=0))

The resulting circuit implements QSVT.

>>> print(qml.draw(example_circuit)())
0: ─╭QSVT─┤  <Z>
1: ─╰QSVT─┤

To see the implementation details, we can expand the circuit:

>>> q_script = qml.tape.QuantumScript(ops=[qml.QSVT(block_encode, shifts)])
>>> print(q_script.expand().draw(decimals=2))
0: ─╭∏_ϕ(0.10)─╭BlockEncode(M0)─╭∏_ϕ(1.10)─╭BlockEncode(M0)†─╭∏_ϕ(2.10)─┤
1: ─╰∏_ϕ(0.10)─╰BlockEncode(M0)─╰∏_ϕ(1.10)─╰BlockEncode(M0)†─╰∏_ϕ(2.10)─┤

When working with this class directly, we can make use of any PennyLane operation
to represent our block-encoding and our phase-shifts.

>>> dev = qml.device("default.qubit", wires=)
>>> block_encoding = qml.Hadamard(wires=0)  # note H is a block encoding of 1/sqrt(2)
>>> phase_shifts = [qml.RZ(-2 * theta, wires=0) for theta in (1.23, -0.5, 4)]  # -2*theta to match convention
>>>
>>> @qml.qnode(dev)
>>> def example_circuit():
...     qml.QSVT(block_encoding, phase_shifts)
...     return qml.expval(qml.PauliZ(wires=0))
>>>
>>> example_circuit()

Once again, we can visualize the circuit as follows:

>>> print(qml.draw(example_circuit)())
0: ──QSVT─┤  <Z>

To see the implementation details, we can expand the circuit:

>>> q_script = qml.tape.QuantumScript(ops=[qml.QSVT(block_encoding, phase_shifts)])
>>> print(q_script.expand().draw(decimals=2))
0: ──RZ(-2.46)──H──RZ(1.00)──H†──RZ(-8.00)─┤
"""

num_wires = AnyWires
"""int: Number of wires that the operator acts on."""

def _flatten(self):
data = (self.hyperparameters["UA"], self.hyperparameters["projectors"])
return data, tuple()

@classmethod
def _unflatten(cls, data, _) -> "QSVT":
return cls(*data)

def __init__(self, UA, projectors, id=None):
if not isinstance(UA, qml.operation.Operator):
raise ValueError("Input block encoding must be an Operator")

self._hyperparameters = {
"UA": UA,
"projectors": projectors,
}

ua_wires = UA.wires.toset()
proj_wires = set.union(*(proj.wires.toset() for proj in projectors))

total_wires = ua_wires.union(proj_wires)
super().__init__(wires=total_wires, id=id)

[docs]    @staticmethod
def compute_decomposition(
UA, projectors, **hyperparameters
):  # pylint: disable=arguments-differ
r"""Representation of the operator as a product of other operators.

The :class:~.QSVT is decomposed into alternating block encoding
and projector-controlled phase shift operators. This is defined by the following
equations, where :math:U is the block encoding operation and both :math:\Pi_\phi and
:math:\tilde{\Pi}_\phi are projector-controlled phase shifts with angle :math:\phi.

When the number of projector-controlled phase shifts is even (:math:d is odd), the QSVT
circuit is defined as:

.. math::

U_{QSVT} = \Pi_{\phi_1}U\left[\prod^{(d-1)/2}_{k=1}\Pi_{\phi_{2k}}U^\dagger
\tilde{\Pi}_{\phi_{2k+1}}U\right]\Pi_{\phi_{d+1}}.

And when the number of projector-controlled phase shifts is odd (:math:d is even):

.. math::

U_{QSVT} = \left[\prod^{d/2}_{k=1}\Pi_{\phi_{2k-1}}U^\dagger\tilde{\Pi}_{\phi_{2k}}U\right]
\Pi_{\phi_{d+1}}.

.. seealso:: :meth:~.QSVT.decomposition.

Args:
UA (Operator): the block encoding circuit, specified as a :class:~.Operator
projectors (list[Operator]): a list of projector-controlled phase
shift circuits that implement the desired polynomial
wires (Iterable): wires that the template acts on

Returns:
list[.Operator]: decomposition of the operator
"""

op_list = []

for idx, op in enumerate(projectors[:-1]):
if qml.QueuingManager.recording():
qml.apply(op)
op_list.append(op)

if idx % 2 == 0:
if qml.QueuingManager.recording():
qml.apply(UA)
op_list.append(UA)

else:

if qml.QueuingManager.recording():
qml.apply(projectors[-1])
op_list.append(projectors[-1])

return op_list

[docs]    def label(self, decimals=None, base_label=None, cache=None):
op_label = base_label or self.__class__.__name__
return op_label

[docs]    def queue(self, context=QueuingManager):
context.remove(self._hyperparameters["UA"])
for op in self._hyperparameters["projectors"]:
context.remove(op)
context.append(self)
return self

[docs]    @staticmethod
def compute_matrix(*args, **kwargs):
r"""Representation of the operator as a canonical matrix in the computational basis (static method).

The canonical matrix is the textbook matrix representation that does not consider wires.
Implicitly, this assumes that the wires of the operator correspond to the global wire order.

.. seealso:: :meth:~.Operator.matrix and :func:~.matrix

Args:
*params (list): trainable parameters of the operator, as stored in the parameters attribute
**hyperparams (dict): non-trainable hyperparameters of the operator, as stored in the hyperparameters attribute

Returns:
tensor_like: matrix representation
"""
# pylint: disable=unused-argument
op_list = []
UA = kwargs["UA"]
projectors = kwargs["projectors"]

with QueuingManager.stop_recording():  # incase this method is called in a queue context, this prevents
UA_copy = copy.copy(UA)  # us from queuing operators unnecessarily

for idx, op in enumerate(projectors[:-1]):
op_list.append(op)
if idx % 2 == 0:
op_list.append(UA)
else:

op_list.append(projectors[-1])
mat = qml.matrix(qml.prod(*tuple(op_list[::-1])))

return mat

def _qsp_to_qsvt(angles):
r"""Converts qsp angles to qsvt angles."""
new_angles = qml.math.array(copy.copy(angles))
new_angles += 3 * np.pi / 4
new_angles[-1] -= np.pi / 4

new_angles[1:-1] += np.pi / 2
return new_angles


Using PennyLane

Development

API

Internals