Replace pauli expectation value cython with multithreaded rust implem…

…entation (#7702)

* Replace pauli expectation value cython with rust implementation

This commit replaces the cython implementation of the pauli expectation
value functions with a multithreaded rust implementation. This was done
primarily for two reasons, the first and primary reason for this change
is because after #7658 this module was the only cython code left in the
qiskit-terra repository so unifying on a single compiled language will
reduce the maintanence burden in qiskit-terra. The second reason is
similar to the rationale in #7658 around why using rust over cython for
multi-threaded hybrid python module. The difference here though is
unlike in stochastic swap this module isn't as performance critical as
it's not nearly as widely used.

* Tune single threaded performance for rust sum

This commit tunes the sum for the single threaded path. Using the
iterator sum() method is very convienent but for the single threaded
path it doesn't create the most efficient output. This was causing a
regression in performance over the previous cython version. To address
that issue, this commit adds a new tuned function which does a chunked
sum which the compiler can handle better. It more closely models how
we'd do this with vectorized SIMD instructions. As a future step we can
look at using simdeez /~https://github.com/jackmott/simdeez
to further optimize this by doing runtime CPU feature detection and
leveraging SIMD instrinsics (we might want to look at using `fast_sum()`
in the multithreaded path if we do that too).

* Add release notes

* Fix lint

* Add docstring and signature to rust functions

* Define parallel threshold as a constant

* Add attribution comment to fast_sum()

* Rename eval_parallel_env -> getenv_use_multiple_threads

* Use inline literal type for size

Co-authored-by: Kevin Hartman <kevin@hart.mn>

* Add overflow check on num_qubits

The functions only work for at most for number of qubits < usize bits
anything larger would cause an overflow. While rust provides overflow
checking in debug mode it disables this for performance in release mode.
Sice we ship binaries in release mode this commit adds an overflow check
for the num_qubits argument to ensure that we don't overflow and produce
incorrect results.

* Remove unecessary setup_requires field from setup.py

The setup_requires field in the setup.py is deprecated and
has been superseded by the pyproject.toml to define build
system dependencies. Since we're already relying on the
pyproject.toml to install setuptools-rust for us having the
setup_requires line will do nothing but potentially cause
issues as it will use an older install mechanism that will potentially conflict with  people's environments.

* Drop `.iter().take(LANES)`.

* Fix typo.

Co-authored-by: Kevin Hartman <kevin@hart.mn>

Cargo.toml

@@ -15,10 +15,11 @@ numpy = "0.16.1"
 rand = "0.8"
 rand_pcg = "0.3"
 rand_distr = "0.4.3"
+num-complex = "0.4"
 
 [dependencies.pyo3]
 version = "0.16.1"
-features = ["extension-module", "hashbrown"]
+features = ["extension-module", "hashbrown", "num-complex"]
 
 [dependencies.ndarray]
 version = "^0.15.0"

qiskit/__init__.py

@@ -25,7 +25,7 @@
 # manually define them on import so people can directly import
 # qiskit._accelerate.* submodules and not have to rely on attribute access
 sys.modules["qiskit._accelerate.stochastic_swap"] = qiskit._accelerate.stochastic_swap
-
+sys.modules["qiskit._accelerate.pauli_expval"] = qiskit._accelerate.pauli_expval
 
 # qiskit errors operator
 from qiskit.exceptions import QiskitError, MissingOptionalLibraryError

qiskit/quantum_info/states/densitymatrix.py

@@ -32,8 +32,8 @@
 from qiskit.quantum_info.operators.channel.quantum_channel import QuantumChannel
 from qiskit.quantum_info.operators.channel.superop import SuperOp
 
-# pylint: disable=no-name-in-module
-from .cython.exp_value import density_expval_pauli_no_x, density_expval_pauli_with_x
+# pylint: disable=import-error
+from qiskit._accelerate.pauli_expval import density_expval_pauli_no_x, density_expval_pauli_with_x
 
 
 class DensityMatrix(QuantumState, TolerancesMixin):

qiskit/quantum_info/states/statevector.py

@@ -31,8 +31,11 @@
 from qiskit.quantum_info.operators.op_shape import OpShape
 from qiskit.quantum_info.operators.predicates import matrix_equal
 
-# pylint: disable=no-name-in-module
-from .cython.exp_value import expval_pauli_no_x, expval_pauli_with_x
+# pylint: disable=import-error
+from qiskit._accelerate.pauli_expval import (
+    expval_pauli_no_x,
+    expval_pauli_with_x,
+)
 
 
 class Statevector(QuantumState, TolerancesMixin):

releasenotes/notes/rust-pauli-expval-f2aa06c5bab85768.yaml

@@ -0,0 +1,16 @@
+---
+features:
+  - |
+    The internal computations of :class:`.Statevector.expectation_value` and
+    :class:`.DensityMatrix.expectation_value methods have been reimplemented
+    in the Rust programming language. This new implementation is multithreaded
+    and by default for a :class:`~.Statevector` or :class:`~.DensityMatrix`
+    >= 19 qubits will spawn a thread pool with the number of logical CPUs
+    available on the local system. You can you can control the number of
+    threads used by setting the ``RAYON_NUM_THREADS`` environment variable to
+    an integer value. For example, setting ``RAYON_NUM_THREADS=4`` will only
+    use 4 threads in the thread pool.
+upgrade:
+  - |
+    Cython is no longer a build dependency of Qiskit Terra and is no longer
+    required to be installed when building Qiskit Terra from source.

setup.py

@@ -18,51 +18,10 @@
 from setuptools import setup, find_packages, Extension
 from setuptools_rust import Binding, RustExtension
 
-try:
-    from Cython.Build import cythonize
-except ImportError:
-    import subprocess
-
-    subprocess.call([sys.executable, "-m", "pip", "install", "Cython>=0.27.1"])
-    from Cython.Build import cythonize
-
 
 with open("requirements.txt") as f:
     REQUIREMENTS = f.read().splitlines()
 
-# Add Cython extensions here
-CYTHON_EXTS = {
-    "qiskit/quantum_info/states/cython/exp_value": "qiskit.quantum_info.states.cython.exp_value",
-}
-
-INCLUDE_DIRS = []
-# Extra link args
-LINK_FLAGS = []
-# If on Win and not in MSYS2 (i.e. Visual studio compile)
-if sys.platform == "win32" and os.environ.get("MSYSTEM") is None:
-    COMPILER_FLAGS = ["/O2"]
-# Everything else
-else:
-    COMPILER_FLAGS = ["-O2", "-funroll-loops", "-std=c++11"]
-    if sys.platform == "darwin":
-        # These are needed for compiling on OSX 10.14+
-        COMPILER_FLAGS.append("-mmacosx-version-min=10.9")
-        LINK_FLAGS.append("-mmacosx-version-min=10.9")
-
-
-EXT_MODULES = []
-# Add Cython Extensions
-for src, module in CYTHON_EXTS.items():
-    ext = Extension(
-        module,
-        sources=[src + ".pyx"],
-        include_dirs=INCLUDE_DIRS,
-        extra_compile_args=COMPILER_FLAGS,
-        extra_link_args=LINK_FLAGS,
-        language="c++",
-    )
-    EXT_MODULES.append(ext)
-
 # Read long description from README.
 README_PATH = os.path.join(os.path.abspath(os.path.dirname(__file__)), "README.md")
 with open(README_PATH) as readme_file:
@@ -117,7 +76,6 @@
     keywords="qiskit sdk quantum",
     packages=find_packages(exclude=["test*"]),
     install_requires=REQUIREMENTS,
-    setup_requires=["Cython>=0.27.1"],
     include_package_data=True,
     python_requires=">=3.7",
     extras_require={
@@ -134,7 +92,6 @@
         "Documentation": "https://qiskit.org/documentation/",
         "Source Code": "/~https://github.com/Qiskit/qiskit-terra",
     },
-    ext_modules=cythonize(EXT_MODULES),
     rust_extensions=[RustExtension("qiskit._accelerate", "Cargo.toml", binding=Binding.PyO3)],
     zip_safe=False,
     entry_points={

src/lib.rs

@@ -10,16 +10,33 @@
 // copyright notice, and modified files need to carry a notice indicating
 // that they have been altered from the originals.
 
+use std::env;
+
 use pyo3::prelude::*;
 use pyo3::wrap_pymodule;
 use pyo3::Python;
 
 mod edge_collections;
 mod nlayout;
+mod pauli_exp_val;
 mod stochastic_swap;
 
+#[inline]
+pub fn getenv_use_multiple_threads() -> bool {
+    let parallel_context = env::var("QISKIT_IN_PARALLEL")
+        .unwrap_or_else(|_| "FALSE".to_string())
+        .to_uppercase()
+        == "TRUE";
+    let force_threads = env::var("QISKIT_FORCE_THREADS")
+        .unwrap_or_else(|_| "FALSE".to_string())
+        .to_uppercase()
+        == "TRUE";
+    !parallel_context || force_threads
+}
+
 #[pymodule]
 fn _accelerate(_py: Python<'_>, m: &PyModule) -> PyResult<()> {
     m.add_wrapped(wrap_pymodule!(stochastic_swap::stochastic_swap))?;
+    m.add_wrapped(wrap_pymodule!(pauli_exp_val::pauli_expval))?;
     Ok(())
 }