Bosonic operators

.pylintdict

@@ -58,7 +58,9 @@ boolean
 bools
 bopes
 boson
+bosons
 bosonic
+bosonicop
 bravyi
 calcinfo
 cargs
@@ -336,6 +338,7 @@ metadata
 meth
 methionine
 microsoft
+miessen
 minao
 mixin
 miyazawa
@@ -521,9 +524,11 @@ sklearn
 sl
 slater
 soham
+somma
 sp
 sparsearray
 sparselabelop
+sparsepauliop
 spinop
 springer
 spsa

qiskit_nature/second_q/mappers/__init__.py

@@ -17,7 +17,8 @@
 
 .. currentmodule:: qiskit_nature.second_q.mappers
 
-The classes here are used to convert fermionic, vibrational and spin operators to qubit operators.
+The classes here are used to convert fermionic, bosonic, vibrational and spin operators to qubit
+operators.
 
 .. autosummary::
    :toctree: ../stubs/
@@ -50,6 +51,17 @@
 
    InterleavedQubitMapper
 
+
+BosonicOp Mappers
++++++++++++++++++++
+
+.. autosummary::
+   :toctree: ../stubs/
+   :template: autosummary/class_with_inherited_members.rst
+   :nosignatures:
+
+   BosonicLinearMapper
+
 VibrationalOp Mappers
 +++++++++++++++++++++
 
@@ -102,6 +114,7 @@
 from .jordan_wigner_mapper import JordanWignerMapper
 from .parity_mapper import ParityMapper
 from .linear_mapper import LinearMapper
+from .bosonic_linear_mapper import BosonicLinearMapper
 from .logarithmic_mapper import LogarithmicMapper
 from .direct_mapper import DirectMapper
 from .qubit_mapper import QubitMapper
@@ -116,6 +129,7 @@
     "JordanWignerMapper",
     "ParityMapper",
     "LinearMapper",
+    "BosonicLinearMapper",
     "LogarithmicMapper",
     "QubitConverter",
     "QubitMapper",

qiskit_nature/second_q/mappers/bosonic_linear_mapper.py

@@ -0,0 +1,195 @@
+# This code is part of Qiskit.
+#
+# (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 import settings
+
+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.
+
+        Raises:
+            NotImplementedError: when `qiskit_nature.settings.use_pauli_sum_op` is set to `True`.
+                This value is deprecated and, thus, not supported by this new implementation.
+                Set it to `False` in order to use this mapper.
+        """
+        if settings.use_pauli_sum_op:
+            raise NotImplementedError(
+                "This class does not support the deprecated `settings.use_pauli_sum_op` value."
+                "Please set `settings.use_pauli_sum_op = False` in order to use this mapper."
+            )
+        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

qiskit_nature/second_q/mappers/bosonic_mapper.py

@@ -0,0 +1,61 @@
+# This code is part of Qiskit.
+#
+# (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.
+
+"""Bosonic Mapper."""
+
+from __future__ import annotations
+
+from qiskit.quantum_info import SparsePauliOp
+
+from qiskit_nature.second_q.operators import BosonicOp
+
+from .qubit_mapper import ListOrDictType, QubitMapper
+
+
+class BosonicMapper(QubitMapper):
+    """Mapper of Bosonic Operator to Qubit Operator
+
+    The following attributes can be read and updated once the ``BosonicMapper`` object
+    has been constructed.
+
+    """
+
+    def __init__(self, max_occupation: int) -> None:
+        """
+        Args:
+            max_occupation: defines the excitation space of the k-th bosonic state. Together with the
+                number of modes required to represent the bosonic operator, it defines the minimum length
+                of the qubit register. The minimum value is 1.
+        """
+        super().__init__()
+        self.max_occupation = max_occupation
+
+    @property
+    def max_occupation(self) -> int:
+        """The maximum occupation of any bosonic state."""
+        return self._max_occupation
+
+    @max_occupation.setter
+    def max_occupation(self, max_occupation: int) -> None:
+        if max_occupation < 1:
+            raise ValueError(
+                f"The maximum occupation must be at least 1, and not {max_occupation}."
+            )
+        self._max_occupation = max_occupation
+
+    def map(
+        self,
+        second_q_ops: BosonicOp | ListOrDictType[BosonicOp],
+        *,
+        register_length: int | None = None,
+    ) -> SparsePauliOp | ListOrDictType[SparsePauliOp]:
+        return super().map(second_q_ops, register_length=register_length)

qiskit_nature/second_q/operators/__init__.py

@@ -23,6 +23,7 @@
 
    ElectronicIntegrals
    FermionicOp
+   BosonicOp
    SparseLabelOp
    SpinOp
    VibrationalOp
@@ -43,6 +44,7 @@
 
 from .electronic_integrals import ElectronicIntegrals
 from .fermionic_op import FermionicOp
+from .bosonic_op import BosonicOp
 from .spin_op import SpinOp
 from .vibrational_op import VibrationalOp
 from .vibrational_integrals import VibrationalIntegrals
@@ -53,6 +55,7 @@
 __all__ = [
     "ElectronicIntegrals",
     "FermionicOp",
+    "BosonicOp",
     "SpinOp",
     "VibrationalOp",
     "VibrationalIntegrals",