Source code for qiskit_algorithms.gradients.reverse.reverse_qgt

# This code is part of a Qiskit project.
# (C) Copyright IBM 2023.
# This code is licensed under the Apache License, Version 2.0. You may
# obtain a copy of this license in the LICENSE.txt file in the root directory
# of this source tree or at
# Any modifications or derivative works of this code must retain this
# copyright notice, and modified files need to carry a notice indicating
# that they have been altered from the originals.

"""QGT with the classically efficient reverse mode."""

from __future__ import annotations
from import Sequence
from typing import cast, List
import logging

import numpy as np

from qiskit.circuit import QuantumCircuit, Parameter
from qiskit.quantum_info import Statevector
from qiskit.providers import Options
from qiskit.primitives import Estimator

from ..base.base_qgt import BaseQGT
from ..base.qgt_result import QGTResult
from ..utils import DerivativeType

from .split_circuits import split
from .bind import bind
from .derive_circuit import derive_circuit

logger = logging.getLogger(__name__)

[docs]class ReverseQGT(BaseQGT): """QGT calculation with the classically efficient reverse mode. .. note:: This QGT implementation is based on statevector manipulations and scales exponentially with the number of qubits. However, for small system sizes it can be very fast compared to circuit-based gradients. This class implements the calculation of the QGT as described in [1]. By keeping track of three statevectors and iteratively sweeping through each parameterized gate, this method scales only quadratically with the number of parameters. **References:** [1]: Jones, T. "Efficient classical calculation of the Quantum Natural Gradient" (2020). `arXiv:2011.02991 <>`_. """ SUPPORTED_GATES = ["rx", "ry", "rz", "cp", "crx", "cry", "crz"] def __init__( self, phase_fix: bool = True, derivative_type: DerivativeType = DerivativeType.COMPLEX ): """ Args: phase_fix: Whether or not to include the phase fix. derivative_type: Determines whether the complex QGT or only the real or imaginary parts are calculated. """ dummy_estimator = Estimator() # this method does not need an estimator super().__init__(dummy_estimator, phase_fix, derivative_type) @property def options(self) -> Options: """There are no options for the reverse QGT, returns an empty options dict. Returns: Empty options. """ return Options() def _run( # pylint: disable=arguments-renamed self, circuits: Sequence[QuantumCircuit], parameter_values: Sequence[Sequence[float]], parameters: Sequence[Sequence[Parameter]], **options, ) -> QGTResult: """Compute the QGT on the given circuits.""" g_circuits, g_parameter_values, g_parameter_sets = self._preprocess( circuits, parameter_values, parameters, self.SUPPORTED_GATES ) results = self._run_unique(g_circuits, g_parameter_values, g_parameter_sets, **options) return self._postprocess(results, circuits, parameter_values, parameters) def _run_unique( self, circuits: Sequence[QuantumCircuit], parameter_values: Sequence[Sequence[float]], parameter_sets: Sequence[Sequence[Parameter]], **options, # pylint: disable=unused-argument ) -> QGTResult: num_qgts = len(circuits) qgts = [] metadata = [] for k in range(num_qgts): values = np.asarray(parameter_values[k]) circuit = circuits[k] parameters = list(parameter_sets[k]) num_parameters = len(parameters) original_parameter_order = [p for p in circuit.parameters if p in parameters] metadata.append({"parameters": original_parameter_order}) unitaries, paramlist = split(circuit, parameters=parameters) # initialize the phase fix vector and the hessian part ``metric`` num_parameters = len(unitaries) phase_fixes = np.zeros(num_parameters, dtype=complex) metric = np.zeros((num_parameters, num_parameters), dtype=complex) # initialize the state variables -- naming convention is the same as the paper parameter_binds = dict(zip(circuit.parameters, values)) bound_unitaries = cast(List[QuantumCircuit], bind(unitaries, parameter_binds)) chi = Statevector(bound_unitaries[0]) psi = chi.copy() phi = Statevector.from_int(0, (2,) * circuit.num_qubits) # Get the analytic gradient of the first unitary # Note: We currently only support gates with a single parameter -- which is reflected # in self.SUPPORTED_GATES -- but generally we could also support gates with multiple # parameters per gate. This is the reason for the second 0-index. deriv = derive_circuit(unitaries[0], paramlist[0][0]) for _, gate in deriv: bind(gate, parameter_binds, inplace=True) grad_coeffs = [coeff for coeff, _ in deriv] grad_states = [phi.evolve(gate) for _, gate in deriv] # compute phase fix (optional) and the hessian part if self._phase_fix: phase_fixes[0] = _phasefix_term(chi, grad_coeffs, grad_states) metric[0, 0] = _l_term(grad_coeffs, grad_states, grad_coeffs, grad_states) for j in range(1, num_parameters): lam = psi.copy() phi = psi.copy() # get the analytic gradient d U_j / d p_j and apply it deriv = derive_circuit(unitaries[j], paramlist[j][0]) for _, gate in deriv: bind(gate, parameter_binds, inplace=True) # compute |phi> (in general it's a sum of states and coeffs) grad_coeffs = [coeff for coeff, _ in deriv] grad_states = [phi.evolve(gate) for _, gate in deriv] # compute the diagonal element L_{j, j} metric[j, j] += _l_term(grad_coeffs, grad_states, grad_coeffs, grad_states) # compute the off diagonal elements L_{i, j} for i in reversed(range(j)): # apply U_{i + 1}_dg unitary_ip_inv = bound_unitaries[i + 1].inverse() grad_states = [state.evolve(unitary_ip_inv) for state in grad_states] lam = lam.evolve(bound_unitaries[i].inverse()) # get the gradient d U_i / d p_i and apply it deriv = derive_circuit(unitaries[i], paramlist[i][0]) for _, gate in deriv: bind(gate, parameter_binds, inplace=True) grad_coeffs_mu = [coeff for coeff, _ in deriv] grad_states_mu = [lam.evolve(gate) for _, gate in deriv] metric[i, j] += _l_term( grad_coeffs_mu, grad_states_mu, grad_coeffs, grad_states ) if self._phase_fix: phase_fixes[j] += _phasefix_term(chi, grad_coeffs, grad_states) psi = psi.evolve(bound_unitaries[j]) # The following code stacks the QGT together and maps the values into the # correct original order of parameters # map circuit parameter to global index in the circuit param_to_circuit = { param: index for index, param in enumerate(original_parameter_order) } # map global index to the local index used in the calculation, the new index can # now be accessed by remap[index] remap = { index: param_to_circuit[_extract_parameter(plist[0])] for index, plist in enumerate(paramlist) } qgt = np.zeros((num_parameters, num_parameters), dtype=complex) for i in range(num_parameters): iloc = remap[i] for j in range(num_parameters): jloc = remap[j] if i <= j: qgt[iloc, jloc] += metric[i, j] else: qgt[iloc, jloc] += np.conj(metric[j, i]) qgt[iloc, jloc] -= np.conj(phase_fixes[i]) * phase_fixes[j] # append and cast to real/imag if required qgts.append(self._to_derivtype(qgt)) result = QGTResult(qgts, self.derivative_type, metadata, options=None) return result def _to_derivtype(self, qgt): if self.derivative_type == DerivativeType.REAL: return np.real(qgt) if self.derivative_type == DerivativeType.IMAG: return np.imag(qgt) return qgt
def _l_term(coeffs_i, states_i, coeffs_j, states_j): return sum( sum( np.conj(coeff_i) * coeff_j * np.conj( for coeff_i, state_i in zip(coeffs_i, states_i) ) for coeff_j, state_j in zip(coeffs_j, states_j) ) def _phasefix_term(chi, coeffs, states): return sum( coeff_i * np.conj( for coeff_i, state_i in zip(coeffs, states) ) def _extract_parameter(expression): if isinstance(expression, Parameter): return expression if len(expression.parameters) > 1: raise ValueError("Expression has more than one parameter.") return list(expression.parameters)[0]