From 364fc459108677a0e23b84e8aed8c77d60cb04ee Mon Sep 17 00:00:00 2001
From: Saeed Bohloul <sbohloul@gmail.com>
Date: Mon, 20 Nov 2023 18:31:11 -0500
Subject: [PATCH 01/11] Add inplace omp parallel implementaion of expval

* Add operations functors in ExpValFunctorsLQubit.hpp

* Add method to call functors in Measurement class

* Add enum class for functors

* Add mapping for operation names to functors
---
 .../measurements/ExpValFunctorsLQubit.hpp     | 210 ++++++++++++++++++
 .../measurements/MeasurementsLQubit.hpp       | 167 +++++++++++++-
 2 files changed, 376 insertions(+), 1 deletion(-)
 create mode 100644 pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp

diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
new file mode 100644
index 0000000000..93210226d1
--- /dev/null
+++ b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
@@ -0,0 +1,210 @@
+// Copyright 2018-2023 Xanadu Quantum Technologies Inc.
+
+// Licensed under the Apache License, Version 2.0 (the License);
+// you may not use this file except in compliance with the License.
+// You may obtain a copy of the License at
+
+// http://www.apache.org/licenses/LICENSE-2.0
+
+// Unless required by applicable law or agreed to in writing, software
+// distributed under the License is distributed on an AS IS BASIS,
+// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+// See the License for the specific language governing permissions and
+// limitations under the License.
+#pragma once
+
+#include "BitUtil.hpp"
+
+/// @cond DEV
+namespace {
+using namespace Pennylane::Util;
+using Pennylane::Util::INVSQRT2;
+} // namespace
+/// @endcond
+
+namespace Pennylane::LightningQubit::Functors {
+
+template <class PrecisionT> struct getExpectationValueIdentityFunctor {
+    std::complex<PrecisionT> *arr;
+    size_t num_qubits;
+
+    getExpectationValueIdentityFunctor(
+        std::complex<PrecisionT> *arr_,
+        [[maybe_unused]] std::size_t num_qubits_,
+        [[maybe_unused]] const std::vector<size_t> &wires) {
+        arr = arr_;
+        num_qubits = num_qubits_;
+    }
+
+    inline void operator()(PrecisionT &expval) const {
+        size_t k;
+#if defined(_OPENMP)
+#pragma omp parallel for default(none) shared(num_qubits, arr) private(k)      \
+    reduction(+ : expval)
+#endif
+        for (k = 0; k < exp2(num_qubits - 1); k++) {
+            expval += real(conj(arr[k] * arr[k]));
+        }
+    }
+};
+
+template <class PrecisionT> struct getExpectationValuePauliXFunctor {
+    std::complex<PrecisionT> *arr;
+    size_t num_qubits;
+
+    size_t rev_wire;
+    size_t rev_wire_shift;
+    size_t wire_parity;
+    size_t wire_parity_inv;
+
+    getExpectationValuePauliXFunctor(std::complex<PrecisionT> *arr_,
+                                     std::size_t num_qubits_,
+                                     const std::vector<size_t> &wires) {
+        arr = arr_;
+        num_qubits = num_qubits_;
+        rev_wire = num_qubits - wires[0] - 1;
+        rev_wire_shift = (static_cast<size_t>(1U) << rev_wire);
+        wire_parity = fillTrailingOnes(rev_wire);
+        wire_parity_inv = fillLeadingOnes(rev_wire + 1);
+    }
+
+    inline void operator()(PrecisionT &expval) const {
+        size_t k, i0, i1;
+#if defined(_OPENMP)
+#pragma omp parallel for default(none)                                         \
+    shared(num_qubits, wire_parity_inv, wire_parity, rev_wire_shift,           \
+               arr) private(k, i0, i1) reduction(+ : expval)
+#endif
+        for (k = 0; k < exp2(num_qubits - 1); k++) {
+            i0 = ((k << 1U) & wire_parity_inv) | (wire_parity & k);
+            i1 = i0 | rev_wire_shift;
+
+            expval +=
+                real(conj(arr[i0]) * arr[i1]) + real(conj(arr[i1]) * arr[i0]);
+        }
+    }
+};
+
+template <class PrecisionT> struct getExpectationValuePauliYFunctor {
+    std::complex<PrecisionT> *arr;
+    size_t num_qubits;
+
+    size_t rev_wire;
+    size_t rev_wire_shift;
+    size_t wire_parity;
+    size_t wire_parity_inv;
+
+    getExpectationValuePauliYFunctor(std::complex<PrecisionT> *arr_,
+                                     std::size_t num_qubits_,
+                                     const std::vector<size_t> &wires) {
+        arr = arr_;
+        num_qubits = num_qubits_;
+        rev_wire = num_qubits - wires[0] - 1;
+        rev_wire_shift = (static_cast<size_t>(1U) << rev_wire);
+        wire_parity = fillTrailingOnes(rev_wire);
+        wire_parity_inv = fillLeadingOnes(rev_wire + 1);
+    }
+
+    inline void operator()(PrecisionT &expval) const {
+        size_t k, i0, i1;
+        std::complex<PrecisionT> v0, v1;
+#if defined(_OPENMP)
+#pragma omp parallel for default(none)                                         \
+    shared(num_qubits, wire_parity_inv, wire_parity, rev_wire_shift,           \
+               arr) private(k, i0, i1, v0, v1) reduction(+ : expval)
+#endif
+        for (k = 0; k < exp2(num_qubits - 1); k++) {
+            i0 = ((k << 1U) & wire_parity_inv) | (wire_parity & k);
+            i1 = i0 | rev_wire_shift;
+            v0 = arr[i0];
+            v1 = arr[i1];
+
+            expval += real(conj(arr[i0]) *
+                           std::complex<PrecisionT>{imag(v1), -real(v1)}) +
+                      real(conj(arr[i1]) *
+                           std::complex<PrecisionT>{-imag(v0), real(v0)});
+        }
+    }
+};
+
+template <class PrecisionT> struct getExpectationValuePauliZFunctor {
+    std::complex<PrecisionT> *arr;
+    size_t num_qubits;
+
+    size_t rev_wire;
+    size_t rev_wire_shift;
+    size_t wire_parity;
+    size_t wire_parity_inv;
+
+    getExpectationValuePauliZFunctor(std::complex<PrecisionT> *arr_,
+                                     std::size_t num_qubits_,
+                                     const std::vector<size_t> &wires) {
+        arr = arr_;
+        num_qubits = num_qubits_;
+        rev_wire = num_qubits - wires[0] - 1;
+        rev_wire_shift = (static_cast<size_t>(1U) << rev_wire);
+        wire_parity = fillTrailingOnes(rev_wire);
+        wire_parity_inv = fillLeadingOnes(rev_wire + 1);
+    }
+
+    inline void operator()(PrecisionT &expval) const {
+        size_t k, i0, i1;
+#if defined(_OPENMP)
+#pragma omp parallel for default(none)                                         \
+    shared(num_qubits, wire_parity_inv, wire_parity, rev_wire_shift,           \
+               arr) private(k, i0, i1) reduction(+ : expval)
+#endif
+        for (k = 0; k < exp2(num_qubits - 1); k++) {
+
+            i0 = ((k << 1U) & wire_parity_inv) | (wire_parity & k);
+            i1 = i0 | rev_wire_shift;
+
+            expval += real(conj(arr[i1]) * (-arr[i1])) +
+                      real(conj(arr[i0]) * (arr[i0]));
+        }
+    }
+};
+
+template <class PrecisionT> struct getExpectationValueHadamardFunctor {
+    std::complex<PrecisionT> *arr;
+    size_t num_qubits;
+
+    size_t rev_wire;
+    size_t rev_wire_shift;
+    size_t wire_parity;
+    size_t wire_parity_inv;
+
+    constexpr static auto isqrt2 = INVSQRT2<PrecisionT>();
+
+    getExpectationValueHadamardFunctor(std::complex<PrecisionT> *arr_,
+                                       std::size_t num_qubits_,
+                                       const std::vector<size_t> &wires) {
+        arr = arr_;
+        num_qubits = num_qubits_;
+        rev_wire = num_qubits - wires[0] - 1;
+        rev_wire_shift = (static_cast<size_t>(1U) << rev_wire);
+        wire_parity = fillTrailingOnes(rev_wire);
+        wire_parity_inv = fillLeadingOnes(rev_wire + 1);
+    }
+
+    inline void operator()(PrecisionT &expval) const {
+        size_t k, i0, i1;
+        std::complex<PrecisionT> v0, v1;
+#if defined(_OPENMP)
+#pragma omp parallel for default(none)                                         \
+    shared(num_qubits, wire_parity_inv, wire_parity, rev_wire_shift,           \
+               arr) private(k, i0, i1, v0, v1) reduction(+ : expval)
+#endif
+
+        for (k = 0; k < exp2(num_qubits - 1); k++) {
+            i0 = ((k << 1U) & wire_parity_inv) | (wire_parity & k);
+            i1 = i0 | rev_wire_shift;
+            v0 = arr[i0];
+            v1 = arr[i1];
+
+            expval += real(isqrt2 * (conj(arr[i0]) * (v0 + v1) +
+                                     conj(arr[i1]) * (v0 - v1)));
+        }
+    }
+};
+} // namespace Pennylane::LightningQubit::Functors
\ No newline at end of file
diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
index e1d3ef551d..856a0d4033 100644
--- a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
+++ b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
@@ -29,6 +29,7 @@
 #include <unordered_map>
 #include <vector>
 
+#include "ExpValFunctorsLQubit.hpp"
 #include "LinearAlgebra.hpp"
 #include "MeasurementsBase.hpp"
 #include "Observables.hpp"
@@ -38,12 +39,26 @@
 #include "TransitionKernels.hpp"
 #include "Util.hpp" //transpose_state_tensor, sorting_indices
 
+#include "BitUtil.hpp"
+#include <iostream>
+
 /// @cond DEV
 namespace {
 using namespace Pennylane::Measures;
 using namespace Pennylane::Observables;
 using Pennylane::LightningQubit::StateVectorLQubitManaged;
 using Pennylane::LightningQubit::Util::innerProdC;
+using namespace Pennylane::LightningQubit::Functors;
+using Pennylane::Util::INVSQRT2;
+enum class ExpValFunc : uint32_t {
+    BEGIN = 1,
+    Identity = 1,
+    PauliX,
+    PauliY,
+    PauliZ,
+    Hadamard,
+    END
+};
 } // namespace
 /// @endcond
 
@@ -67,7 +82,63 @@ class Measurements final
 
   public:
     explicit Measurements(const StateVectorT &statevector)
-        : BaseType{statevector} {};
+        : BaseType{statevector} {
+        init_expval_funcs_();
+    };
+
+    /**
+     * @brief Templated method that returns the expectation value of named
+     * observables using in-place operations, without creating extra copy of
+     * the statevector.
+     *
+     * @tparam functor_t Expectation value functor class for in-place
+     * operations.
+     * @tparam nqubits Number of wires.
+     * @param wires Wires to which the observable is applied.
+     * @return expectation value of the observable.
+     */
+    template <template <class> class functor_t, size_t num_wires>
+    PrecisionT applyExpValNamedFunctor(const std::vector<size_t> &wires) {
+        if constexpr (num_wires > 0) {
+            PL_ASSERT(wires.size() == num_wires);
+        }
+
+        const size_t num_qubits = this->_statevector.getNumQubits();
+        auto arr_data = this->_statevector.getData();
+
+        PrecisionT expval = 0.0;
+        functor_t(arr_data, num_qubits, wires)(expval);
+
+        /*
+        size_t rev_wire = num_qubits - wires[0] - 1;
+        size_t rev_wire_shift = (static_cast<size_t>(1U) << rev_wire);
+        size_t wire_parity = fillTrailingOnes(rev_wire);
+        size_t wire_parity_inv = fillLeadingOnes(rev_wire + 1);
+
+        // ======
+        // PauliX
+        // ======
+        size_t k, i0, i1;
+        PrecisionT tmp_val = 0.0;
+        double t_start = omp_get_wtime();
+        #if defined(_OPENMP)
+        #pragma omp parallel for default(none) \
+        shared(num_qubits, wire_parity_inv, wire_parity, rev_wire_shift,
+        arr_data) \ private(k, i0, i1) reduction(+ : tmp_val) #endif for (k = 0;
+        k < exp2(num_qubits - 1); k++) { i0 =
+                ((k << 1U) & wire_parity_inv) | (wire_parity & k);
+            i1 = i0 | rev_wire_shift;
+
+            tmp_val += real(conj(arr_data[i0]) * arr_data[i1]) +
+                       real(conj(arr_data[i1]) * arr_data[i0]);
+        }
+        double t_stop = omp_get_wtime();
+        double t_elapsed = t_stop - t_start;
+        std::cout << "t_kernel_direct: " << t_elapsed << std::endl;
+        */
+
+        return expval;
+    }
 
     /**
      * @brief Probabilities of each computational basis state.
@@ -175,6 +246,83 @@ class Measurements final
         ComplexT expected_value = innerProdC(this->_statevector.getData(),
                                              operator_statevector.getData(),
                                              this->_statevector.getLength());
+
+        std::cout << "Inside expval" << std::endl;
+        std::cout << "operation: " << operation << std::endl;
+        std::cout << "wires.size(): " << wires.size() << std::endl;
+        size_t data_size = this->_statevector.getLength();
+        size_t num_qubits = this->_statevector.getNumQubits();
+        std::cout << "data_size: " << data_size << std::endl;
+        std::cout << "exp2(num_qubits - 1): " << exp2(num_qubits - 1)
+                  << std::endl;
+
+        // In-place evaluation
+        PrecisionT expval_inplace = 0;
+
+        double t_start = omp_get_wtime();
+        switch (expval_funcs_[operation]) {
+        case ExpValFunc::Identity:
+            // return
+            // applyExpValNamedFunctor<getExpectationValueIdentityFunctor,
+            //                                0>(wires);
+
+            expval_inplace =
+                applyExpValNamedFunctor<getExpectationValueIdentityFunctor, 0>(
+                    wires);
+            break;
+
+        case ExpValFunc::PauliX:
+            // return applyExpValNamedFunctor<getExpectationValuePauliXFunctor,
+            // 1>(
+            //     wires);
+            expval_inplace =
+                applyExpValNamedFunctor<getExpectationValuePauliXFunctor, 1>(
+                    wires);
+            break;
+
+        case ExpValFunc::PauliY:
+            // return applyExpValNamedFunctor<getExpectationValuePauliYFunctor,
+            // 1>(
+            //     wires);
+            expval_inplace =
+                applyExpValNamedFunctor<getExpectationValuePauliYFunctor, 1>(
+                    wires);
+            break;
+
+        case ExpValFunc::PauliZ:
+            // return applyExpValNamedFunctor<getExpectationValuePauliZFunctor,
+            // 1>(
+            //     wires);
+
+            expval_inplace =
+                applyExpValNamedFunctor<getExpectationValuePauliZFunctor, 1>(
+                    wires);
+            break;
+
+        case ExpValFunc::Hadamard:
+            // return
+            // applyExpValNamedFunctor<getExpectationValueHadamardFunctor,
+            //                                1>(wires);
+
+            expval_inplace =
+                applyExpValNamedFunctor<getExpectationValueHadamardFunctor, 1>(
+                    wires);
+            break;
+
+        default:
+            PL_ABORT(
+                std::string("Expval does not exist for named observable ") +
+                operation);
+            break;
+        }
+        double t_stop = omp_get_wtime();
+        double duration = t_stop - t_start;
+        std::cout << "duration: " << duration << std::endl;
+
+        PrecisionT err_expval =
+            std::abs(std::real(expected_value) - expval_inplace);
+        std::cout << "err_expval: " << err_expval << std::endl;
+
         return std::real(expected_value);
     };
 
@@ -661,5 +809,22 @@ class Measurements final
         }
         return init_idx;
     }
+
+    std::unordered_map<std::string, ExpValFunc> expval_funcs_;
+
+    // clang-format off
+        /**
+        * @brief Register generator operations in the generators_indices_ attribute:
+        *        an unordered_map mapping strings to GateOperation enumeration keywords.
+        */
+        void init_expval_funcs_() {
+            expval_funcs_["Identity"] = ExpValFunc::Identity;
+            expval_funcs_["PauliX"]   = ExpValFunc::PauliX;
+            expval_funcs_["PauliY"]   = ExpValFunc::PauliY;
+            expval_funcs_["PauliZ"]   = ExpValFunc::PauliZ;
+            expval_funcs_["Hadamard"] = ExpValFunc::Hadamard;
+        }
+    // clang-format on
+
 }; // class Measurements
 } // namespace Pennylane::LightningQubit::Measures
\ No newline at end of file

From 756033b9699a7d889c3d5f8afa4b36d61f27dc96 Mon Sep 17 00:00:00 2001
From: Saeed Bohloul <sbohloul@gmail.com>
Date: Tue, 21 Nov 2023 14:32:05 -0500
Subject: [PATCH 02/11] Add minor fixes to pass coverage and cpp-test

---
 .../measurements/ExpValFunctorsLQubit.hpp     | 46 +++++++++++--------
 .../measurements/MeasurementsLQubit.hpp       |  6 +--
 2 files changed, 31 insertions(+), 21 deletions(-)

diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
index 93210226d1..a3a157b1ed 100644
--- a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
+++ b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
@@ -14,22 +14,22 @@
 #pragma once
 
 #include "BitUtil.hpp"
+#include "Util.hpp"
 
 /// @cond DEV
 namespace {
 using namespace Pennylane::Util;
-using Pennylane::Util::INVSQRT2;
 } // namespace
 /// @endcond
 
 namespace Pennylane::LightningQubit::Functors {
 
 template <class PrecisionT> struct getExpectationValueIdentityFunctor {
-    std::complex<PrecisionT> *arr;
+    const std::complex<PrecisionT> *arr;
     size_t num_qubits;
 
     getExpectationValueIdentityFunctor(
-        std::complex<PrecisionT> *arr_,
+        const std::complex<PrecisionT> *arr_,
         [[maybe_unused]] std::size_t num_qubits_,
         [[maybe_unused]] const std::vector<size_t> &wires) {
         arr = arr_;
@@ -49,7 +49,7 @@ template <class PrecisionT> struct getExpectationValueIdentityFunctor {
 };
 
 template <class PrecisionT> struct getExpectationValuePauliXFunctor {
-    std::complex<PrecisionT> *arr;
+    const std::complex<PrecisionT> *arr;
     size_t num_qubits;
 
     size_t rev_wire;
@@ -57,7 +57,7 @@ template <class PrecisionT> struct getExpectationValuePauliXFunctor {
     size_t wire_parity;
     size_t wire_parity_inv;
 
-    getExpectationValuePauliXFunctor(std::complex<PrecisionT> *arr_,
+    getExpectationValuePauliXFunctor(const std::complex<PrecisionT> *arr_,
                                      std::size_t num_qubits_,
                                      const std::vector<size_t> &wires) {
         arr = arr_;
@@ -69,7 +69,9 @@ template <class PrecisionT> struct getExpectationValuePauliXFunctor {
     }
 
     inline void operator()(PrecisionT &expval) const {
-        size_t k, i0, i1;
+        size_t k;
+        size_t i0;
+        size_t i1;
 #if defined(_OPENMP)
 #pragma omp parallel for default(none)                                         \
     shared(num_qubits, wire_parity_inv, wire_parity, rev_wire_shift,           \
@@ -86,7 +88,7 @@ template <class PrecisionT> struct getExpectationValuePauliXFunctor {
 };
 
 template <class PrecisionT> struct getExpectationValuePauliYFunctor {
-    std::complex<PrecisionT> *arr;
+    const std::complex<PrecisionT> *arr;
     size_t num_qubits;
 
     size_t rev_wire;
@@ -94,7 +96,7 @@ template <class PrecisionT> struct getExpectationValuePauliYFunctor {
     size_t wire_parity;
     size_t wire_parity_inv;
 
-    getExpectationValuePauliYFunctor(std::complex<PrecisionT> *arr_,
+    getExpectationValuePauliYFunctor(const std::complex<PrecisionT> *arr_,
                                      std::size_t num_qubits_,
                                      const std::vector<size_t> &wires) {
         arr = arr_;
@@ -106,8 +108,11 @@ template <class PrecisionT> struct getExpectationValuePauliYFunctor {
     }
 
     inline void operator()(PrecisionT &expval) const {
-        size_t k, i0, i1;
-        std::complex<PrecisionT> v0, v1;
+        size_t k;
+        size_t i0;
+        size_t i1;
+        std::complex<PrecisionT> v0;
+        std::complex<PrecisionT> v1;
 #if defined(_OPENMP)
 #pragma omp parallel for default(none)                                         \
     shared(num_qubits, wire_parity_inv, wire_parity, rev_wire_shift,           \
@@ -128,7 +133,7 @@ template <class PrecisionT> struct getExpectationValuePauliYFunctor {
 };
 
 template <class PrecisionT> struct getExpectationValuePauliZFunctor {
-    std::complex<PrecisionT> *arr;
+    const std::complex<PrecisionT> *arr;
     size_t num_qubits;
 
     size_t rev_wire;
@@ -136,7 +141,7 @@ template <class PrecisionT> struct getExpectationValuePauliZFunctor {
     size_t wire_parity;
     size_t wire_parity_inv;
 
-    getExpectationValuePauliZFunctor(std::complex<PrecisionT> *arr_,
+    getExpectationValuePauliZFunctor(const std::complex<PrecisionT> *arr_,
                                      std::size_t num_qubits_,
                                      const std::vector<size_t> &wires) {
         arr = arr_;
@@ -148,7 +153,9 @@ template <class PrecisionT> struct getExpectationValuePauliZFunctor {
     }
 
     inline void operator()(PrecisionT &expval) const {
-        size_t k, i0, i1;
+        size_t k;
+        size_t i0;
+        size_t i1;
 #if defined(_OPENMP)
 #pragma omp parallel for default(none)                                         \
     shared(num_qubits, wire_parity_inv, wire_parity, rev_wire_shift,           \
@@ -166,7 +173,7 @@ template <class PrecisionT> struct getExpectationValuePauliZFunctor {
 };
 
 template <class PrecisionT> struct getExpectationValueHadamardFunctor {
-    std::complex<PrecisionT> *arr;
+    const std::complex<PrecisionT> *arr;
     size_t num_qubits;
 
     size_t rev_wire;
@@ -174,9 +181,9 @@ template <class PrecisionT> struct getExpectationValueHadamardFunctor {
     size_t wire_parity;
     size_t wire_parity_inv;
 
-    constexpr static auto isqrt2 = INVSQRT2<PrecisionT>();
+    PrecisionT isqrt2 = INVSQRT2<PrecisionT>();
 
-    getExpectationValueHadamardFunctor(std::complex<PrecisionT> *arr_,
+    getExpectationValueHadamardFunctor(const std::complex<PrecisionT> *arr_,
                                        std::size_t num_qubits_,
                                        const std::vector<size_t> &wires) {
         arr = arr_;
@@ -188,8 +195,11 @@ template <class PrecisionT> struct getExpectationValueHadamardFunctor {
     }
 
     inline void operator()(PrecisionT &expval) const {
-        size_t k, i0, i1;
-        std::complex<PrecisionT> v0, v1;
+        size_t k;
+        size_t i0;
+        size_t i1;
+        std::complex<PrecisionT> v0;
+        std::complex<PrecisionT> v1;
 #if defined(_OPENMP)
 #pragma omp parallel for default(none)                                         \
     shared(num_qubits, wire_parity_inv, wire_parity, rev_wire_shift,           \
diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
index 856a0d4033..fe39c1d611 100644
--- a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
+++ b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
@@ -103,11 +103,11 @@ class Measurements final
             PL_ASSERT(wires.size() == num_wires);
         }
 
-        const size_t num_qubits = this->_statevector.getNumQubits();
-        auto arr_data = this->_statevector.getData();
+        size_t num_qubits = this->_statevector.getNumQubits();
+        const std::complex<PrecisionT> *arr_data = this->_statevector.getData();
 
         PrecisionT expval = 0.0;
-        functor_t(arr_data, num_qubits, wires)(expval);
+        functor_t<PrecisionT>(arr_data, num_qubits, wires)(expval);
 
         /*
         size_t rev_wire = num_qubits - wires[0] - 1;

From 8c18563cba37f5c5fd5e086ab18682a1460cc6bc Mon Sep 17 00:00:00 2001
From: Saeed Bohloul <sbohloul@gmail.com>
Date: Tue, 21 Nov 2023 18:33:23 -0500
Subject: [PATCH 03/11] Minor fixes to pass test

* Fix types and const to pass tidy, cpp, format, coverage tests
---
 .../measurements/MeasurementsLQubit.hpp       | 83 +++----------------
 1 file changed, 12 insertions(+), 71 deletions(-)

diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
index fe39c1d611..1a3d0109e9 100644
--- a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
+++ b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
@@ -236,94 +236,35 @@ class Measurements final
      */
     PrecisionT expval(const std::string &operation,
                       const std::vector<size_t> &wires) {
-        // Copying the original state vector, for the application of the
-        // observable operator.
-        StateVectorLQubitManaged<PrecisionT> operator_statevector(
-            this->_statevector);
-
-        operator_statevector.applyOperation(operation, wires);
-
-        ComplexT expected_value = innerProdC(this->_statevector.getData(),
-                                             operator_statevector.getData(),
-                                             this->_statevector.getLength());
-
-        std::cout << "Inside expval" << std::endl;
-        std::cout << "operation: " << operation << std::endl;
-        std::cout << "wires.size(): " << wires.size() << std::endl;
-        size_t data_size = this->_statevector.getLength();
-        size_t num_qubits = this->_statevector.getNumQubits();
-        std::cout << "data_size: " << data_size << std::endl;
-        std::cout << "exp2(num_qubits - 1): " << exp2(num_qubits - 1)
-                  << std::endl;
-
-        // In-place evaluation
-        PrecisionT expval_inplace = 0;
-
-        double t_start = omp_get_wtime();
+        // In-place calculation of expval without creating duplicate of the
+        // statevector.
         switch (expval_funcs_[operation]) {
         case ExpValFunc::Identity:
-            // return
-            // applyExpValNamedFunctor<getExpectationValueIdentityFunctor,
-            //                                0>(wires);
-
-            expval_inplace =
-                applyExpValNamedFunctor<getExpectationValueIdentityFunctor, 0>(
-                    wires);
+            return applyExpValNamedFunctor<getExpectationValueIdentityFunctor,
+                                           0>(wires);
             break;
-
         case ExpValFunc::PauliX:
-            // return applyExpValNamedFunctor<getExpectationValuePauliXFunctor,
-            // 1>(
-            //     wires);
-            expval_inplace =
-                applyExpValNamedFunctor<getExpectationValuePauliXFunctor, 1>(
-                    wires);
+            return applyExpValNamedFunctor<getExpectationValuePauliXFunctor, 1>(
+                wires);
             break;
-
         case ExpValFunc::PauliY:
-            // return applyExpValNamedFunctor<getExpectationValuePauliYFunctor,
-            // 1>(
-            //     wires);
-            expval_inplace =
-                applyExpValNamedFunctor<getExpectationValuePauliYFunctor, 1>(
-                    wires);
+            return applyExpValNamedFunctor<getExpectationValuePauliYFunctor, 1>(
+                wires);
             break;
-
         case ExpValFunc::PauliZ:
-            // return applyExpValNamedFunctor<getExpectationValuePauliZFunctor,
-            // 1>(
-            //     wires);
-
-            expval_inplace =
-                applyExpValNamedFunctor<getExpectationValuePauliZFunctor, 1>(
-                    wires);
+            return applyExpValNamedFunctor<getExpectationValuePauliZFunctor, 1>(
+                wires);
             break;
-
         case ExpValFunc::Hadamard:
-            // return
-            // applyExpValNamedFunctor<getExpectationValueHadamardFunctor,
-            //                                1>(wires);
-
-            expval_inplace =
-                applyExpValNamedFunctor<getExpectationValueHadamardFunctor, 1>(
-                    wires);
+            return applyExpValNamedFunctor<getExpectationValueHadamardFunctor,
+                                           1>(wires);
             break;
-
         default:
             PL_ABORT(
                 std::string("Expval does not exist for named observable ") +
                 operation);
             break;
         }
-        double t_stop = omp_get_wtime();
-        double duration = t_stop - t_start;
-        std::cout << "duration: " << duration << std::endl;
-
-        PrecisionT err_expval =
-            std::abs(std::real(expected_value) - expval_inplace);
-        std::cout << "err_expval: " << err_expval << std::endl;
-
-        return std::real(expected_value);
     };
 
     /**

From edb97a554a7c11654d24b44e2810e4006e0b9b80 Mon Sep 17 00:00:00 2001
From: Saeed Bohloul <sbohloul@gmail.com>
Date: Tue, 21 Nov 2023 19:33:08 -0500
Subject: [PATCH 04/11] Add cpp test

* Added Hadamard and Identity tests to 'Expval - NamedObs'
---
 .../core/src/measurements/tests/Test_MeasurementsBase.cpp  | 7 +++++--
 1 file changed, 5 insertions(+), 2 deletions(-)

diff --git a/pennylane_lightning/core/src/measurements/tests/Test_MeasurementsBase.cpp b/pennylane_lightning/core/src/measurements/tests/Test_MeasurementsBase.cpp
index 2f6d8adce8..1cf9a37ad9 100644
--- a/pennylane_lightning/core/src/measurements/tests/Test_MeasurementsBase.cpp
+++ b/pennylane_lightning/core/src/measurements/tests/Test_MeasurementsBase.cpp
@@ -162,12 +162,15 @@ template <typename TypeList> void testNamedObsExpval() {
         Measurements<StateVectorT> Measurer(statevector);
 
         std::vector<std::vector<size_t>> wires_list = {{0}, {1}, {2}};
-        std::vector<std::string> obs_name = {"PauliX", "PauliY", "PauliZ"};
+        std::vector<std::string> obs_name = {"PauliX", "PauliY", "PauliZ",
+                                             "Hadamard", "Identity"};
         // Expected results calculated with Pennylane default.qubit:
         std::vector<std::vector<PrecisionT>> exp_values_ref = {
             {0.49272486, 0.42073549, 0.28232124},
             {-0.64421768, -0.47942553, -0.29552020},
-            {0.58498357, 0.77015115, 0.91266780}};
+            {0.58498357, 0.77015115, 0.91266780},
+            {0.76205494, 0.84208402, 0.84498485},
+            {1.0, 1.0, 1.0}};
 
         for (size_t ind_obs = 0; ind_obs < obs_name.size(); ind_obs++) {
             DYNAMIC_SECTION(obs_name[ind_obs]

From 44d8c17733ea116434ce192031c3c5d4e0349f4f Mon Sep 17 00:00:00 2001
From: Saeed Bohloul <sbohloul@gmail.com>
Date: Thu, 23 Nov 2023 14:08:13 -0500
Subject: [PATCH 05/11] Clean up dev related comments and lines

---
 .../measurements/MeasurementsLQubit.hpp       | 36 -------------------
 1 file changed, 36 deletions(-)

diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
index 1a3d0109e9..5b2d2cab32 100644
--- a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
+++ b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
@@ -39,9 +39,6 @@
 #include "TransitionKernels.hpp"
 #include "Util.hpp" //transpose_state_tensor, sorting_indices
 
-#include "BitUtil.hpp"
-#include <iostream>
-
 /// @cond DEV
 namespace {
 using namespace Pennylane::Measures;
@@ -109,34 +106,6 @@ class Measurements final
         PrecisionT expval = 0.0;
         functor_t<PrecisionT>(arr_data, num_qubits, wires)(expval);
 
-        /*
-        size_t rev_wire = num_qubits - wires[0] - 1;
-        size_t rev_wire_shift = (static_cast<size_t>(1U) << rev_wire);
-        size_t wire_parity = fillTrailingOnes(rev_wire);
-        size_t wire_parity_inv = fillLeadingOnes(rev_wire + 1);
-
-        // ======
-        // PauliX
-        // ======
-        size_t k, i0, i1;
-        PrecisionT tmp_val = 0.0;
-        double t_start = omp_get_wtime();
-        #if defined(_OPENMP)
-        #pragma omp parallel for default(none) \
-        shared(num_qubits, wire_parity_inv, wire_parity, rev_wire_shift,
-        arr_data) \ private(k, i0, i1) reduction(+ : tmp_val) #endif for (k = 0;
-        k < exp2(num_qubits - 1); k++) { i0 =
-                ((k << 1U) & wire_parity_inv) | (wire_parity & k);
-            i1 = i0 | rev_wire_shift;
-
-            tmp_val += real(conj(arr_data[i0]) * arr_data[i1]) +
-                       real(conj(arr_data[i1]) * arr_data[i0]);
-        }
-        double t_stop = omp_get_wtime();
-        double t_elapsed = t_stop - t_start;
-        std::cout << "t_kernel_direct: " << t_elapsed << std::endl;
-        */
-
         return expval;
     }
 
@@ -242,23 +211,18 @@ class Measurements final
         case ExpValFunc::Identity:
             return applyExpValNamedFunctor<getExpectationValueIdentityFunctor,
                                            0>(wires);
-            break;
         case ExpValFunc::PauliX:
             return applyExpValNamedFunctor<getExpectationValuePauliXFunctor, 1>(
                 wires);
-            break;
         case ExpValFunc::PauliY:
             return applyExpValNamedFunctor<getExpectationValuePauliYFunctor, 1>(
                 wires);
-            break;
         case ExpValFunc::PauliZ:
             return applyExpValNamedFunctor<getExpectationValuePauliZFunctor, 1>(
                 wires);
-            break;
         case ExpValFunc::Hadamard:
             return applyExpValNamedFunctor<getExpectationValueHadamardFunctor,
                                            1>(wires);
-            break;
         default:
             PL_ABORT(
                 std::string("Expval does not exist for named observable ") +

From 200d5b1db9f0880e2f9feb06d30940e9d7734ec2 Mon Sep 17 00:00:00 2001
From: Saeed Bohloul <sbohloul@gmail.com>
Date: Fri, 24 Nov 2023 01:23:16 -0500
Subject: [PATCH 06/11] Fix codefactor issues

---
 .../lightning_qubit/measurements/ExpValFunctorsLQubit.hpp        | 1 -
 .../lightning_qubit/measurements/MeasurementsLQubit.hpp          | 1 -
 2 files changed, 2 deletions(-)

diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
index a3a157b1ed..1261c1a0eb 100644
--- a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
+++ b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
@@ -162,7 +162,6 @@ template <class PrecisionT> struct getExpectationValuePauliZFunctor {
                arr) private(k, i0, i1) reduction(+ : expval)
 #endif
         for (k = 0; k < exp2(num_qubits - 1); k++) {
-
             i0 = ((k << 1U) & wire_parity_inv) | (wire_parity & k);
             i1 = i0 | rev_wire_shift;
 
diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
index 5b2d2cab32..df6dbc5d5b 100644
--- a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
+++ b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
@@ -730,6 +730,5 @@ class Measurements final
             expval_funcs_["Hadamard"] = ExpValFunc::Hadamard;
         }
     // clang-format on
-
 }; // class Measurements
 } // namespace Pennylane::LightningQubit::Measures
\ No newline at end of file

From ba48a7f7c1e594657c5bc80499f4a303907fa689 Mon Sep 17 00:00:00 2001
From: Saeed Bohloul <sbohloul@gmail.com>
Date: Fri, 24 Nov 2023 08:39:32 -0500
Subject: [PATCH 07/11] Update changelog

---
 .github/CHANGELOG.md | 3 +++
 1 file changed, 3 insertions(+)

diff --git a/.github/CHANGELOG.md b/.github/CHANGELOG.md
index 0bc27c0897..afa228aa96 100644
--- a/.github/CHANGELOG.md
+++ b/.github/CHANGELOG.md
@@ -15,6 +15,9 @@
 
 ### Improvements
 
+* Modify expval of named operators in LightningQubit for in-place computation of expectation value, to avoid creating an intermediate statevector
+  [(##565)] (https://github.com/PennyLaneAI/pennylane-lightning/pull/565)
+
 * Modify `setup.py` to use backend-specific build directory (`f"build_{backend}"`) to accelerate rebuilding backends in alternance.
   [(#540)] (https://github.com/PennyLaneAI/pennylane-lightning/pull/540)
 

From 803bc42f317d104e6479460679606ad750741b1e Mon Sep 17 00:00:00 2001
From: Saeed Bohloul <sbohloul@gmail.com>
Date: Tue, 28 Nov 2023 08:11:25 -0500
Subject: [PATCH 08/11] Add fixes and expand codecov

* Fix identity operator
* Move tests to 'Test_MeasurementsLQubit.cpp' and add Identity and Hadamard tests
---
 .../src/measurements/tests/Test_MeasurementsBase.cpp   |  7 ++-----
 .../measurements/ExpValFunctorsLQubit.hpp              |  4 ++--
 .../measurements/tests/Test_MeasurementsLQubit.cpp     | 10 ++++++++++
 3 files changed, 14 insertions(+), 7 deletions(-)

diff --git a/pennylane_lightning/core/src/measurements/tests/Test_MeasurementsBase.cpp b/pennylane_lightning/core/src/measurements/tests/Test_MeasurementsBase.cpp
index 1cf9a37ad9..2f6d8adce8 100644
--- a/pennylane_lightning/core/src/measurements/tests/Test_MeasurementsBase.cpp
+++ b/pennylane_lightning/core/src/measurements/tests/Test_MeasurementsBase.cpp
@@ -162,15 +162,12 @@ template <typename TypeList> void testNamedObsExpval() {
         Measurements<StateVectorT> Measurer(statevector);
 
         std::vector<std::vector<size_t>> wires_list = {{0}, {1}, {2}};
-        std::vector<std::string> obs_name = {"PauliX", "PauliY", "PauliZ",
-                                             "Hadamard", "Identity"};
+        std::vector<std::string> obs_name = {"PauliX", "PauliY", "PauliZ"};
         // Expected results calculated with Pennylane default.qubit:
         std::vector<std::vector<PrecisionT>> exp_values_ref = {
             {0.49272486, 0.42073549, 0.28232124},
             {-0.64421768, -0.47942553, -0.29552020},
-            {0.58498357, 0.77015115, 0.91266780},
-            {0.76205494, 0.84208402, 0.84498485},
-            {1.0, 1.0, 1.0}};
+            {0.58498357, 0.77015115, 0.91266780}};
 
         for (size_t ind_obs = 0; ind_obs < obs_name.size(); ind_obs++) {
             DYNAMIC_SECTION(obs_name[ind_obs]
diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
index 1261c1a0eb..cac8f99241 100644
--- a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
+++ b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
@@ -42,8 +42,8 @@ template <class PrecisionT> struct getExpectationValueIdentityFunctor {
 #pragma omp parallel for default(none) shared(num_qubits, arr) private(k)      \
     reduction(+ : expval)
 #endif
-        for (k = 0; k < exp2(num_qubits - 1); k++) {
-            expval += real(conj(arr[k] * arr[k]));
+        for (k = 0; k < exp2(num_qubits); k++) {
+            expval += real(conj(arr[k]) * arr[k]);
         }
     }
 };
diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/tests/Test_MeasurementsLQubit.cpp b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/tests/Test_MeasurementsLQubit.cpp
index b4be45e99a..35aa2989d0 100644
--- a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/tests/Test_MeasurementsLQubit.cpp
+++ b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/tests/Test_MeasurementsLQubit.cpp
@@ -100,6 +100,11 @@ TEMPLATE_PRODUCT_TEST_CASE("Expected Values", "[Measurements]",
         std::vector<std::vector<size_t>> wires_list = {{0}, {1}, {2}};
         std::vector<std::string> operations_list;
 
+        operations_list = {"Identity", "Identity", "Identity"};
+        exp_values = Measurer.expval(operations_list, wires_list);
+        exp_values_ref = {1.0, 1.0, 1.0};
+        REQUIRE_THAT(exp_values, Catch::Approx(exp_values_ref).margin(1e-6));
+
         operations_list = {"PauliX", "PauliX", "PauliX"};
         exp_values = Measurer.expval(operations_list, wires_list);
         exp_values_ref = {0.49272486, 0.42073549, 0.28232124};
@@ -114,6 +119,11 @@ TEMPLATE_PRODUCT_TEST_CASE("Expected Values", "[Measurements]",
         exp_values = Measurer.expval(operations_list, wires_list);
         exp_values_ref = {0.58498357, 0.77015115, 0.91266780};
         REQUIRE_THAT(exp_values, Catch::Approx(exp_values_ref).margin(1e-6));
+
+        operations_list = {"Hadamard", "Hadamard", "Hadamard"};
+        exp_values = Measurer.expval(operations_list, wires_list);
+        exp_values_ref = {0.7620549436, 0.8420840225, 0.8449848566};
+        REQUIRE_THAT(exp_values, Catch::Approx(exp_values_ref).margin(1e-6));
     }
 }
 

From 815faa1710b39fafc54d738d1354a9711ec840a9 Mon Sep 17 00:00:00 2001
From: Saeed Bohloul <sbohloul@gmail.com>
Date: Tue, 28 Nov 2023 22:17:27 -0500
Subject: [PATCH 09/11] Refactor for complying with coding standards

* Add docstrings for functors
* Return expval by value instead of passing by reference
* Make fuctors data members private
* Use initializer list for functors constructors
* Make isqrt2
* Make expval_funcs_ const and initialize with list initializer
---
 .../measurements/ExpValFunctorsLQubit.hpp     | 187 ++++++++++++++----
 .../measurements/MeasurementsLQubit.hpp       |  35 +---
 2 files changed, 154 insertions(+), 68 deletions(-)

diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
index cac8f99241..8502ee9bfe 100644
--- a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
+++ b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
@@ -11,6 +11,13 @@
 // WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 // See the License for the specific language governing permissions and
 // limitations under the License.
+
+/**
+ * @file ExpValFunctorsQubit.hpp
+ * Define functors for in-place computation of expectation value of
+ * Identity, PauliX, PauliY, PauliZ, and Hadamard operators.
+ */
+
 #pragma once
 
 #include "BitUtil.hpp"
@@ -23,21 +30,37 @@ using namespace Pennylane::Util;
 /// @endcond
 
 namespace Pennylane::LightningQubit::Functors {
-
+/**
+ * @brief Functor for computing expectation value of Identity operator.
+ *
+ * This functor provides an in-place implementation for computing the
+ * expectation value of Identity operator for a given statevector.
+ *
+ * @tparam PrecisionT precision data type of the statevector.
+ */
 template <class PrecisionT> struct getExpectationValueIdentityFunctor {
+
+  private:
     const std::complex<PrecisionT> *arr;
     size_t num_qubits;
 
+  public:
+    /**
+     * @brief Get the Expectation Value Identity Functor object
+     *
+     * @param arr_ Pointer to the statevector
+     * @param num_qubits_ Number of qubits
+     * @param wires Holds index of the wire to apply the operation on
+     */
     getExpectationValueIdentityFunctor(
         const std::complex<PrecisionT> *arr_,
         [[maybe_unused]] std::size_t num_qubits_,
-        [[maybe_unused]] const std::vector<size_t> &wires) {
-        arr = arr_;
-        num_qubits = num_qubits_;
-    }
+        [[maybe_unused]] const std::vector<size_t> &wires)
+        : arr(arr_), num_qubits(num_qubits_) {}
 
-    inline void operator()(PrecisionT &expval) const {
+    inline PrecisionT operator()() const {
         size_t k;
+        PrecisionT expval = 0.0;
 #if defined(_OPENMP)
 #pragma omp parallel for default(none) shared(num_qubits, arr) private(k)      \
     reduction(+ : expval)
@@ -45,10 +68,22 @@ template <class PrecisionT> struct getExpectationValueIdentityFunctor {
         for (k = 0; k < exp2(num_qubits); k++) {
             expval += real(conj(arr[k]) * arr[k]);
         }
+
+        return expval;
     }
 };
 
+/**
+ * @brief Functor for computing expectation value of PauliX operator.
+ *
+ * This functor provides an in-place implementation for computing the
+ * expectation value of PauliX operator for a given statevector.
+ *
+ * @tparam PrecisionT precision data type of the statevector.
+ */
 template <class PrecisionT> struct getExpectationValuePauliXFunctor {
+
+  private:
     const std::complex<PrecisionT> *arr;
     size_t num_qubits;
 
@@ -57,21 +92,28 @@ template <class PrecisionT> struct getExpectationValuePauliXFunctor {
     size_t wire_parity;
     size_t wire_parity_inv;
 
+  public:
+    /**
+     * @brief Get the Expectation Value PauliX Functor object
+     *
+     * @param arr_ Pointer to the statevector
+     * @param num_qubits_ Number of qubits
+     * @param wires Holds index of the wire to apply the operation on
+     */
     getExpectationValuePauliXFunctor(const std::complex<PrecisionT> *arr_,
                                      std::size_t num_qubits_,
-                                     const std::vector<size_t> &wires) {
-        arr = arr_;
-        num_qubits = num_qubits_;
-        rev_wire = num_qubits - wires[0] - 1;
-        rev_wire_shift = (static_cast<size_t>(1U) << rev_wire);
-        wire_parity = fillTrailingOnes(rev_wire);
-        wire_parity_inv = fillLeadingOnes(rev_wire + 1);
-    }
+                                     const std::vector<size_t> &wires)
+        : arr(arr_), num_qubits(num_qubits_),
+          rev_wire(num_qubits - wires[0] - 1),
+          rev_wire_shift((static_cast<size_t>(1U) << rev_wire)),
+          wire_parity(fillTrailingOnes(rev_wire)),
+          wire_parity_inv(fillLeadingOnes(rev_wire + 1)) {}
 
-    inline void operator()(PrecisionT &expval) const {
+    inline PrecisionT operator()() const {
         size_t k;
         size_t i0;
         size_t i1;
+        PrecisionT expval = 0.0;
 #if defined(_OPENMP)
 #pragma omp parallel for default(none)                                         \
     shared(num_qubits, wire_parity_inv, wire_parity, rev_wire_shift,           \
@@ -84,10 +126,22 @@ template <class PrecisionT> struct getExpectationValuePauliXFunctor {
             expval +=
                 real(conj(arr[i0]) * arr[i1]) + real(conj(arr[i1]) * arr[i0]);
         }
+
+        return expval;
     }
 };
 
+/**
+ * @brief Functor for computing expectation value of PauliY operator.
+ *
+ * This functor provides an in-place implementation for computing the
+ * expectation value of PauliY operator for a given statevector.
+ *
+ * @tparam PrecisionT precision data type of the statevector.
+ */
 template <class PrecisionT> struct getExpectationValuePauliYFunctor {
+
+  private:
     const std::complex<PrecisionT> *arr;
     size_t num_qubits;
 
@@ -96,23 +150,30 @@ template <class PrecisionT> struct getExpectationValuePauliYFunctor {
     size_t wire_parity;
     size_t wire_parity_inv;
 
+  public:
+    /**
+     * @brief Get the Expectation Value PauliY Functor object
+     *
+     * @param arr_ Pointer to the statevector
+     * @param num_qubits_ Number of qubits
+     * @param wires Holds index of the wire to apply the operation on
+     */
     getExpectationValuePauliYFunctor(const std::complex<PrecisionT> *arr_,
                                      std::size_t num_qubits_,
-                                     const std::vector<size_t> &wires) {
-        arr = arr_;
-        num_qubits = num_qubits_;
-        rev_wire = num_qubits - wires[0] - 1;
-        rev_wire_shift = (static_cast<size_t>(1U) << rev_wire);
-        wire_parity = fillTrailingOnes(rev_wire);
-        wire_parity_inv = fillLeadingOnes(rev_wire + 1);
-    }
+                                     const std::vector<size_t> &wires)
+        : arr(arr_), num_qubits(num_qubits_),
+          rev_wire(num_qubits - wires[0] - 1),
+          rev_wire_shift((static_cast<size_t>(1U) << rev_wire)),
+          wire_parity(fillTrailingOnes(rev_wire)),
+          wire_parity_inv(fillLeadingOnes(rev_wire + 1)) {}
 
-    inline void operator()(PrecisionT &expval) const {
+    inline PrecisionT operator()() const {
         size_t k;
         size_t i0;
         size_t i1;
         std::complex<PrecisionT> v0;
         std::complex<PrecisionT> v1;
+        PrecisionT expval = 0.0;
 #if defined(_OPENMP)
 #pragma omp parallel for default(none)                                         \
     shared(num_qubits, wire_parity_inv, wire_parity, rev_wire_shift,           \
@@ -129,10 +190,22 @@ template <class PrecisionT> struct getExpectationValuePauliYFunctor {
                       real(conj(arr[i1]) *
                            std::complex<PrecisionT>{-imag(v0), real(v0)});
         }
+
+        return expval;
     }
 };
 
+/**
+ * @brief Functor for computing expectation value of PauliZ operator.
+ *
+ * This functor provides an in-place implementation for computing the
+ * expectation value of PauliZ operator for a given statevector.
+ *
+ * @tparam PrecisionT precision data type of the statevector.
+ */
 template <class PrecisionT> struct getExpectationValuePauliZFunctor {
+
+  private:
     const std::complex<PrecisionT> *arr;
     size_t num_qubits;
 
@@ -141,21 +214,28 @@ template <class PrecisionT> struct getExpectationValuePauliZFunctor {
     size_t wire_parity;
     size_t wire_parity_inv;
 
+  public:
+    /**
+     * @brief Get the Expectation Value PauliZ Functor object
+     *
+     * @param arr_ Pointer to the statevector
+     * @param num_qubits_ Number of qubits
+     * @param wires Holds index of the wire to apply the operation on
+     */
     getExpectationValuePauliZFunctor(const std::complex<PrecisionT> *arr_,
                                      std::size_t num_qubits_,
-                                     const std::vector<size_t> &wires) {
-        arr = arr_;
-        num_qubits = num_qubits_;
-        rev_wire = num_qubits - wires[0] - 1;
-        rev_wire_shift = (static_cast<size_t>(1U) << rev_wire);
-        wire_parity = fillTrailingOnes(rev_wire);
-        wire_parity_inv = fillLeadingOnes(rev_wire + 1);
-    }
+                                     const std::vector<size_t> &wires)
+        : arr(arr_), num_qubits(num_qubits_),
+          rev_wire(num_qubits - wires[0] - 1),
+          rev_wire_shift((static_cast<size_t>(1U) << rev_wire)),
+          wire_parity(fillTrailingOnes(rev_wire)),
+          wire_parity_inv(fillLeadingOnes(rev_wire + 1)) {}
 
-    inline void operator()(PrecisionT &expval) const {
+    inline PrecisionT operator()() const {
         size_t k;
         size_t i0;
         size_t i1;
+        PrecisionT expval = 0.0;
 #if defined(_OPENMP)
 #pragma omp parallel for default(none)                                         \
     shared(num_qubits, wire_parity_inv, wire_parity, rev_wire_shift,           \
@@ -168,10 +248,22 @@ template <class PrecisionT> struct getExpectationValuePauliZFunctor {
             expval += real(conj(arr[i1]) * (-arr[i1])) +
                       real(conj(arr[i0]) * (arr[i0]));
         }
+
+        return expval;
     }
 };
 
+/**
+ * @brief Functor for computing expectation value of Hadamard operator.
+ *
+ * This functor provides an in-place implementation for computing the
+ * expectation value of Hadamard operator for a given statevector.
+ *
+ * @tparam PrecisionT precision data type of the statevector.
+ */
 template <class PrecisionT> struct getExpectationValueHadamardFunctor {
+
+  private:
     const std::complex<PrecisionT> *arr;
     size_t num_qubits;
 
@@ -180,25 +272,32 @@ template <class PrecisionT> struct getExpectationValueHadamardFunctor {
     size_t wire_parity;
     size_t wire_parity_inv;
 
-    PrecisionT isqrt2 = INVSQRT2<PrecisionT>();
+    constexpr static auto isqrt2 = INVSQRT2<PrecisionT>();
 
+  public:
+    /**
+     * @brief Get the Expectation Value Hadamard Functor object
+     *
+     * @param arr_ Pointer to the statevector
+     * @param num_qubits_ Number of qubits
+     * @param wires Holds index of the wire to apply the operation on
+     */
     getExpectationValueHadamardFunctor(const std::complex<PrecisionT> *arr_,
                                        std::size_t num_qubits_,
-                                       const std::vector<size_t> &wires) {
-        arr = arr_;
-        num_qubits = num_qubits_;
-        rev_wire = num_qubits - wires[0] - 1;
-        rev_wire_shift = (static_cast<size_t>(1U) << rev_wire);
-        wire_parity = fillTrailingOnes(rev_wire);
-        wire_parity_inv = fillLeadingOnes(rev_wire + 1);
-    }
+                                       const std::vector<size_t> &wires)
+        : arr(arr_), num_qubits(num_qubits_),
+          rev_wire(num_qubits - wires[0] - 1),
+          rev_wire_shift((static_cast<size_t>(1U) << rev_wire)),
+          wire_parity(fillTrailingOnes(rev_wire)),
+          wire_parity_inv(fillLeadingOnes(rev_wire + 1)) {}
 
-    inline void operator()(PrecisionT &expval) const {
+    inline PrecisionT operator()() const {
         size_t k;
         size_t i0;
         size_t i1;
         std::complex<PrecisionT> v0;
         std::complex<PrecisionT> v1;
+        PrecisionT expval = 0.0;
 #if defined(_OPENMP)
 #pragma omp parallel for default(none)                                         \
     shared(num_qubits, wire_parity_inv, wire_parity, rev_wire_shift,           \
@@ -214,6 +313,8 @@ template <class PrecisionT> struct getExpectationValueHadamardFunctor {
             expval += real(isqrt2 * (conj(arr[i0]) * (v0 + v1) +
                                      conj(arr[i1]) * (v0 - v1)));
         }
+
+        return expval;
     }
 };
 } // namespace Pennylane::LightningQubit::Functors
\ No newline at end of file
diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
index df6dbc5d5b..92bc3d05e6 100644
--- a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
+++ b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
@@ -46,7 +46,6 @@ using namespace Pennylane::Observables;
 using Pennylane::LightningQubit::StateVectorLQubitManaged;
 using Pennylane::LightningQubit::Util::innerProdC;
 using namespace Pennylane::LightningQubit::Functors;
-using Pennylane::Util::INVSQRT2;
 enum class ExpValFunc : uint32_t {
     BEGIN = 1,
     Identity = 1,
@@ -79,9 +78,12 @@ class Measurements final
 
   public:
     explicit Measurements(const StateVectorT &statevector)
-        : BaseType{statevector} {
-        init_expval_funcs_();
-    };
+        : BaseType{statevector},
+          expval_funcs_{{"Identity", ExpValFunc::Identity},
+                        {"PauliX", ExpValFunc::PauliX},
+                        {"PauliY", ExpValFunc::PauliY},
+                        {"PauliZ", ExpValFunc::PauliZ},
+                        {"Hadamard", ExpValFunc::Hadamard}} {};
 
     /**
      * @brief Templated method that returns the expectation value of named
@@ -103,10 +105,7 @@ class Measurements final
         size_t num_qubits = this->_statevector.getNumQubits();
         const std::complex<PrecisionT> *arr_data = this->_statevector.getData();
 
-        PrecisionT expval = 0.0;
-        functor_t<PrecisionT>(arr_data, num_qubits, wires)(expval);
-
-        return expval;
+        return functor_t<PrecisionT>(arr_data, num_qubits, wires)();
     }
 
     /**
@@ -207,7 +206,7 @@ class Measurements final
                       const std::vector<size_t> &wires) {
         // In-place calculation of expval without creating duplicate of the
         // statevector.
-        switch (expval_funcs_[operation]) {
+        switch (expval_funcs_.at(operation)) {
         case ExpValFunc::Identity:
             return applyExpValNamedFunctor<getExpectationValueIdentityFunctor,
                                            0>(wires);
@@ -634,6 +633,8 @@ class Measurements final
     }
 
   private:
+    const std::unordered_map<std::string, ExpValFunc> expval_funcs_;
+
     /**
      * @brief Support function that calculates <bra|obs|ket> to obtain the
      * observable's expectation value.
@@ -714,21 +715,5 @@ class Measurements final
         }
         return init_idx;
     }
-
-    std::unordered_map<std::string, ExpValFunc> expval_funcs_;
-
-    // clang-format off
-        /**
-        * @brief Register generator operations in the generators_indices_ attribute:
-        *        an unordered_map mapping strings to GateOperation enumeration keywords.
-        */
-        void init_expval_funcs_() {
-            expval_funcs_["Identity"] = ExpValFunc::Identity;
-            expval_funcs_["PauliX"]   = ExpValFunc::PauliX;
-            expval_funcs_["PauliY"]   = ExpValFunc::PauliY;
-            expval_funcs_["PauliZ"]   = ExpValFunc::PauliZ;
-            expval_funcs_["Hadamard"] = ExpValFunc::Hadamard;
-        }
-    // clang-format on
 }; // class Measurements
 } // namespace Pennylane::LightningQubit::Measures
\ No newline at end of file

From d1a0e32f0d24d8aa5d113fc8220fa3ba0c2cbba9 Mon Sep 17 00:00:00 2001
From: Saeed Bohloul <sbohloul@gmail.com>
Date: Tue, 28 Nov 2023 22:40:01 -0500
Subject: [PATCH 10/11] Fix CodeFactor issues

---
 .../lightning_qubit/measurements/ExpValFunctorsLQubit.hpp    | 5 -----
 1 file changed, 5 deletions(-)

diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
index 8502ee9bfe..c9a1e84795 100644
--- a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
+++ b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
@@ -39,7 +39,6 @@ namespace Pennylane::LightningQubit::Functors {
  * @tparam PrecisionT precision data type of the statevector.
  */
 template <class PrecisionT> struct getExpectationValueIdentityFunctor {
-
   private:
     const std::complex<PrecisionT> *arr;
     size_t num_qubits;
@@ -82,7 +81,6 @@ template <class PrecisionT> struct getExpectationValueIdentityFunctor {
  * @tparam PrecisionT precision data type of the statevector.
  */
 template <class PrecisionT> struct getExpectationValuePauliXFunctor {
-
   private:
     const std::complex<PrecisionT> *arr;
     size_t num_qubits;
@@ -140,7 +138,6 @@ template <class PrecisionT> struct getExpectationValuePauliXFunctor {
  * @tparam PrecisionT precision data type of the statevector.
  */
 template <class PrecisionT> struct getExpectationValuePauliYFunctor {
-
   private:
     const std::complex<PrecisionT> *arr;
     size_t num_qubits;
@@ -204,7 +201,6 @@ template <class PrecisionT> struct getExpectationValuePauliYFunctor {
  * @tparam PrecisionT precision data type of the statevector.
  */
 template <class PrecisionT> struct getExpectationValuePauliZFunctor {
-
   private:
     const std::complex<PrecisionT> *arr;
     size_t num_qubits;
@@ -262,7 +258,6 @@ template <class PrecisionT> struct getExpectationValuePauliZFunctor {
  * @tparam PrecisionT precision data type of the statevector.
  */
 template <class PrecisionT> struct getExpectationValueHadamardFunctor {
-
   private:
     const std::complex<PrecisionT> *arr;
     size_t num_qubits;

From 20133669767a965ea5c172af89fd55dee9cd7a60 Mon Sep 17 00:00:00 2001
From: Saeed Bohloul <sbohloul@gmail.com>
Date: Wed, 29 Nov 2023 18:34:24 -0500
Subject: [PATCH 11/11] Minor refactor

* Change functors from struct to class
---
 .../measurements/ExpValFunctorsLQubit.hpp              | 10 +++++-----
 .../measurements/MeasurementsLQubit.hpp                |  4 ++--
 2 files changed, 7 insertions(+), 7 deletions(-)

diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
index c9a1e84795..8fec14ce8e 100644
--- a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
+++ b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/ExpValFunctorsLQubit.hpp
@@ -38,7 +38,7 @@ namespace Pennylane::LightningQubit::Functors {
  *
  * @tparam PrecisionT precision data type of the statevector.
  */
-template <class PrecisionT> struct getExpectationValueIdentityFunctor {
+template <class PrecisionT> class getExpectationValueIdentityFunctor {
   private:
     const std::complex<PrecisionT> *arr;
     size_t num_qubits;
@@ -80,7 +80,7 @@ template <class PrecisionT> struct getExpectationValueIdentityFunctor {
  *
  * @tparam PrecisionT precision data type of the statevector.
  */
-template <class PrecisionT> struct getExpectationValuePauliXFunctor {
+template <class PrecisionT> class getExpectationValuePauliXFunctor {
   private:
     const std::complex<PrecisionT> *arr;
     size_t num_qubits;
@@ -137,7 +137,7 @@ template <class PrecisionT> struct getExpectationValuePauliXFunctor {
  *
  * @tparam PrecisionT precision data type of the statevector.
  */
-template <class PrecisionT> struct getExpectationValuePauliYFunctor {
+template <class PrecisionT> class getExpectationValuePauliYFunctor {
   private:
     const std::complex<PrecisionT> *arr;
     size_t num_qubits;
@@ -200,7 +200,7 @@ template <class PrecisionT> struct getExpectationValuePauliYFunctor {
  *
  * @tparam PrecisionT precision data type of the statevector.
  */
-template <class PrecisionT> struct getExpectationValuePauliZFunctor {
+template <class PrecisionT> class getExpectationValuePauliZFunctor {
   private:
     const std::complex<PrecisionT> *arr;
     size_t num_qubits;
@@ -257,7 +257,7 @@ template <class PrecisionT> struct getExpectationValuePauliZFunctor {
  *
  * @tparam PrecisionT precision data type of the statevector.
  */
-template <class PrecisionT> struct getExpectationValueHadamardFunctor {
+template <class PrecisionT> class getExpectationValueHadamardFunctor {
   private:
     const std::complex<PrecisionT> *arr;
     size_t num_qubits;
diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
index 92bc3d05e6..034b01ec6d 100644
--- a/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
+++ b/pennylane_lightning/core/src/simulators/lightning_qubit/measurements/MeasurementsLQubit.hpp
@@ -76,6 +76,8 @@ class Measurements final
     using ComplexT = typename StateVectorT::ComplexT;
     using BaseType = MeasurementsBase<StateVectorT, Measurements<StateVectorT>>;
 
+    const std::unordered_map<std::string, ExpValFunc> expval_funcs_;
+
   public:
     explicit Measurements(const StateVectorT &statevector)
         : BaseType{statevector},
@@ -633,8 +635,6 @@ class Measurements final
     }
 
   private:
-    const std::unordered_map<std::string, ExpValFunc> expval_funcs_;
-
     /**
      * @brief Support function that calculates <bra|obs|ket> to obtain the
      * observable's expectation value.