qiskit_nature.second_q.mappers.bosonic_linear_mapper のソースコード

# 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 http://www.apache.org/licenses/LICENSE-2.0.
# 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.

"""The Linear Mapper for Bosons."""

from __future__ import annotations
import operator

from functools import reduce, lru_cache

import numpy as np

from qiskit.quantum_info import Pauli, SparsePauliOp

from qiskit_nature.second_q.operators import BosonicOp
from .bosonic_mapper import BosonicMapper

[ドキュメント]class BosonicLinearMapper(BosonicMapper): """The Linear boson-to-qubit mapping. This mapper generates a linear encoding of the Bosonic operator :math:`b_k^\\dagger, b_k` to qubit operators (linear combinations of pauli strings). In this linear encoding each bosonic mode is represented via :math:`n_k^{max} + 1` qubits, where :math:`n_k^{max}` is the max occupation of the mode (meaning the number of states used in the expansion of the mode, or equivalently the state at which the maximum excitation can take place). The mode :math:`|k\\rangle` is then mapped to the occupation number vector :math:`|0_{n_k^{max}}, 0_{n_k^{max} - 1},..., 0_{n_k + 1}, 1_{n_k}, 0_{n_k - 1},..., 0_{0_k}\\rangle` It implements the formula in Section II.C of Reference [1]: .. math:: b_k^\\dagger = \\sum_{n_k =0}^{n_k^{max}-1}(\\sqrt{n_k +1}\\sigma_{n_k}^{+}\\sigma_{n_k + 1}^{-}) from :math:`n_k = 0` to :math:`n_k^{max} + 1` where :math:`n_k^{max}` is the maximum occupation (defined by the user). In the following implementation, we explicit the operators :math:`\\sigma^+` and :math:`\\sigma^-` with the Pauli matrices: .. math:: \\sigma_{n_k}^+ := S_j^+ = 0.5 * (X_j + \\textit{i}Y_j) \\sigma_{n_k}^- := S_j^- = 0.5 * (X_j - \\textit{i}Y_j) The length of the qubit register is: .. code-block:: python BosonicOp.num_modes * (BosonicLinearMapper.max_occupation + 1) To use this mapper one can for example: .. code-block:: python from qiskit_nature.second_q.mappers import BosonicLinearMapper from qiskit_nature.second_q.operators import BosonicOp mapper = BosonicLinearMapper(max_occupation=1) qubit_op = mapper.map(BosonicOp({'+_0 -_0': 1}, num_modes=1)) .. note:: Since this mapper truncates the maximum occupation of a bosonic state as represented in the qubit register, the commutation relation after the mapping differ from the standard ones. Please refer to Section 4, equation 22 of Reference [2] for more details References: [1] A. Miessen et al., Quantum algorithms for quantum dynamics: A performance study on the spin-boson model, Phys. Rev. Research 3, 043212. https://link.aps.org/doi/10.1103/PhysRevResearch.3.043212 [2] R. Somma et al., Quantum Simulations of Physics Problems, Arxiv https://doi.org/10.48550/arXiv.quant-ph/0304063 """ def _map_single( self, second_q_op: BosonicOp, *, register_length: int | None = None ) -> SparsePauliOp: """Maps a :class:`~qiskit_nature.second_q.operators.SparseLabelOp` to a``SparsePauliOp``. Args: second_q_op: the ``SparseLabelOp`` to be mapped. register_length: when provided, this will be used to overwrite the ``register_length`` attribute of the operator being mapped. This is possible because the ``register_length`` is considered a lower bound in a ``SparseLabelOp``. Returns: The qubit operator corresponding to the problem-Hamiltonian in the qubit space. """ if register_length is None: register_length = second_q_op.num_modes qubit_register_length = register_length * (self.max_occupation + 1) # Create a Pauli operator, which we will fill in this method pauli_op: list[SparsePauliOp] = [] # Then we loop over all the terms of the bosonic operator for terms, coeff in second_q_op.terms(): # Then loop over each term (terms -> List[Tuple[string, int]]) bos_op_to_pauli_op = SparsePauliOp(["I" * qubit_register_length], coeffs=[1.0]) for op, idx in terms: if op not in ("+", "-"): break pauli_expansion: list[SparsePauliOp] = [] # Now we are dealing with a single bosonic operator. We have to perform the linear mapper for n_k in range(self.max_occupation): prefactor = np.sqrt(n_k + 1) / 4.0 # Define the actual index in the qubit register. It is given by n_k plus the shift # due to the mode onto which the operator is acting register_index = n_k + idx * (self.max_occupation + 1) # Now build the Pauli operators XX, XY, YX, YY, which arise from S_i^+ S_j^- x_x, x_y, y_x, y_y = self._get_ij_pauli_matrix( register_index, qubit_register_length ) tmp_op = SparsePauliOp(x_x) + SparsePauliOp(y_y) if op == "+": tmp_op += -1j * SparsePauliOp(x_y) + 1j * SparsePauliOp(y_x) else: tmp_op += +1j * SparsePauliOp(x_y) - 1j * SparsePauliOp(y_x) pauli_expansion.append(prefactor * tmp_op) # Add the Pauli expansion for a single n_k to map of the bosonic operator bos_op_to_pauli_op = reduce(operator.add, pauli_expansion).compose( bos_op_to_pauli_op ) # Add the map of the single boson op (e.g. +_0) to the map of the full bosonic operator pauli_op.append(coeff * reduce(operator.add, bos_op_to_pauli_op.simplify())) # return the lookup table for the transformed XYZI operators bos_op_encoding = reduce(operator.add, pauli_op) return bos_op_encoding @classmethod @lru_cache(maxsize=32) def _get_ij_pauli_matrix( cls, register_index: int, register_length: int ) -> tuple[Pauli, Pauli, Pauli, Pauli]: """This method builds the Qiskit Pauli operators of the operators XX, YY, XY and YX Args: register_index: the index of the qubit register where the mapped operator should be placed. register_length: the length of the qubit register. Returns: Four Pauli operators that represent XX, XY, YX and YY at the specified index in the current qubit register. """ # Define recurrent variables prefix_zeros = [0] * register_index suffix_zeros = [0] * (register_length - 2 - register_index) # Build the Pauli strings x_x = Pauli( ( [0] * register_length, prefix_zeros + [1, 1] + suffix_zeros, ) ) x_y = Pauli( ( prefix_zeros + [0, 1] + suffix_zeros, prefix_zeros + [1, 1] + suffix_zeros, ) ) y_x = Pauli( ( prefix_zeros + [1, 0] + suffix_zeros, prefix_zeros + [1, 1] + suffix_zeros, ) ) y_y = Pauli( ( prefix_zeros + [1, 1] + suffix_zeros, prefix_zeros + [1, 1] + suffix_zeros, ) ) return x_x, x_y, y_x, y_y