diff --git a/.github/CHANGELOG.md b/.github/CHANGELOG.md index b32d65eb8..05a8371b5 100644 --- a/.github/CHANGELOG.md +++ b/.github/CHANGELOG.md @@ -5,6 +5,9 @@ * Add native N-controlled generators and adjoint support to `lightning.gpu`'s single-GPU backend. [(#970)](https://github.com/PennyLaneAI/pennylane-lightning/pull/970) +* Add a Catalyst-specific wrapping class for Lightning Qubit. + [(#960)](https://github.com/PennyLaneAI/pennylane-lightning/pull/960) + * Add native N-controlled gates support to `lightning.gpu`'s single-GPU backend. [(#938)](https://github.com/PennyLaneAI/pennylane-lightning/pull/938) @@ -35,7 +38,7 @@ This release contains contributions from (in alphabetical order): -Ali Asadi, Joseph Lee, Luis Alfredo Nuñez Meneses, Shuli Shu +Ali Asadi, Joseph Lee, Luis Alfredo Nuñez Meneses, Shuli Shu, Raul Torres, Haochen Paul Wang --- diff --git a/.github/workflows/tests_lqcpu_python.yml b/.github/workflows/tests_lqcpu_python.yml index e6734ac20..b89056a4c 100644 --- a/.github/workflows/tests_lqcpu_python.yml +++ b/.github/workflows/tests_lqcpu_python.yml @@ -197,10 +197,11 @@ jobs: - name: Run PennyLane-Lightning unit tests run: | DEVICENAME=`echo ${{ matrix.pl_backend }} | sed "s/_/./g"` - OMP_NUM_THREADS=1 PL_DEVICE=${DEVICENAME} python -m pytest -n auto tests/ -k "not unitary_correct" \ + # Remove `python -m` to avoid running tests with relative modules + OMP_NUM_THREADS=1 PL_DEVICE=${DEVICENAME} pytest -n auto tests/ -k "not unitary_correct" \ $COVERAGE_FLAGS --splits 4 --group ${{ matrix.group }} \ --durations-path='.github/workflows/python_lightning_qubit_test_durations.json' --splitting-algorithm=least_duration - PL_DEVICE=${DEVICENAME} python -m pytest tests/ -k "unitary_correct" $COVERAGE_FLAGS --cov-append + PL_DEVICE=${DEVICENAME} pytest tests/ -k "unitary_correct" $COVERAGE_FLAGS --cov-append pl-device-test --device ${DEVICENAME} --skip-ops --shots=20000 $COVERAGE_FLAGS --cov-append pl-device-test --device ${DEVICENAME} --shots=None --skip-ops $COVERAGE_FLAGS --cov-append mv .coverage ${{ github.workspace }}/.coverage-${{ github.job }}-${{ matrix.pl_backend }}-${{ matrix.blas }}-${{ matrix.group }} diff --git a/cmake/support_catalyst.cmake b/cmake/support_catalyst.cmake index 1f5c55622..9e8c8b673 100644 --- a/cmake/support_catalyst.cmake +++ b/cmake/support_catalyst.cmake @@ -42,7 +42,7 @@ macro(FindCatalyst target_name) ${HEADER_NAME} URL https://raw.githubusercontent.com/PennyLaneAI/catalyst/${CATALYST_GIT_TAG}/runtime/lib/backend/common/${HEADER} DOWNLOAD_NO_EXTRACT True - SOURCE_DIR include + SOURCE_DIR ../../include ) FetchContent_MakeAvailable(${HEADER_NAME}) @@ -62,13 +62,13 @@ macro(FindCatalyst target_name) ${HEADER_NAME} URL https://raw.githubusercontent.com/PennyLaneAI/catalyst/${CATALYST_GIT_TAG}/runtime/include/${HEADER} DOWNLOAD_NO_EXTRACT True - SOURCE_DIR include + SOURCE_DIR ../../include ) FetchContent_MakeAvailable(${HEADER_NAME}) endforeach() - target_include_directories(${target_name} PUBLIC ${CMAKE_CURRENT_BINARY_DIR}/include) + target_include_directories(${target_name} PUBLIC ${CMAKE_CURRENT_BINARY_DIR}/../../include) endif() endmacro() diff --git a/pennylane_lightning/core/_version.py b/pennylane_lightning/core/_version.py index 84b4957c2..c625ad341 100644 --- a/pennylane_lightning/core/_version.py +++ b/pennylane_lightning/core/_version.py @@ -16,4 +16,4 @@ Version number (major.minor.patch[-label]) """ -__version__ = "0.40.0-dev6" +__version__ = "0.40.0-dev7" diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/CMakeLists.txt b/pennylane_lightning/core/src/simulators/lightning_qubit/CMakeLists.txt index da457fd80..b11e3f70e 100644 --- a/pennylane_lightning/core/src/simulators/lightning_qubit/CMakeLists.txt +++ b/pennylane_lightning/core/src/simulators/lightning_qubit/CMakeLists.txt @@ -89,6 +89,13 @@ set(COMPONENT_SUBDIRS algorithms utils ) +# Do not build the LQ plugin for Windows +if (NOT WIN32) + set(COMPONENT_SUBDIRS ${COMPONENT_SUBDIRS} + catalyst + ) +endif() + foreach(COMP ${COMPONENT_SUBDIRS}) add_subdirectory(${COMP}) endforeach() diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/CMakeLists.txt b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/CMakeLists.txt new file mode 100644 index 000000000..81ee885bb --- /dev/null +++ b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/CMakeLists.txt @@ -0,0 +1,24 @@ +cmake_minimum_required(VERSION 3.20) + +project(lightning_qubit_catalyst LANGUAGES CXX) + +set(LQ_CATALYST_FILES LightningSimulator.cpp CACHE INTERNAL "") +add_library(lightning_qubit_catalyst SHARED ${LQ_CATALYST_FILES}) + +include(FetchContent) + +include("${pennylane_lightning_SOURCE_DIR}/cmake/support_catalyst.cmake") +FindCatalyst(lightning_qubit_catalyst) + +target_include_directories(lightning_qubit_catalyst INTERFACE ${CMAKE_CURRENT_SOURCE_DIR}) +target_link_libraries(lightning_qubit_catalyst PUBLIC lightning_compile_options + lightning_external_libs + lightning_qubit_algorithms + lightning_qubit_measurements + lightning_qubit_observables +) + +if (BUILD_TESTS) + enable_testing() + add_subdirectory("tests") +endif() diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/LightningObsManager.hpp b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/LightningObsManager.hpp new file mode 100644 index 000000000..6a92a220e --- /dev/null +++ b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/LightningObsManager.hpp @@ -0,0 +1,205 @@ +// Copyright 2022 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 +#include +#include +#include + +#include "Exception.hpp" +#include "Types.h" +#include "Utils.hpp" + +#include "ObservablesLQubit.hpp" +#include "StateVectorLQubitDynamic.hpp" + +namespace { +using Pennylane::LightningQubit::StateVectorLQubitDynamic; +using Pennylane::LightningQubit::Observables::HermitianObs; +using Pennylane::LightningQubit::Observables::NamedObs; +using Pennylane::LightningQubit::Observables::TensorProdObs; +using Pennylane::Observables::Observable; +} // namespace + +namespace Catalyst::Runtime::Simulator { + +/** + * @brief The LightningObsManager caches observables of a program at runtime + * and maps each one to a const unique index (`int64_t`) in the scope + * of the global context manager. + */ +template class LightningObsManager { + private: + using VectorStateT = StateVectorLQubitDynamic; + using ComplexT = typename VectorStateT::ComplexT; + using ObservablePairType = + std::pair>, ObsType>; + std::vector observables_{}; + + public: + LightningObsManager() = default; + ~LightningObsManager() = default; + + LightningObsManager(const LightningObsManager &) = delete; + LightningObsManager &operator=(const LightningObsManager &) = delete; + LightningObsManager(LightningObsManager &&) = delete; + LightningObsManager &operator=(LightningObsManager &&) = delete; + + /** + * @brief A helper function to clear constructed observables in the program. + */ + void clear() { observables_.clear(); } + + /** + * @brief Check the validity of observable keys. + * + * @param obsKeys The vector of observable keys + * @return bool + */ + [[nodiscard]] auto + isValidObservables(const std::vector &obsKeys) const -> bool { + return std::all_of(obsKeys.begin(), obsKeys.end(), [this](auto i) { + return (i >= 0 && static_cast(i) < observables_.size()); + }); + } + + /** + * @brief Get the constructed observable instance. + * + * @param key The observable key + * @return std::shared_ptr> + */ + [[nodiscard]] auto getObservable(ObsIdType key) + -> std::shared_ptr> { + RT_FAIL_IF(!isValidObservables({key}), "Invalid observable key"); + return std::get<0>(observables_[key]); + } + + /** + * @brief Get the number of observables. + * + * @return size_t + */ + [[nodiscard]] auto numObservables() const -> size_t { + return observables_.size(); + } + + /** + * @brief Create and cache a new NamedObs instance. + * + * @param obsId The named observable id of type ObsId + * @param wires The vector of wires the observable acts on + * @return ObsIdType + */ + [[nodiscard]] auto createNamedObs(ObsId obsId, + const std::vector &wires) + -> ObsIdType { + auto &&obs_str = std::string( + Lightning::lookup_obs( + Lightning::simulator_observable_support, obsId)); + + observables_.push_back(std::make_pair( + std::make_shared>(obs_str, wires), + ObsType::Basic)); + return static_cast(observables_.size() - 1); + } + + /** + * @brief Create and cache a new HermitianObs instance. + * + * @param matrix The row-wise Hermitian matrix + * @param wires The vector of wires the observable acts on + * @return ObsIdType + */ + [[nodiscard]] auto createHermitianObs(const std::vector &matrix, + const std::vector &wires) + -> ObsIdType { + observables_.push_back( + std::make_pair(std::make_shared>( + HermitianObs{matrix, wires}), + ObsType::Basic)); + + return static_cast(observables_.size() - 1); + } + + /** + * @brief Create and cache a new TensorProd instance. + * + * @param obsKeys The vector of observable keys + * @return ObsIdType + */ + [[nodiscard]] auto + createTensorProdObs(const std::vector &obsKeys) -> ObsIdType { + const auto key_size = obsKeys.size(); + const auto obs_size = observables_.size(); + + std::vector>> obs_vec; + obs_vec.reserve(key_size); + + for (const auto &key : obsKeys) { + RT_FAIL_IF(static_cast(key) >= obs_size || key < 0, + "Invalid observable key"); + + auto &&[obs, type] = observables_[key]; + obs_vec.push_back(obs); + } + + observables_.push_back(std::make_pair( + TensorProdObs::create(obs_vec), ObsType::TensorProd)); + + return static_cast(obs_size); + } + + /** + * @brief Create and cache a new HamiltonianObs instance. + * + * @param coeffs The vector of coefficients + * @param obsKeys The vector of observable keys + * @return ObsIdType + */ + [[nodiscard]] auto + createHamiltonianObs(const std::vector &coeffs, + const std::vector &obsKeys) -> ObsIdType { + const auto key_size = obsKeys.size(); + const auto obs_size = observables_.size(); + + RT_FAIL_IF( + key_size != coeffs.size(), + "Incompatible list of observables and coefficients; " + "Number of observables and number of coefficients must be equal"); + + std::vector>> obs_vec; + obs_vec.reserve(key_size); + + for (auto key : obsKeys) { + RT_FAIL_IF(static_cast(key) >= obs_size || key < 0, + "Invalid observable key"); + + auto &&[obs, type] = observables_[key]; + obs_vec.push_back(obs); + } + + observables_.push_back(std::make_pair( + std::make_shared>( + Pennylane::LightningQubit::Observables::Hamiltonian< + VectorStateT>(coeffs, std::move(obs_vec))), + ObsType::Hamiltonian)); + + return static_cast(obs_size); + } +}; +} // namespace Catalyst::Runtime::Simulator diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/LightningSimulator.cpp b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/LightningSimulator.cpp new file mode 100644 index 000000000..cec08e162 --- /dev/null +++ b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/LightningSimulator.cpp @@ -0,0 +1,564 @@ +// Copyright 2022 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. + +#include "LightningSimulator.hpp" + +#include "AdjointJacobianLQubit.hpp" +#include "JacobianData.hpp" +#include "LinearAlgebra.hpp" +#include "MeasurementsLQubit.hpp" + +namespace Catalyst::Runtime::Simulator { + +auto LightningSimulator::AllocateQubit() -> QubitIdType { + size_t sv_id = this->device_sv->allocateWire(); + return this->qubit_manager.Allocate(sv_id); +} + +auto LightningSimulator::AllocateQubits(size_t num_qubits) + -> std::vector { + if (num_qubits == 0U) { + return {}; + } + + // at the first call when num_qubits == 0 + if (this->GetNumQubits() == 0U) { + this->device_sv = std::make_unique(num_qubits); + return this->qubit_manager.AllocateRange(0, num_qubits); + } + + std::vector result(num_qubits); + std::generate_n(result.begin(), num_qubits, + [this]() { return AllocateQubit(); }); + return result; +} + +void LightningSimulator::ReleaseAllQubits() { + this->device_sv->clearData(); + this->qubit_manager.ReleaseAll(); +} + +void LightningSimulator::ReleaseQubit(QubitIdType q) { + if (this->qubit_manager.isValidQubitId(q)) { + this->device_sv->releaseWire(this->qubit_manager.getDeviceId(q)); + } + this->qubit_manager.Release(q); +} + +auto LightningSimulator::GetNumQubits() const -> size_t { + return this->device_sv->getNumQubits(); +} + +void LightningSimulator::StartTapeRecording() { + RT_FAIL_IF(this->tape_recording, "Cannot re-activate the cache manager"); + this->tape_recording = true; + this->cache_manager.Reset(); +} + +void LightningSimulator::StopTapeRecording() { + RT_FAIL_IF(!this->tape_recording, + "Cannot stop an already stopped cache manager"); + this->tape_recording = false; +} + +auto LightningSimulator::CacheManagerInfo() + -> std::tuple, + std::vector> { + return {this->cache_manager.getNumOperations(), + this->cache_manager.getNumObservables(), + this->cache_manager.getNumParams(), + this->cache_manager.getOperationsNames(), + this->cache_manager.getObservablesKeys()}; +} + +void LightningSimulator::SetDeviceShots(size_t shots) { + this->device_shots = shots; +} + +auto LightningSimulator::GetDeviceShots() const -> size_t { + return this->device_shots; +} + +void LightningSimulator::SetDevicePRNG(std::mt19937 *gen) { this->gen = gen; } + +void LightningSimulator::PrintState() { + using std::cout; + using std::endl; + + const size_t num_qubits = this->device_sv->getNumQubits(); + const size_t size = Pennylane::Util::exp2(num_qubits); + size_t idx = 0; + cout << "*** State-Vector of Size " << size << " ***\n"; + cout << "["; + auto &&state = this->device_sv->getDataVector(); + for (; idx < size - 1; idx++) { + cout << state[idx] << ", "; + } + cout << state[idx] << "]" << endl; +} + +void LightningSimulator::SetState(DataView, 1> &state, + std::vector &wires) { + std::vector> data_vector(state.begin(), state.end()); + this->device_sv->setStateVector(data_vector, getDeviceWires(wires)); +} + +void LightningSimulator::SetBasisState(DataView &n, + std::vector &wires) { + std::vector data_vector(n.begin(), n.end()); + this->device_sv->setBasisState(data_vector, getDeviceWires(wires)); +} + +auto LightningSimulator::Zero() const -> Result { + return const_cast(&GLOBAL_RESULT_FALSE_CONST); +} + +auto LightningSimulator::One() const -> Result { + return const_cast(&GLOBAL_RESULT_TRUE_CONST); +} + +void LightningSimulator::NamedOperation( + const std::string &name, const std::vector ¶ms, + const std::vector &wires, bool inverse, + const std::vector &controlled_wires, + const std::vector &controlled_values) { + // Check the validity of number of qubits and parameters + RT_FAIL_IF(controlled_wires.size() != controlled_values.size(), + "Controlled wires/values size mismatch"); + RT_FAIL_IF(!isValidQubits(wires), "Given wires do not refer to qubits"); + RT_FAIL_IF(!isValidQubits(controlled_wires), + "Given controlled wires do not refer to qubits"); + + // Convert wires to device wires + auto &&dev_wires = getDeviceWires(wires); + auto &&dev_controlled_wires = getDeviceWires(controlled_wires); + + // Update the state-vector + if (controlled_wires.empty()) { + this->device_sv->applyOperation(name, dev_wires, inverse, params); + } else { + this->device_sv->applyOperation(name, dev_controlled_wires, + controlled_values, dev_wires, inverse, + params); + } + + // Update tape caching if required + if (this->tape_recording) { + this->cache_manager.addOperation(name, params, dev_wires, inverse, {}, + dev_controlled_wires, + controlled_values); + } +} + +void LightningSimulator::MatrixOperation( + const std::vector> &matrix, + const std::vector &wires, bool inverse, + const std::vector &controlled_wires, + const std::vector &controlled_values) { + RT_FAIL_IF(controlled_wires.size() != controlled_values.size(), + "Controlled wires/values size mismatch"); + RT_FAIL_IF(!isValidQubits(wires), "Given wires do not refer to qubits"); + RT_FAIL_IF(!isValidQubits(controlled_wires), + "Given controlled wires do not refer to qubits"); + + // Convert wires to device wires + // with checking validity of wires + auto &&dev_wires = getDeviceWires(wires); + auto &&dev_controlled_wires = getDeviceWires(controlled_wires); + + // Update the state-vector + if (controlled_wires.empty()) { + this->device_sv->applyMatrix(matrix.data(), dev_wires, inverse); + } else { + this->device_sv->applyControlledMatrix( + matrix.data(), dev_controlled_wires, controlled_values, dev_wires, + inverse); + } + + // Update tape caching if required + if (this->tape_recording) { + this->cache_manager.addOperation("QubitUnitary", {}, dev_wires, inverse, + matrix, dev_controlled_wires, + controlled_values); + } +} + +auto LightningSimulator::Observable( + ObsId id, const std::vector> &matrix, + const std::vector &wires) -> ObsIdType { + RT_FAIL_IF(wires.size() > this->GetNumQubits(), "Invalid number of wires"); + RT_FAIL_IF(!isValidQubits(wires), "Invalid given wires"); + + auto &&dev_wires = getDeviceWires(wires); + + if (id == ObsId::Hermitian) { + return this->obs_manager.createHermitianObs(matrix, dev_wires); + } + + return this->obs_manager.createNamedObs(id, dev_wires); +} + +auto LightningSimulator::TensorObservable(const std::vector &obs) + -> ObsIdType { + return this->obs_manager.createTensorProdObs(obs); +} + +auto LightningSimulator::HamiltonianObservable( + const std::vector &coeffs, const std::vector &obs) + -> ObsIdType { + return this->obs_manager.createHamiltonianObs(coeffs, obs); +} + +auto LightningSimulator::Expval(ObsIdType obsKey) -> double { + RT_FAIL_IF(!this->obs_manager.isValidObservables({obsKey}), + "Invalid key for cached observables"); + auto &&obs = this->obs_manager.getObservable(obsKey); + + // update tape caching + if (this->tape_recording) { + this->cache_manager.addObservable(obsKey, MeasurementsT::Expval); + } + + Pennylane::LightningQubit::Measures::Measurements m{ + *(this->device_sv)}; + + return (device_shots != 0U) ? m.expval(*obs, device_shots, {}) + : m.expval(*obs); +} + +auto LightningSimulator::Var(ObsIdType obsKey) -> double { + RT_FAIL_IF(!this->obs_manager.isValidObservables({obsKey}), + "Invalid key for cached observables"); + auto &&obs = this->obs_manager.getObservable(obsKey); + + // update tape caching + if (this->tape_recording) { + this->cache_manager.addObservable(obsKey, MeasurementsT::Var); + } + + Pennylane::LightningQubit::Measures::Measurements m{ + *(this->device_sv)}; + + return (device_shots != 0U) ? m.var(*obs, device_shots) : m.var(*obs); +} + +void LightningSimulator::State(DataView, 1> &state) { + auto &&dv_state = this->device_sv->getDataVector(); + RT_FAIL_IF(state.size() != dv_state.size(), + "Invalid size for the pre-allocated state vector"); + + std::move(dv_state.begin(), dv_state.end(), state.begin()); +} + +void LightningSimulator::Probs(DataView &probs) { + Pennylane::LightningQubit::Measures::Measurements m{ + *(this->device_sv)}; + auto &&dv_probs = (device_shots != 0U) ? m.probs(device_shots) : m.probs(); + + RT_FAIL_IF(probs.size() != dv_probs.size(), + "Invalid size for the pre-allocated probabilities"); + + std::move(dv_probs.begin(), dv_probs.end(), probs.begin()); +} + +void LightningSimulator::PartialProbs(DataView &probs, + const std::vector &wires) { + const size_t numWires = wires.size(); + const size_t numQubits = this->GetNumQubits(); + + RT_FAIL_IF(numWires > numQubits, "Invalid number of wires"); + RT_FAIL_IF(!isValidQubits(wires), "Invalid given wires to measure"); + + auto dev_wires = getDeviceWires(wires); + Pennylane::LightningQubit::Measures::Measurements m{ + *(this->device_sv)}; + auto &&dv_probs = (device_shots != 0U) ? m.probs(dev_wires, device_shots) + : m.probs(dev_wires); + + RT_FAIL_IF(probs.size() != dv_probs.size(), + "Invalid size for the pre-allocated partial-probabilities"); + + std::move(dv_probs.begin(), dv_probs.end(), probs.begin()); +} + +std::vector +LightningSimulator::GenerateSamplesMetropolis(size_t shots) { + // generate_samples_metropolis is a member function of the Measures class. + Pennylane::LightningQubit::Measures::Measurements m{ + *(this->device_sv)}; + + // PL-Lightning generates samples using the alias method. + // Reference: https://en.wikipedia.org/wiki/Alias_method + // Given the number of samples, returns 1-D vector of samples + // in binary, each sample is separated by a stride equal to + // the number of qubits. + // + // Return Value Optimization (RVO) + return m.generate_samples_metropolis(this->kernel_name, this->num_burnin, + shots); +} + +std::vector LightningSimulator::GenerateSamples(size_t shots) { + if (this->mcmc) { + return this->GenerateSamplesMetropolis(shots); + } + // generate_samples is a member function of the Measures class. + Pennylane::LightningQubit::Measures::Measurements m{ + *(this->device_sv)}; + + // PL-Lightning generates samples using the alias method. + // Reference: https://en.wikipedia.org/wiki/Alias_method + // Given the number of samples, returns 1-D vector of samples + // in binary, each sample is separated by a stride equal to + // the number of qubits. + // + // Return Value Optimization (RVO) + if (this->gen != nullptr) { + return m.generate_samples(shots, (*(this->gen))()); + } + return m.generate_samples(shots); +} + +void LightningSimulator::Sample(DataView &samples, size_t shots) { + auto li_samples = this->GenerateSamples(shots); + + RT_FAIL_IF(samples.size() != li_samples.size(), + "Invalid size for the pre-allocated samples"); + + const size_t numQubits = this->GetNumQubits(); + + // The lightning samples are layed out as a single vector of size + // shots*qubits, where each element represents a single bit. The + // corresponding shape is (shots, qubits). Gather the desired bits + // corresponding to the input wires into a bitstring. + auto samplesIter = samples.begin(); + for (size_t shot = 0; shot < shots; shot++) { + for (size_t wire = 0; wire < numQubits; wire++) { + *(samplesIter++) = + static_cast(li_samples[shot * numQubits + wire]); + } + } +} + +void LightningSimulator::PartialSample(DataView &samples, + const std::vector &wires, + size_t shots) { + const size_t numWires = wires.size(); + const size_t numQubits = this->GetNumQubits(); + + RT_FAIL_IF(numWires > numQubits, "Invalid number of wires"); + RT_FAIL_IF(!isValidQubits(wires), "Invalid given wires to measure"); + RT_FAIL_IF(samples.size() != shots * numWires, + "Invalid size for the pre-allocated partial-samples"); + + // get device wires + auto &&dev_wires = getDeviceWires(wires); + + auto li_samples = this->GenerateSamples(shots); + + // The lightning samples are layed out as a single vector of size + // shots*qubits, where each element represents a single bit. The + // corresponding shape is (shots, qubits). Gather the desired bits + // corresponding to the input wires into a bitstring. + auto samplesIter = samples.begin(); + for (size_t shot = 0; shot < shots; shot++) { + for (auto wire : dev_wires) { + *(samplesIter++) = + static_cast(li_samples[shot * numQubits + wire]); + } + } +} + +void LightningSimulator::Counts(DataView &eigvals, + DataView &counts, size_t shots) { + const size_t numQubits = this->GetNumQubits(); + const size_t numElements = 1U << numQubits; + + RT_FAIL_IF(eigvals.size() != numElements || counts.size() != numElements, + "Invalid size for the pre-allocated counts"); + + auto li_samples = this->GenerateSamples(shots); + + // Fill the eigenvalues with the integer representation of the corresponding + // computational basis bitstring. In the future, eigenvalues can also be + // obtained from an observable, hence the bitstring integer is stored as a + // double. + std::iota(eigvals.begin(), eigvals.end(), 0); + std::fill(counts.begin(), counts.end(), 0); + + // The lightning samples are layed out as a single vector of size + // shots*qubits, where each element represents a single bit. The + // corresponding shape is (shots, qubits). Gather the bits of all qubits + // into a bitstring. + for (size_t shot = 0; shot < shots; shot++) { + std::bitset basisState; + size_t idx = numQubits; + for (size_t wire = 0; wire < numQubits; wire++) { + basisState[--idx] = (li_samples[shot * numQubits + wire] != 0U); + } + counts(static_cast(basisState.to_ulong())) += 1; + } +} + +void LightningSimulator::PartialCounts(DataView &eigvals, + DataView &counts, + const std::vector &wires, + size_t shots) { + const size_t numWires = wires.size(); + const size_t numQubits = this->GetNumQubits(); + const size_t numElements = 1U << numWires; + + RT_FAIL_IF(numWires > numQubits, "Invalid number of wires"); + RT_FAIL_IF(!isValidQubits(wires), "Invalid given wires to measure"); + RT_FAIL_IF((eigvals.size() != numElements || counts.size() != numElements), + "Invalid size for the pre-allocated partial-counts"); + + // get device wires + auto &&dev_wires = getDeviceWires(wires); + + auto li_samples = this->GenerateSamples(shots); + + // Fill the eigenvalues with the integer representation of the corresponding + // computational basis bitstring. In the future, eigenvalues can also be + // obtained from an observable, hence the bitstring integer is stored as a + // double. + std::iota(eigvals.begin(), eigvals.end(), 0); + std::fill(counts.begin(), counts.end(), 0); + + // The lightning samples are layed out as a single vector of size + // shots*qubits, where each element represents a single bit. The + // corresponding shape is (shots, qubits). Gather the desired bits + // corresponding to the input wires into a bitstring. + for (size_t shot = 0; shot < shots; shot++) { + std::bitset basisState; + size_t idx = dev_wires.size(); + for (auto wire : dev_wires) { + basisState[--idx] = (li_samples[shot * numQubits + wire] != 0U); + } + counts(static_cast(basisState.to_ulong())) += 1; + } +} + +auto LightningSimulator::Measure(QubitIdType wire, + std::optional postselect) -> Result { + // get a measurement + std::vector wires = {reinterpret_cast(wire)}; + + std::vector probs(1U << wires.size()); + DataView buffer_view(probs); + auto device_shots = GetDeviceShots(); + SetDeviceShots(0); + PartialProbs(buffer_view, wires); + SetDeviceShots(device_shots); + + // It represents the measured result, true for 1, false for 0 + bool mres = Lightning::simulateDraw(probs, postselect, this->gen); + auto dev_wires = getDeviceWires(wires); + this->device_sv->collapse(dev_wires[0], mres); + return mres ? this->One() : this->Zero(); +} + +// Gradient +void LightningSimulator::Gradient(std::vector> &gradients, + const std::vector &trainParams) { + const bool tp_empty = trainParams.empty(); + const size_t num_observables = this->cache_manager.getNumObservables(); + const size_t num_params = this->cache_manager.getNumParams(); + const size_t num_train_params = tp_empty ? num_params : trainParams.size(); + const size_t jac_size = + num_train_params * this->cache_manager.getNumObservables(); + + if (jac_size == 0U) { + return; + } + + RT_FAIL_IF(gradients.size() != num_observables, + "Invalid number of pre-allocated gradients"); + + auto &&obs_callees = this->cache_manager.getObservablesCallees(); + bool is_valid_measurements = + std::all_of(obs_callees.begin(), obs_callees.end(), + [](const auto &m) { return m == MeasurementsT::Expval; }); + RT_FAIL_IF( + !is_valid_measurements, + "Unsupported measurements to compute gradient; " + "Adjoint differentiation method only supports expectation return type"); + + auto &&state = this->device_sv->getDataVector(); + + // create OpsData + auto &&ops_names = this->cache_manager.getOperationsNames(); + auto &&ops_params = this->cache_manager.getOperationsParameters(); + auto &&ops_wires = this->cache_manager.getOperationsWires(); + auto &&ops_inverses = this->cache_manager.getOperationsInverses(); + auto &&ops_matrices = this->cache_manager.getOperationsMatrices(); + auto &&ops_controlled_wires = + this->cache_manager.getOperationsControlledWires(); + auto &&ops_controlled_values = + this->cache_manager.getOperationsControlledValues(); + + const auto &&ops = Pennylane::Algorithms::OpsData( + ops_names, ops_params, ops_wires, ops_inverses, ops_matrices, + ops_controlled_wires, ops_controlled_values); + + // create the vector of observables + auto &&obs_keys = this->cache_manager.getObservablesKeys(); + std::vector< + std::shared_ptr>> + obs_vec; + obs_vec.reserve(obs_keys.size()); + for (auto idx : obs_keys) { + obs_vec.emplace_back(this->obs_manager.getObservable(idx)); + } + + std::vector all_params; + if (tp_empty) { + all_params.reserve(num_params); + for (size_t i = 0; i < num_params; i++) { + all_params.push_back(i); + } + } + + // construct the Jacobian data + Pennylane::Algorithms::JacobianData tape{ + num_params, state.size(), state.data(), + obs_vec, ops, tp_empty ? all_params : trainParams}; + + Pennylane::LightningQubit::Algorithms::AdjointJacobian adj; + std::vector jacobian(jac_size, 0); + adj.adjointJacobian(std::span{jacobian}, tape, + /* ref_data */ *this->device_sv, + /* apply_operations */ false); + + // convert jacobians to a list of lists for each observable + std::vector jacobian_t = Pennylane::LightningQubit::Util::Transpose( + jacobian, num_train_params, num_observables); + + std::vector cur_buffer(num_train_params); + auto begin_loc_iter = jacobian_t.begin(); + for (size_t obs_idx = 0; obs_idx < num_observables; obs_idx++) { + RT_ASSERT(begin_loc_iter != jacobian_t.end()); + RT_ASSERT(num_train_params <= gradients[obs_idx].size()); + std::move(begin_loc_iter, begin_loc_iter + num_train_params, + cur_buffer.begin()); + std::move(cur_buffer.begin(), cur_buffer.end(), + gradients[obs_idx].begin()); + begin_loc_iter += num_train_params; + } +} + +} // namespace Catalyst::Runtime::Simulator + +GENERATE_DEVICE_FACTORY(LightningSimulator, + Catalyst::Runtime::Simulator::LightningSimulator); diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/LightningSimulator.hpp b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/LightningSimulator.hpp new file mode 100644 index 000000000..948f78147 --- /dev/null +++ b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/LightningSimulator.hpp @@ -0,0 +1,128 @@ +// Copyright 2022 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 +#include +#include +#include +#include +#include +#include + +#include "StateVectorLQubitDynamic.hpp" + +#include "CacheManager.hpp" +#include "Exception.hpp" +#include "LightningObsManager.hpp" +#include "QuantumDevice.hpp" +#include "QubitManager.hpp" +#include "Utils.hpp" + +namespace Catalyst::Runtime::Simulator { +class LightningSimulator final : public Catalyst::Runtime::QuantumDevice { + private: + using StateVectorT = + Pennylane::LightningQubit::StateVectorLQubitDynamic; + + // static constants for RESULT values + static constexpr bool GLOBAL_RESULT_TRUE_CONST = true; + static constexpr bool GLOBAL_RESULT_FALSE_CONST = false; + + static constexpr size_t default_num_burnin{ + 100}; // tidy: readability-magic-numbers + static constexpr std::string_view default_kernel_name{ + "Local"}; // tidy: readability-magic-numbers + + Catalyst::Runtime::QubitManager qubit_manager{}; + Catalyst::Runtime::CacheManager> cache_manager{}; + bool tape_recording{false}; + size_t device_shots; + + std::mt19937 *gen = nullptr; + + bool mcmc{false}; + size_t num_burnin{0}; + std::string kernel_name; + + std::unique_ptr device_sv = std::make_unique(0); + LightningObsManager obs_manager{}; + + inline auto isValidQubit(QubitIdType wire) -> bool { + return this->qubit_manager.isValidQubitId(wire); + } + + inline auto isValidQubits(const std::vector &wires) -> bool { + return std::all_of(wires.begin(), wires.end(), [this](QubitIdType w) { + return this->isValidQubit(w); + }); + } + + inline auto isValidQubits(size_t numWires, const QubitIdType *wires) + -> bool { + return std::all_of(wires, wires + numWires, [this](QubitIdType w) { + return this->isValidQubit(w); + }); + } + + inline auto getDeviceWires(const std::vector &wires) + -> std::vector { + std::vector res; + res.reserve(wires.size()); + std::transform( + wires.begin(), wires.end(), std::back_inserter(res), + [this](auto w) { return this->qubit_manager.getDeviceId(w); }); + return res; + } + + public: + explicit LightningSimulator( + const std::string &kwargs = "{}") // NOLINT(hicpp-member-init) + { + auto &&args = Catalyst::Runtime::parse_kwargs(kwargs); + device_shots = args.contains("shots") + ? static_cast(std::stoll(args["shots"])) + : 0; + mcmc = args.contains("mcmc") ? args["mcmc"] == "True" : false; + num_burnin = args.contains("num_burnin") + ? static_cast(std::stoll(args["num_burnin"])) + : default_num_burnin; + kernel_name = args.contains("kernel_name") ? args["kernel_name"] + : default_kernel_name; + } + ~LightningSimulator() override = default; + + void SetDevicePRNG(std::mt19937 *gen) override; + void SetState(DataView, 1> &state, + std::vector &wires) override; + void SetBasisState(DataView &n, + std::vector &wires) override; + + QUANTUM_DEVICE_DEL_DECLARATIONS(LightningSimulator); + + // TODO: properly refactor the common device methods, + // instead of using #define macros in + // runtime/lib/backend/common/Utils.hpp + // When done, remove the NOLINT(hicpp-member-init) on the class constructor + QUANTUM_DEVICE_RT_DECLARATIONS; + QUANTUM_DEVICE_QIS_DECLARATIONS; + + auto CacheManagerInfo() + -> std::tuple, + std::vector>; + auto GenerateSamplesMetropolis(size_t shots) -> std::vector; + auto GenerateSamples(size_t shots) -> std::vector; +}; +} // namespace Catalyst::Runtime::Simulator diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/StateVectorLQubitDynamic.cpp b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/StateVectorLQubitDynamic.cpp new file mode 100644 index 000000000..4516e67c4 --- /dev/null +++ b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/StateVectorLQubitDynamic.cpp @@ -0,0 +1,19 @@ +// Copyright 2022 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. + +#include "StateVectorLQubitDynamic.hpp" + +// explicit instantiation +template class Pennylane::LightningQubit::StateVectorLQubitDynamic; +template class Pennylane::LightningQubit::StateVectorLQubitDynamic; diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/StateVectorLQubitDynamic.hpp b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/StateVectorLQubitDynamic.hpp new file mode 100644 index 000000000..dc967c235 --- /dev/null +++ b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/StateVectorLQubitDynamic.hpp @@ -0,0 +1,388 @@ +// Copyright 2022 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 +#include +#include + +#include "BitUtil.hpp" +#include "LinearAlgebra.hpp" +#include "Util.hpp" + +#include "Error.hpp" + +#include + +#include + +namespace { +using Pennylane::Util::AlignedAllocator; +using Pennylane::Util::bestCPUMemoryModel; +using Pennylane::Util::exp2; +using Pennylane::Util::isPerfectPowerOf2; +using Pennylane::Util::log2PerfectPower; +using Pennylane::Util::ONE; +using Pennylane::Util::squaredNorm; +using Pennylane::Util::ZERO; +} // namespace + +namespace Pennylane::LightningQubit { +/** + * @brief State-vector dynamic class. + * + * This class allocates and deallocates qubits/wires dynamically, + * and defines all operations to manipulate the statevector data for + * quantum circuit simulation. + * + */ +template +class StateVectorLQubitDynamic + : public StateVectorLQubit> { + public: + using PrecisionT = fp_t; + using ComplexT = std::complex; + using CFP_t = ComplexT; + using MemoryStorageT = Pennylane::Util::MemoryStorageLocation::Internal; + using AlignedAllocComplexT = AlignedAllocator; + + private: + using BaseType = + StateVectorLQubit>; + std::vector> data_; + + static constexpr PrecisionT epsilon_ = + std::numeric_limits::epsilon() * 100; + + template + inline OIter _move_data_elements(IIter first, size_t distance, + OIter second) { + *second++ = std::move(*first); + for (size_t i = 1; i < distance; i++) { + *second++ = std::move(*++first); + } + return second; + } + + template + inline OIter _shallow_move_data_elements(IIter first, size_t distance, + OIter second) { + for (size_t i = 0; i < distance; i++) { + *second++ = std::move(*first); + *first = ZERO(); + first++; + } + return second; + } + + inline void + _scalar_mul_data(std::vector &data, + ComplexT scalar) { + std::transform( + data.begin(), data.end(), data.begin(), + [scalar](const ComplexT &elem) { return elem * scalar; }); + } + + inline void + _normalize_data(std::vector &data) { + _scalar_mul_data(data, + ONE() / + std::sqrt(squaredNorm(data.data(), data.size()))); + } + + public: + /** + * @brief Create a new statevector + * + * @param num_qubits Number of qubits + * @param threading Threading option the statevector to use + * @param memory_model Memory model the statevector will use + */ + explicit StateVectorLQubitDynamic( + size_t num_qubits, Threading threading = Threading::SingleThread, + CPUMemoryModel memory_model = bestCPUMemoryModel()) + : BaseType{num_qubits, threading, memory_model}, + data_{exp2(num_qubits), ZERO(), + getAllocator( // LCOV_EXCL_LINE + this->memory_model_)} { + data_[0] = ONE(); + } + + /** + * @brief Construct a statevector from another statevector + * + * @tparam OtherDerived A derived type of StateVectorLQubit to use for + * construction. + * @param other Another statevector to construct the statevector from + */ + template + explicit StateVectorLQubitDynamic( + const StateVectorLQubit &other) + : BaseType(other.getNumQubits(), other.threading(), + other.memoryModel()), + data_{other.getData(), other.getData() + other.getLength(), + getAllocator(this->memory_model_)} {} + + /** + * @brief Construct a statevector from data pointer + * + * @param other_data Data pointer to construct the statvector from. + * @param other_size Size of the data + * @param threading Threading option the statevector to use + * @param memory_model Memory model the statevector will use + */ + StateVectorLQubitDynamic(const ComplexT *other_data, size_t other_size, + Threading threading = Threading::SingleThread, + CPUMemoryModel memory_model = bestCPUMemoryModel()) + : BaseType(log2PerfectPower(other_size), threading, memory_model), + data_{other_data, other_data + other_size, + getAllocator(this->memory_model_)} { + PL_ABORT_IF_NOT(isPerfectPowerOf2(other_size), + "The size of provided data must be a power of 2."); + } + + /** + * @brief Construct a statevector from a data vector + * + * @tparam Alloc Allocator type of std::vector to use for constructing + * statevector. + * @param other Data to construct the statevector from + * @param threading Threading option the statevector to use + * @param memory_model Memory model the statevector will use + */ + template + explicit StateVectorLQubitDynamic( + const std::vector &other, + Threading threading = Threading::SingleThread, + CPUMemoryModel memory_model = bestCPUMemoryModel()) + : StateVectorLQubitDynamic(other.data(), other.size(), threading, + memory_model) {} + + StateVectorLQubitDynamic(const StateVectorLQubitDynamic &rhs) = default; + StateVectorLQubitDynamic(StateVectorLQubitDynamic &&) noexcept = default; + + StateVectorLQubitDynamic & + operator=(const StateVectorLQubitDynamic &) = default; + StateVectorLQubitDynamic & + operator=(StateVectorLQubitDynamic &&) noexcept = default; + + ~StateVectorLQubitDynamic() = default; + + /** + * @brief Set the number of qubits represented by the statevector data. + * + * @return std::size_t + */ + void setNumQubits(size_t qubits) noexcept { this->num_qubits_ = qubits; } + + /** + * @brief Get underlying C-style data of the state-vector. + */ + [[nodiscard]] auto getData() -> ComplexT * { return data_.data(); } + + /** + * @brief Get underlying C-style data of the state-vector. + */ + [[nodiscard]] auto getData() const -> const ComplexT * { + return data_.data(); + } + + /** + * @brief Get underlying data vector. + */ + [[nodiscard]] auto getDataVector() + -> std::vector & { + return data_; + } + + /** + * @brief Get underlying data vector. + */ + [[nodiscard]] auto getDataVector() const + -> const std::vector & { + return data_; + } + + /** + * @brief Update data of the class to new_data + * + * @param new_data data pointer to new data. + * @param new_size size of underlying data storage. + */ + void updateData(const ComplexT *new_data, size_t new_size) { + PL_ASSERT(data_.size() == new_size); + std::copy(new_data, new_data + new_size, data_.data()); + } + + /** + * @brief Update data of the class to new_data + * + * @tparam Alloc Allocator type of std::vector to use for updating data. + * @param new_data std::vector contains data. + */ + template + void updateData(const std::vector &new_data) { + updateData(new_data.data(), new_data.size()); + } + + [[nodiscard]] AlignedAllocComplexT allocator() const { + return data_.get_allocator(); + } + + [[nodiscard]] auto isValidWire(size_t wire) -> bool { + return wire < this->getNumQubits(); + } + + /** + * @brief Compute the purity of the system after releasing (a qubit) `wire`. + * + * This traces out the complement of the wire for a more efficient + * computation of the purity in O(N) with calculating the reduced density + * matrix after tracing out the complement of qubit `wire`. + * + * @param wire Index of the wire. + * @return ComplexT + */ + auto getSubsystemPurity(size_t wire) -> ComplexT { + PL_ABORT_IF_NOT(isValidWire(wire), + "Invalid wire: The wire must be in the range of wires"); + + const size_t sv_size = data_.size(); + + // With `k` indexing the subsystem on n-1 qubits, we need to insert an + // addtional bit into the index of the full state-vector at position + // `wire`. These masks enable us to split the bits of the index `k` into + // those above and below `wire`. + const size_t lower_mask = (1UL << wire) - 1; + const size_t upper_mask = sv_size - lower_mask - 1; + + // The resulting 2x2 reduced density matrix of the complement system to + // qubit `wire`. + std::vector rho(4, {0, 0}); + + for (uint8_t i = 0; i < 2; i++) { + for (uint8_t j = 0; j < 2; j++) { + ComplexT sum{0, 0}; + for (size_t k = 0; k < (sv_size / 2); k++) { + size_t idx_wire_0 = + (/* upper_bits: */ (upper_mask & k) << 1UL) + + /* lower_bits: */ (lower_mask & k); + size_t idx_i = idx_wire_0 + (i << wire); + size_t idx_j = idx_wire_0 + (j << wire); + + // This computes <00..i..00|psi> on the first + // iteration, with the last iteration computing + // <11..i..11|psi>. + sum += data_[idx_i] * std::conj(data_[idx_j]); + } + rho[2 * i + j] = sum; + } + } + + // Compute/Return the trace of rho**2 + return (rho[0] * rho[0]) + (ComplexT{2, 0} * rho[1] * rho[2]) + + (rho[3] * rho[3]); + } + + /** + * @brief Check the purity of a system after releasing/disabling `wire`. + * + * @param wire Index of the wire. + * @param eps The comparing precision threshold. + * @return bool + */ + [[nodiscard]] auto checkSubsystemPurity(size_t wire, double eps = epsilon_) + -> bool { + ComplexT purity = getSubsystemPurity(wire); + return (std::abs(1.0 - purity.real()) < eps) && (purity.imag() < eps); + } + + /** + * @brief Allocate a new wire. + * + * @return It updates the state-vector and the number of qubits, + * and returns index of the activated wire. + */ + auto allocateWire() -> size_t { + const size_t next_idx = this->getNumQubits(); + const size_t dsize = data_.size(); + data_.resize(dsize << 1UL); + + auto src = data_.begin(); + std::advance(src, dsize - 1); + for (auto dst = data_.end() - 2; src != data_.begin(); + std::advance(src, -1), std::advance(dst, -2)) { + _shallow_move_data_elements(src, 1, dst); + } + + this->setNumQubits(next_idx + 1); + return next_idx; + } + + /** + * @brief Release a given wire. + * + * @param wire Index of the wire to be released + */ + void releaseWire(size_t wire) { + PL_ABORT_IF_NOT( + checkSubsystemPurity(wire), + "Invalid wire: " + "The state-vector must remain pure after releasing a wire") + + const size_t distance = 1UL << wire; + + auto dst = data_.begin(); + + // Check if the reduced state-vector is the first-half + bool is_first_half = false; + for (auto src = dst; src < data_.end(); + std::advance(src, 2 * distance)) { + is_first_half = std::any_of( + src, src + static_cast(distance), + [](ComplexT &e) { return e != ZERO(); }); + if (is_first_half) { + break; + } + } + + auto src = dst; + if (!is_first_half) { + std::advance(src, distance); + } + + for (; src < data_.end(); + std::advance(src, 2 * distance), std::advance(dst, distance)) { + _move_data_elements(src, distance, dst); + } + data_.resize(data_.size() / 2); + // normalize state-vector + _normalize_data(data_); + + this->setNumQubits(this->getNumQubits() - 1); + } + + /** + * @brief Update the state-vector to the initial state with 0-qubit. + */ + void clearData() { + data_.clear(); + this->setNumQubits(0); + + // the init state-vector + data_.push_back(ONE()); + } +}; + +} // namespace Pennylane::LightningQubit diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/CMakeLists.txt b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/CMakeLists.txt new file mode 100644 index 000000000..51e7f0f1f --- /dev/null +++ b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/CMakeLists.txt @@ -0,0 +1,42 @@ +cmake_minimum_required(VERSION 3.20) + +project(lightning_qubit_catalyst_tests) + +# Default build type for test code is Debug +if(NOT CMAKE_BUILD_TYPE) + set(CMAKE_BUILD_TYPE Debug) +endif() + +include("${pennylane_lightning_SOURCE_DIR}/cmake/support_tests.cmake") +FetchAndIncludeCatch() + +################################################################################ +# Define library +################################################################################ + +add_library(lightning_qubit_catalyst_tests INTERFACE) +target_link_libraries(lightning_qubit_catalyst_tests INTERFACE Catch2::Catch2 + lightning_qubit_catalyst +) + +ProcessTestOptions(lightning_qubit_catalyst_tests) + +target_sources(lightning_qubit_catalyst_tests INTERFACE runner_lightning_qubit_catalyst.cpp) + +################################################################################ +# Define targets +################################################################################ +set(TEST_SOURCES Test_LightningSimulator.cpp + Test_LightningMeasures.cpp + Test_LightningGradient.cpp + Test_LightningDriver.cpp + Test_SVDynamicCPU_Allocation.cpp + Test_SVDynamicCPU_Core.cpp +) + +add_executable(lightning_qubit_catalyst_tests_runner ${TEST_SOURCES}) +target_link_libraries(lightning_qubit_catalyst_tests_runner PRIVATE lightning_qubit_catalyst_tests) + +catch_discover_tests(lightning_qubit_catalyst_tests_runner) + +install(TARGETS lightning_qubit_catalyst_tests_runner DESTINATION bin) diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/Test_LightningDriver.cpp b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/Test_LightningDriver.cpp new file mode 100644 index 000000000..961eff6e9 --- /dev/null +++ b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/Test_LightningDriver.cpp @@ -0,0 +1,267 @@ +// Copyright 2022 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. + +#include +#include + +#include "LightningSimulator.hpp" +#include "catch2/catch.hpp" + +using namespace Catalyst::Runtime; +using namespace Catalyst::Runtime::Simulator; +using LQSimulator = LightningSimulator; + +TEST_CASE("Test parse_kwargs coverage", "[Utils]") { + std::string case1; + CHECK(parse_kwargs(case1).empty()); + + std::string case2{"{shots : 1000}"}; + std::string case3{"shots : 1000"}; + std::string case4{"'shots':'1000'"}; + CHECK(parse_kwargs(case2) == parse_kwargs(case3)); + CHECK(parse_kwargs(case3) == parse_kwargs(case4)); + + std::string case5{"{'A':'B', 'C':'D', 'E':'F'}"}; + auto res5 = parse_kwargs(case5); + CHECK(res5.size() == 3); + CHECK((res5.contains("A") && res5["A"] == "B")); + CHECK((res5.contains("C") && res5["C"] == "D")); + CHECK((res5.contains("E") && res5["E"] == "F")); + + std::string case6{ + "device_type : braket.aws.qubit,{'shots': 0, 'device_arn': 'sv1', " + "'s3_destination_folder': \"('catalyst-op3-s3', 'prefix')\"}"}; + auto res6 = parse_kwargs(case6); + CHECK(res6.size() == 4); + CHECK((res6.contains("device_type") && + res6["device_type"] == "braket.aws.qubit")); + CHECK((res6.contains("shots") && res6["shots"] == "0")); + CHECK((res6.contains("device_arn") && res6["device_arn"] == "sv1")); + CHECK((res6.contains("s3_destination_folder") && + res6["s3_destination_folder"] == "('catalyst-op3-s3', 'prefix')")); +} + +TEST_CASE("lightning Basis vector", "[Driver]") { + std::unique_ptr sim = std::make_unique(); + + [[maybe_unused]] QubitIdType q1 = sim->AllocateQubit(); + [[maybe_unused]] QubitIdType q2 = sim->AllocateQubit(); + QubitIdType q3 = sim->AllocateQubit(); + + sim->ReleaseQubit(q3); + + std::vector> state(1U << sim->GetNumQubits()); + DataView, 1> view(state); + sim->State(view); + + CHECK(view(0).real() == Approx(1.0).epsilon(1e-5)); + CHECK(view(0).imag() == Approx(0.0).epsilon(1e-5)); + CHECK(view(1).real() == Approx(0.0).epsilon(1e-5)); + CHECK(view(1).imag() == Approx(0.0).epsilon(1e-5)); + CHECK(view(2).real() == Approx(0.0).epsilon(1e-5)); + CHECK(view(2).imag() == Approx(0.0).epsilon(1e-5)); + CHECK(view(3).real() == Approx(0.0).epsilon(1e-5)); + CHECK(view(3).imag() == Approx(0.0).epsilon(1e-5)); +} + +TEST_CASE("Qubit allocatation and deallocation", "[Driver]") { + std::unique_ptr sim = std::make_unique(); + + constexpr size_t n = 1; + constexpr size_t sz = (1UL << n); + + QubitIdType q; + for (size_t i = 0; i < n; i++) { + q = sim->AllocateQubit(); + } + + CHECK(n == static_cast(q) + 1); + + std::vector> state(1U << sim->GetNumQubits()); + DataView, 1> view(state); + sim->State(view); + + CHECK(state.size() == (1UL << n)); + CHECK(state[0].real() == Approx(1.0).epsilon(1e-5)); + CHECK(state[0].imag() == Approx(0.0).epsilon(1e-5)); + + std::complex sum{0, 0}; + for (size_t i = 1; i < sz; i++) { + sum += state[i]; + } + + CHECK(sum.real() == Approx(0.0).epsilon(1e-5)); + CHECK(sum.imag() == Approx(0.0).epsilon(1e-5)); + + for (size_t i = n; i > 0; i--) { + CHECK(state.size() == sz); + + sim->ReleaseQubit(i - 1); + sim->AllocateQubit(); + + DataView, 1> view(state); + sim->State(view); + } +} + +TEST_CASE("test AllocateQubits", "[Driver]") { + std::unique_ptr sim = std::make_unique(); + + CHECK(sim->AllocateQubits(0).empty()); + + auto &&q = sim->AllocateQubits(2); + + sim->ReleaseQubit(q[0]); + + std::vector> state(1U << sim->GetNumQubits()); + DataView, 1> view(state); + sim->State(view); + + CHECK(state[0].real() == Approx(1.0).epsilon(1e-5)); +} + +TEST_CASE("test multiple AllocateQubits", "[Driver]") { + std::unique_ptr sim = std::make_unique(); + + auto &&q1 = sim->AllocateQubits(2); + CHECK(q1[0] == 0); + CHECK(q1[1] == 1); + + auto &&q2 = sim->AllocateQubits(3); + CHECK(q2.size() == 3); + CHECK(q2[0] == 2); + CHECK(q2[2] == 4); +} + +TEST_CASE("test DeviceShots", "[Driver]") { + std::unique_ptr sim = std::make_unique(); + + CHECK(sim->GetDeviceShots() == 0); + + sim->SetDeviceShots(500); + + CHECK(sim->GetDeviceShots() == 500); +} + +TEST_CASE("compute register tests", "[Driver]") { + std::unique_ptr sim = std::make_unique(); + + constexpr size_t n = 10; + std::vector Qs; + Qs.reserve(n); + + // allocate a few qubits + for (size_t i = 0; i < n; i++) { + Qs[i] = sim->AllocateQubit(); + } + + // release some of them + sim->ReleaseQubit(n - 1); + sim->ReleaseQubit(n - 2); + + const size_t new_n = n - 2; + + // check the correctness + std::vector Qs_expected(new_n); + std::iota(Qs_expected.begin(), Qs_expected.end(), + static_cast(0)); + + for (size_t i = 0; i < new_n; i++) { + CHECK(Qs_expected[i] == Qs[i]); + } +} + +TEST_CASE("Check an unsupported operation", "[Driver]") { + REQUIRE_THROWS_WITH( + Lightning::lookup_gates(Lightning::simulator_gate_info, + "UnsupportedGateName"), + Catch::Contains( + "The given operation is not supported by the simulator")); +} + +TEST_CASE("QuantumDevice object test [lightning.qubit]", "[Driver]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr size_t n = 10; + std::vector Qs; + Qs.reserve(n); + for (size_t i = 0; i < n; i++) { + Qs[i] = sim->AllocateQubit(); + } + + sim->NamedOperation("Identity", {}, {Qs[0]}, false); + sim->NamedOperation("Identity", {}, {Qs[2]}, false); + sim->NamedOperation("Identity", {}, {Qs[4]}, false); + sim->NamedOperation("Identity", {}, {Qs[6]}, false); + sim->NamedOperation("Identity", {}, {Qs[8]}, false); + + std::vector> state(1U << sim->GetNumQubits()); + DataView, 1> view(state); + sim->State(view); + + CHECK(state[0].real() == Approx(1.0).epsilon(1e-5)); + CHECK(state[0].imag() == Approx(0.0).epsilon(1e-5)); + + std::complex sum{0, 0}; + for (size_t i = 1; i < state.size(); i++) { + sum += state[i]; + } + + CHECK(sum.real() == Approx(0.0).epsilon(1e-5)); + CHECK(sum.imag() == Approx(0.0).epsilon(1e-5)); + + for (size_t i = 0; i < n; i++) { + sim->ReleaseQubit(static_cast(i)); + // 0, 1, 2, ..., 9 + } + + for (size_t i = 10; i < n + 10; i++) { + CHECK(static_cast(i) == sim->AllocateQubit()); + // 10, 11, ..., 19 + } + + for (size_t i = 10; i < n + 10; i++) { + sim->ReleaseQubit(static_cast(i)); + // 10, 11, ..., 19 + } + + for (size_t i = 20; i < n + 20; i++) { + CHECK(static_cast(i) == sim->AllocateQubit()); + // 20, 21, ..., 29 + } +} + +TEST_CASE("Check re-AllocateQubit", "[Driver]") { + std::unique_ptr sim = std::make_unique(); + + sim->AllocateQubit(); + sim->NamedOperation("Hadamard", {}, {0}, false); + + std::vector> state(1U << sim->GetNumQubits()); + DataView, 1> view(state); + sim->State(view); + CHECK(state[0].real() == Approx(0.707107).epsilon(1e-5)); + CHECK(state[1].real() == Approx(0.707107).epsilon(1e-5)); + + sim->AllocateQubit(); + sim->AllocateQubit(); + sim->AllocateQubit(); + + state = std::vector>(1U << sim->GetNumQubits()); + view = DataView, 1>(state); + sim->State(view); + CHECK(state[0].real() == Approx(0.707107).epsilon(1e-5)); + CHECK(state[8].real() == Approx(0.707107).epsilon(1e-5)); +} diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/Test_LightningGradient.cpp b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/Test_LightningGradient.cpp new file mode 100644 index 000000000..87e9cd6bf --- /dev/null +++ b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/Test_LightningGradient.cpp @@ -0,0 +1,298 @@ +// Copyright 2022 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. + +#include "LightningSimulator.hpp" +#include "catch2/catch.hpp" + +/// @cond DEV +namespace { +// MemRef type definition (Helper) +template struct MemRefT { + T *data_allocated; + T *data_aligned; + std::size_t offset; + std::size_t sizes[R]; + std::size_t strides[R]; +}; +using namespace Catalyst::Runtime::Simulator; +using LQSimulator = LightningSimulator; +} // namespace +/// @endcond + +TEST_CASE("Zero qubits. Zero parameters", "[Gradient]") { + std::unique_ptr LQsim = std::make_unique(); + + std::vector> gradients; + std::vector Qs = LQsim->AllocateQubits(0); + REQUIRE_NOTHROW(LQsim->Gradient(gradients, {})); +} + +TEST_CASE("Test Gradient with zero number of obs", "[Gradient]") { + std::unique_ptr sim = std::make_unique(); + + std::vector buffer(1); + std::vector> gradients; + gradients.emplace_back(buffer); + + const std::vector trainParams{0}; + + const auto q = sim->AllocateQubit(); + + sim->StartTapeRecording(); + + sim->NamedOperation("S", {}, {q}, false); + sim->NamedOperation("T", {}, {q}, false); + + REQUIRE_NOTHROW(sim->Gradient(gradients, trainParams)); + + sim->StopTapeRecording(); +} + +TEST_CASE("Test Gradient with Var", "[Gradient]") { + std::unique_ptr sim = std::make_unique(); + + std::vector buffer(1); + std::vector> gradients; + gradients.emplace_back(buffer); + + const std::vector trainParams{0}; + + const auto q = sim->AllocateQubit(); + + sim->StartTapeRecording(); + + sim->NamedOperation("RX", {-M_PI / 7}, {q}, false); + auto pz = sim->Observable(ObsId::PauliZ, {}, {q}); + sim->Var(pz); + + REQUIRE_THROWS_WITH( + sim->Gradient(gradients, trainParams), + Catch::Contains("Unsupported measurements to compute gradient")); + + REQUIRE_THROWS_WITH( + sim->Gradient(gradients, {}), + Catch::Contains("Unsupported measurements to compute gradient")); + + sim->StopTapeRecording(); +} + +TEST_CASE("Test Gradient with Op=RX, Obs=Z", "[Gradient]") { + std::unique_ptr sim = std::make_unique(); + + std::vector buffer(1); + std::vector> gradients; + gradients.emplace_back(buffer); + + const std::vector trainParams{0}; + + const auto q = sim->AllocateQubit(); + + sim->StartTapeRecording(); + + sim->NamedOperation("RX", {-M_PI / 7}, {q}, false); + auto obs = sim->Observable(ObsId::PauliZ, {}, {q}); + sim->Expval(obs); + + sim->Gradient(gradients, trainParams); + CHECK(-sin(-M_PI / 7) == Approx(buffer[0]).margin(1e-5)); + + // Update buffer + buffer[0] = 0.0; + + sim->Gradient(gradients, {}); + CHECK(-sin(-M_PI / 7) == Approx(buffer[0]).margin(1e-5)); + + sim->StopTapeRecording(); +} + +TEST_CASE("Test Gradient with Op=RX, Obs=Hermitian", "[Gradient]") { + std::unique_ptr sim = std::make_unique(); + + std::vector buffer(1); + std::vector> gradients; + gradients.emplace_back(buffer); + + const std::vector trainParams{0}; + + constexpr double expected{0.2169418696}; + + const auto q = sim->AllocateQubit(); + + sim->StartTapeRecording(); + + sim->NamedOperation("RX", {-M_PI / 7}, {q}, false); + + std::vector> mat{ + {1.0, 0.0}, {0.0, 0.0}, {2.0, 0.0}, {0.0, 0.0}}; + + auto obs = sim->Observable(ObsId::Hermitian, mat, {q}); + + sim->Expval(obs); + + sim->Gradient(gradients, trainParams); + CHECK(expected == Approx(buffer[0]).margin(1e-5)); + + // Update buffer + buffer[0] = 0.0; + + sim->Gradient(gradients, {}); + CHECK(expected == Approx(buffer[0]).margin(1e-5)); + + sim->StopTapeRecording(); +} + +TEST_CASE("Test Gradient with Op=[RX,RX,RX,CZ], Obs=[Z,Z,Z]", "[Gradient]") { + std::unique_ptr sim = std::make_unique(); + + constexpr std::size_t num_parms = 3; + + std::vector buffer_p0(num_parms); + std::vector buffer_p1(num_parms); + std::vector buffer_p2(num_parms); + std::vector> gradients; + gradients.emplace_back(buffer_p0); + gradients.emplace_back(buffer_p1); + gradients.emplace_back(buffer_p2); + + const std::vector trainParams{0, 1, 2}; + + const std::vector param{-M_PI / 7, M_PI / 5, 2 * M_PI / 3}; + const std::vector expected{-sin(param[0]), -sin(param[1]), + -sin(param[2])}; + + const auto Qs = sim->AllocateQubits(num_parms); + + sim->StartTapeRecording(); + + sim->NamedOperation("RX", {param[0]}, {Qs[0]}, false); + sim->NamedOperation("RX", {param[1]}, {Qs[1]}, false); + sim->NamedOperation("RX", {param[2]}, {Qs[2]}, false); + sim->NamedOperation("CZ", {}, {Qs[0], Qs[2]}, false); + + std::vector> mat{ + {1.0, 0.0}, {0.0, 0.0}, {2.0, 0.0}, {0.0, 0.0}}; + + auto obs0 = sim->Observable(ObsId::PauliZ, {}, {Qs[0]}); + auto obs1 = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + auto obs2 = sim->Observable(ObsId::PauliZ, {}, {Qs[2]}); + + sim->Expval(obs0); + sim->Expval(obs1); + sim->Expval(obs2); + + sim->Gradient(gradients, trainParams); + CHECK(expected[0] == Approx(buffer_p0[0]).margin(1e-5)); + CHECK(expected[1] == Approx(buffer_p1[1]).margin(1e-5)); + CHECK(expected[2] == Approx(buffer_p2[2]).margin(1e-5)); + + sim->StopTapeRecording(); +} + +TEST_CASE("Test Gradient with Op=Mixed, Obs=Hamiltonian([Z@Z, H], {0.2, 0.6})", + "[Gradient]") { + std::unique_ptr sim = std::make_unique(); + + constexpr std::size_t num_parms = 6; + + std::vector buffer(num_parms); + std::vector> gradients; + gradients.emplace_back(buffer); + + const std::vector trainParams{0, 1, 2, 3, 4, 5}; + + const std::vector param{-M_PI / 7, M_PI / 5, 2 * M_PI / 3}; + const std::vector expected{0.0, -0.2493761627, 0.0, + 0.0, -0.1175570505, 0.0}; + + const auto Qs = sim->AllocateQubits(3); + + sim->StartTapeRecording(); + + sim->NamedOperation("RZ", {param[0]}, {Qs[0]}, false); + sim->NamedOperation("RY", {param[1]}, {Qs[0]}, false); + sim->NamedOperation("RZ", {param[2]}, {Qs[0]}, false); + sim->NamedOperation("CNOT", {}, {Qs[0], Qs[1]}, false); + sim->NamedOperation("CNOT", {}, {Qs[1], Qs[2]}, false); + sim->NamedOperation("RZ", {param[0]}, {Qs[1]}, false); + sim->NamedOperation("RY", {param[1]}, {Qs[1]}, false); + sim->NamedOperation("RZ", {param[2]}, {Qs[1]}, false); + + std::vector> mat{ + {1.0, 0.0}, {0.0, 0.0}, {2.0, 0.0}, {0.0, 0.0}}; + + auto obs0 = sim->Observable(ObsId::PauliZ, {}, {Qs[0]}); + auto obs1 = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + auto obs2 = sim->TensorObservable({obs0, obs1}); + auto obs3 = sim->Observable(ObsId::Hadamard, {}, {Qs[2]}); + auto obs4 = sim->HamiltonianObservable({0.2, 0.6}, {obs2, obs3}); + + sim->Expval(obs4); + + sim->Gradient(gradients, trainParams); + + for (std::size_t i = 0; i < num_parms; i++) { + CAPTURE(i); + CHECK(expected[i] == Approx(buffer[i]).margin(1e-5)); + buffer[i] = 0.0; + } + + sim->Gradient(gradients, {}); + + for (std::size_t i = 0; i < num_parms; i++) { + CAPTURE(i); + CHECK(expected[i] == Approx(buffer[i]).margin(1e-5)); + } + + sim->StopTapeRecording(); +} + +TEST_CASE("Test Gradient with QubitUnitary", "[Gradient]") { + std::unique_ptr sim = std::make_unique(); + + std::vector buffer(1); + std::vector> gradients; + gradients.emplace_back(buffer); + + const std::vector trainParams{0}; + + constexpr double expected{-0.8611041863}; + + const std::vector> matrix{ + {-0.6709485262524046, -0.6304426335363695}, + {-0.14885403153998722, 0.3608498832392019}, + {-0.2376311670004963, 0.3096798175687841}, + {-0.8818365947322423, -0.26456390390903695}, + }; + + const auto Qs = sim->AllocateQubits(1); + + sim->StartTapeRecording(); + + sim->NamedOperation("RX", {-M_PI / 7}, {Qs[0]}, false); + sim->MatrixOperation(matrix, {Qs[0]}, false); + + auto obs = sim->Observable(ObsId::PauliY, {}, {Qs[0]}); + sim->Expval(obs); + + sim->Gradient(gradients, trainParams); + CHECK(expected == Approx(buffer[0]).margin(1e-5)); + + // Update buffer + buffer[0] = 0.0; + + sim->Gradient(gradients, {}); + CHECK(expected == Approx(buffer[0]).margin(1e-5)); + + sim->StopTapeRecording(); +} diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/Test_LightningMeasures.cpp b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/Test_LightningMeasures.cpp new file mode 100644 index 000000000..3189fbef7 --- /dev/null +++ b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/Test_LightningMeasures.cpp @@ -0,0 +1,1833 @@ +// Copyright 2022 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. + +#include + +#include "CacheManager.hpp" +#include "LightningSimulator.hpp" +#include "QuantumDevice.hpp" +#include "Types.h" +#include "Utils.hpp" +#include "catch2/catch.hpp" +#include "cmath" + +/// @cond DEV +namespace { +// MemRef type definition (Helper) +template struct MemRefT { + T *data_allocated; + T *data_aligned; + std::size_t offset; + std::size_t sizes[R]; + std::size_t strides[R]; +}; +using namespace Catalyst::Runtime::Simulator; +using LQSimulator = LightningSimulator; +} // namespace +/// @endcond + +TEST_CASE("NameObs test with invalid number of wires", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + REQUIRE_THROWS_WITH(sim->Observable(ObsId::PauliX, {}, {1}), + Catch::Contains("Invalid number of wires")); +} + +TEST_CASE("NameObs test with invalid given wires for NamedObs", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + sim->AllocateQubit(); + + REQUIRE_THROWS_WITH(sim->Observable(ObsId::PauliX, {}, {1}), + Catch::Contains("Invalid given wires")); +} + +TEST_CASE("HermitianObs test with invalid number of wires", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + REQUIRE_THROWS_WITH(sim->Observable(ObsId::Hermitian, {}, {1}), + Catch::Contains("Invalid number of wires")); +} + +TEST_CASE("HermitianObs test with invalid given wires for HermitianObs", + "[Measures]") { + std::unique_ptr sim = std::make_unique(); + sim->AllocateQubit(); + + REQUIRE_THROWS_WITH(sim->Observable(ObsId::Hermitian, {}, {1}), + Catch::Contains("Invalid given wires")); +} + +TEST_CASE("Check an unsupported observable", "[Measures]") { + REQUIRE_THROWS_WITH( + Lightning::lookup_obs( + Lightning::simulator_observable_support, static_cast(10)), + Catch::Contains( + "The given observable is not supported by the simulator")); +} + +TEST_CASE("Measurement collapse test with 2 wires", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + constexpr std::size_t n = 2; + std::vector Qs = sim->AllocateQubits(n); + + sim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + auto *m = sim->Measure(Qs[0]); + std::vector> state(1U << sim->GetNumQubits()); + DataView, 1> view(state); + sim->State(view); + + // LCOV_EXCL_START + // This is conditional over the measurement result + if (*m) { + CHECK(pow(std::abs(std::real(state[2])), 2) + + pow(std::abs(std::imag(state[2])), 2) == + Approx(1.0).margin(1e-5)); + } else { + CHECK(pow(std::abs(std::real(state[0])), 2) + + pow(std::abs(std::imag(state[0])), 2) == + Approx(1.0).margin(1e-5)); + } + // LCOV_EXCL_STOP +} + +TEST_CASE("Measurement collapse concrete logical qubit difference", + "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + constexpr std::size_t n = 1; + // The first time an array is allocated, logical and concrete qubits + // are the same. + std::vector Qs = sim->AllocateQubits(n); + sim->ReleaseAllQubits(); + + // Now in this the concrete qubits are shifted by n. + Qs = sim->AllocateQubits(n); + + sim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + sim->Measure(Qs[0]); + std::vector> state(1U << sim->GetNumQubits()); + DataView, 1> view(state); + sim->State(view); + + // LCOV_EXCL_START + bool is_zero = pow(std::abs(std::real(state[0])), 2) + + pow(std::abs(std::imag(state[0])), 2) == + Approx(1.0).margin(1e-5); + bool is_one = pow(std::abs(std::real(state[1])), 2) + + pow(std::abs(std::imag(state[1])), 2) == + Approx(1.0).margin(1e-5); + bool is_valid = is_zero ^ is_one; + CHECK(is_valid); + // LCOV_EXCL_STOP +} + +TEST_CASE("Mid-circuit measurement naive test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + intptr_t q; + + q = sim->AllocateQubit(); + + sim->NamedOperation("PauliX", {}, {q}, false); + + auto *m = sim->Measure(q); + + CHECK(*m); +} + +TEST_CASE("Mid-circuit measurement test with postselect = 0", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + intptr_t q; + + q = sim->AllocateQubit(); + + sim->NamedOperation("Hadamard", {}, {q}, false); + + auto *m = sim->Measure(q, 0); + + CHECK(*m == 0); +} + +TEST_CASE("Mid-circuit measurement test with postselect = 1", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + intptr_t q; + + q = sim->AllocateQubit(); + + sim->NamedOperation("Hadamard", {}, {q}, false); + + auto *m = sim->Measure(q, 1); + + CHECK(*m == 1); +} + +TEST_CASE("Mid-circuit measurement test with invalid postselect value", + "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + intptr_t q; + + q = sim->AllocateQubit(); + + sim->NamedOperation("Hadamard", {}, {q}, false); + + REQUIRE_THROWS_WITH(sim->Measure(q, 2), + Catch::Contains("Invalid postselect value")); +} + +TEST_CASE("Expval(ObsT) test with invalid key for cached observables", + "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + REQUIRE_THROWS_WITH(sim->Expval(0), + Catch::Contains("Invalid key for cached observables")); +} + +TEST_CASE("Expval(NamedObs) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[1]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + + CHECK(sim->Expval(px) == Approx(1.0).margin(1e-5)); + CHECK(sim->Expval(py) == Approx(.0).margin(1e-5)); + CHECK(sim->Expval(pz) == Approx(-1.0).margin(1e-5)); +} + +TEST_CASE("Expval(NamedObs) shots test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[1]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + + constexpr std::size_t num_shots = 10000; + sim->SetDeviceShots(num_shots); + + CHECK(sim->Expval(px) == Approx(0.0).margin(5e-2)); + CHECK(sim->Expval(py) == Approx(0.0).margin(5e-2)); + CHECK(sim->Expval(pz) == Approx(-1.0).margin(5e-2)); +} + +TEST_CASE("Expval(HermitianObs) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 2; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + + std::vector> mat1(16, {0, 0}); + std::vector> mat2{ + {1.0, 0.0}, {0.0, 0.0}, {-1.0, 0.0}, {0.0, 0.0}}; + + ObsIdType h1 = sim->Observable(ObsId::Hermitian, mat1, {Qs[0], Qs[1]}); + ObsIdType h2 = sim->Observable(ObsId::Hermitian, mat2, {Qs[0]}); + + CHECK(sim->Expval(h1) == Approx(.0).margin(1e-5)); + CHECK(sim->Expval(h2) == Approx(.0).margin(1e-5)); +} + +TEST_CASE("Expval(HermitianObs) shots test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 2; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + constexpr std::size_t num_shots = 10000; + sim->SetDeviceShots(num_shots); + + sim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + + std::vector> mat1(16, {0, 0}); + + ObsIdType h1 = sim->Observable(ObsId::Hermitian, mat1, {Qs[0], Qs[1]}); + +#ifndef PL_USE_LAPACK + REQUIRE_THROWS_WITH( + sim->Expval(h1), + Catch::Contains( + "Hermitian observables with shot measurement are not supported")); +#else + CHECK(sim->Expval(h1) == Approx(0.0).margin(1e-5)); +#endif +} + +TEST_CASE("Var(HermitianObs) shots test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 2; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + constexpr std::size_t num_shots = 10000; + sim->SetDeviceShots(num_shots); + + sim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + + std::vector> mat1(16, {0, 0}); + + ObsIdType h1 = sim->Observable(ObsId::Hermitian, mat1, {Qs[0], Qs[1]}); +#ifndef PL_USE_LAPACK + REQUIRE_THROWS_WITH( + sim->Var(h1), + Catch::Contains( + "Hermitian observables with shot measurement are not supported")); +#else + CHECK(sim->Var(h1) == Approx(0.0).margin(1e-5)); +#endif +} + +TEST_CASE("Expval(TensorProd(NamedObs)) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[1]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + ObsIdType tpx = sim->TensorObservable({px}); + ObsIdType tpy = sim->TensorObservable({py}); + ObsIdType tpz = sim->TensorObservable({pz}); + + CHECK(sim->Expval(tpx) == Approx(1.0).margin(1e-5)); + CHECK(sim->Expval(tpy) == Approx(.0).margin(1e-5)); + CHECK(sim->Expval(tpz) == Approx(-1.0).margin(1e-5)); +} + +TEST_CASE("Expval(TensorProd(NamedObs)) shots test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[1]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + ObsIdType tpx = sim->TensorObservable({px}); + ObsIdType tpy = sim->TensorObservable({py}); + ObsIdType tpz = sim->TensorObservable({pz}); + + constexpr std::size_t num_shots = 10000; + sim->SetDeviceShots(num_shots); + + CHECK(sim->Expval(tpx) == Approx(1.0).margin(5e-2)); + CHECK(sim->Expval(tpy) == Approx(.0).margin(5e-2)); + CHECK(sim->Expval(tpz) == Approx(-1.0).margin(5e-2)); +} + +TEST_CASE("Expval(TensorProd(NamedObs[])) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[1]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + ObsIdType tpxy = sim->TensorObservable({px, py}); + ObsIdType tpxz = sim->TensorObservable({px, pz}); + + REQUIRE_THROWS_WITH( + sim->TensorObservable({px, py, pz}), + Catch::Contains("All wires in observables must be disjoint.")); + + CHECK(sim->Expval(tpxy) == Approx(0.0).margin(1e-5)); + CHECK(sim->Expval(tpxz) == Approx(-1.0).margin(1e-5)); +} + +TEST_CASE("Expval(TensorProd(NamedObs[])) shots test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + std::unique_ptr sim0 = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + Qs.push_back(sim0->AllocateQubit()); + } + + constexpr std::size_t num_shots = 10000; + sim->SetDeviceShots(num_shots); + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + sim0->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim0->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim0->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[1]}); + ObsIdType tpxy = sim->TensorObservable({px, py}); + + ObsIdType px0 = sim0->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py0 = sim0->Observable(ObsId::PauliY, {}, {Qs[1]}); + ObsIdType tpxy0 = sim0->TensorObservable({px0, py0}); + + CHECK(sim->Expval(tpxy) == Approx(sim0->Expval(tpxy0)).margin(5e-2)); +} + +TEST_CASE("Expval(TensorProd(HermitianObs))", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 2; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + + std::vector> mat1(16, {0, 0}); + std::vector> mat2{ + {1.0, 0.0}, {0.0, 0.0}, {-1.0, 0.0}, {0.0, 0.0}}; + + ObsIdType h1 = sim->Observable(ObsId::Hermitian, mat1, {Qs[0], Qs[1]}); + ObsIdType h2 = sim->Observable(ObsId::Hermitian, mat2, {Qs[0]}); + ObsIdType tph1 = sim->TensorObservable({h1}); + ObsIdType tph2 = sim->TensorObservable({h2}); + + CHECK(sim->Expval(tph1) == Approx(.0).margin(1e-5)); + CHECK(sim->Expval(tph2) == Approx(.0).margin(1e-5)); +} + +TEST_CASE("Expval(TensorProd(HermitianObs[]))", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 2; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + + std::vector> mat1(4, {1.0, 0}); + std::vector> mat2{ + {1.0, 0.0}, {0.0, 0.0}, {-1.0, 0.0}, {0.0, 0.0}}; + + ObsIdType h1 = sim->Observable(ObsId::Hermitian, mat1, {Qs[1]}); + ObsIdType h2 = sim->Observable(ObsId::Hermitian, mat2, {Qs[0]}); + ObsIdType tp = sim->TensorObservable({h1, h2}); + + CHECK(sim->Expval(tp) == Approx(.0).margin(1e-5)); +} + +TEST_CASE("Expval(TensorProd(Obs[]))", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + + std::vector> mat2{ + {1.0, 0.0}, {2.0, 0.0}, {-1.0, 0.0}, {3.0, 0.0}}; + + ObsIdType h = sim->Observable(ObsId::Hermitian, mat2, {Qs[0]}); + ObsIdType tp = sim->TensorObservable({px, h, pz}); + + CHECK(sim->Expval(tp) == Approx(-3.0).margin(1e-5)); +} + +TEST_CASE("Expval(Tensor(Hamiltonian(NamedObs[]), NamedObs)) test", + "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[1]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[0]}); + ObsIdType hxy = sim->HamiltonianObservable({0.4, 0.8}, {px, py}); + ObsIdType thz = sim->TensorObservable({hxy, pz}); + + CHECK(sim->Expval(thz) == Approx(-0.4).margin(1e-5)); +} + +TEST_CASE("Expval(Tensor(HermitianObs, Hamiltonian()) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 3; + std::vector Qs = sim->AllocateQubits(n); + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + + std::vector> mat2{ + {1.0, 0.0}, {0.0, 0.0}, {-1.0, 0.0}, {0.0, 0.0}}; + + ObsIdType her = sim->Observable(ObsId::Hermitian, mat2, {Qs[0]}); + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[1]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[2]}); + ObsIdType hxy = sim->HamiltonianObservable({0.4, 0.8}, {px, py}); + ObsIdType ten = sim->TensorObservable({her, hxy}); + + CHECK(sim->Expval(ten) == Approx(0.0).margin(1e-5)); +} + +TEST_CASE("Expval(Hamiltonian(NamedObs[])) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[1]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + ObsIdType hxyz = sim->HamiltonianObservable({0.4, 0.8, 0.2}, {px, py, pz}); + + CHECK(sim->Expval(hxyz) == Approx(0.2).margin(1e-5)); +} + +TEST_CASE("Expval(Hamiltonian(NamedObs[])) shots test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[1]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + ObsIdType hxyz = sim->HamiltonianObservable({0.4, 0.8, 0.2}, {px, py, pz}); + + constexpr std::size_t num_shots = 10000; + sim->SetDeviceShots(num_shots); + + CHECK(sim->Expval(hxyz) == Approx(0.2).margin(5e-2)); +} + +TEST_CASE("Expval(Hamiltonian(TensorObs[])) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[1]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + ObsIdType tpxy = sim->TensorObservable({px, py}); + ObsIdType tpxz = sim->TensorObservable({px, pz}); + ObsIdType hxyz = sim->HamiltonianObservable({0.2, 0.6}, {tpxy, tpxz}); + + CHECK(sim->Expval(hxyz) == Approx(-.6).margin(1e-5)); +} + +TEST_CASE("Expval(Hamiltonian(Hermitian[])) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + + std::vector> mat2{ + {1.0, 0.0}, {2.0, 0.0}, {-1.0, 0.0}, {3.0, 0.0}}; + ObsIdType h = sim->Observable(ObsId::Hermitian, mat2, {Qs[0]}); + ObsIdType hxhz = sim->HamiltonianObservable({0.2, 0.3, 0.6}, {px, h, pz}); + + CHECK(sim->Expval(hxhz) == Approx(0.5).margin(1e-5)); +} + +TEST_CASE("Expval(Hamiltonian({TensorProd, Hermitian}[])) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + ObsIdType tp = sim->TensorObservable({px, pz}); + + std::vector> mat2{ + {1.0, 0.0}, {2.0, 0.0}, {-1.0, 0.0}, {3.0, 0.0}}; + ObsIdType h = sim->Observable(ObsId::Hermitian, mat2, {Qs[0]}); + ObsIdType hhtp = sim->HamiltonianObservable({0.5, 0.3}, {h, tp}); + + CHECK(sim->Expval(hhtp) == Approx(1.2).margin(1e-5)); +} + +TEST_CASE("Expval(Hamiltonian({Hamiltonian, Hermitian}[])) test", + "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + ObsIdType hp = sim->HamiltonianObservable({0.2, 0.6}, {px, pz}); + + std::vector> mat2{ + {1.0, 0.0}, {2.0, 0.0}, {-1.0, 0.0}, {3.0, 0.0}}; + ObsIdType h = sim->Observable(ObsId::Hermitian, mat2, {Qs[0]}); + ObsIdType hhtp = sim->HamiltonianObservable({0.5, 0.3}, {h, hp}); + + CHECK(sim->Expval(hhtp) == Approx(1.38).margin(1e-5)); +} + +TEST_CASE("Expval(Hamiltonian({Hamiltonian(Hamiltonian), Hermitian}[])) test", + "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + ObsIdType hp = sim->HamiltonianObservable({0.2, 0.6}, {px, pz}); + ObsIdType hhp = sim->HamiltonianObservable({1}, {hp}); + + std::vector> mat2{ + {1.0, 0.0}, {2.0, 0.0}, {-1.0, 0.0}, {3.0, 0.0}}; + ObsIdType h = sim->Observable(ObsId::Hermitian, mat2, {Qs[0]}); + ObsIdType hhtp = sim->HamiltonianObservable({0.5, 0.3}, {hhp, h}); + + CHECK(sim->Expval(hhtp) == Approx(0.7).margin(1e-5)); +} + +TEST_CASE("Var(NamedObs) test with numWires=4", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[0]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[3]}); + + CHECK(sim->Var(px) == Approx(.0).margin(1e-5)); + CHECK(sim->Var(py) == Approx(1.0).margin(1e-5)); + CHECK(sim->Var(pz) == Approx(.0).margin(1e-5)); +} + +TEST_CASE("Var(NamedObs) shots test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 2; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + constexpr std::size_t num_shots = 5000; + sim->SetDeviceShots(num_shots); + + sim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[0]}); + + CHECK(sim->Var(py) == Approx(1.0).margin(5e-2)); +} + +TEST_CASE("Var(HermitianObs) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 2; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + + std::vector> mat1(16, {0, 0}); + std::vector> mat2{ + {1.0, 0.0}, {0.0, 0.0}, {-1.0, 0.0}, {0.0, 0.0}}; + + ObsIdType h1 = sim->Observable(ObsId::Hermitian, mat1, {Qs[0], Qs[1]}); + ObsIdType h2 = sim->Observable(ObsId::Hermitian, mat2, {Qs[0]}); + + CHECK(sim->Var(h1) == Approx(.0).margin(1e-5)); + CHECK(sim->Var(h2) == Approx(1.0).margin(1e-5)); +} + +TEST_CASE("Var(TensorProd(NamedObs)) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[1]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + ObsIdType tpx = sim->TensorObservable({px}); + ObsIdType tpy = sim->TensorObservable({py}); + ObsIdType tpz = sim->TensorObservable({pz}); + + CHECK(sim->Var(tpx) == Approx(.0).margin(1e-5)); + CHECK(sim->Var(tpy) == Approx(1.0).margin(1e-5)); + CHECK(sim->Var(tpz) == Approx(.0).margin(1e-5)); +} + +TEST_CASE("Var(TensorProd(NamedObs)) shots test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + constexpr std::size_t num_shots = 10000; + sim->SetDeviceShots(num_shots); + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + ObsIdType tpx = sim->TensorObservable({px}); + ObsIdType tpz = sim->TensorObservable({pz}); + + CHECK(sim->Var(tpx) == Approx(.0).margin(5e-2)); + CHECK(sim->Var(tpz) == Approx(.0).margin(5e-2)); +} + +TEST_CASE("Var(TensorProd(NamedObs[])) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[1]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + ObsIdType tpxy = sim->TensorObservable({px, py}); + ObsIdType tpxz = sim->TensorObservable({px, pz}); + + CHECK(sim->Var(tpxy) == Approx(1.0).margin(1e-5)); + CHECK(sim->Var(tpxz) == Approx(0.0).margin(1e-5)); +} + +TEST_CASE("Var(TensorProd(HermitianObs)) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 2; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + + std::vector> mat1(16, {0, 0}); + std::vector> mat2{ + {1.0, 0.0}, {0.0, 0.0}, {-1.0, 0.0}, {0.0, 0.0}}; + + ObsIdType h1 = sim->Observable(ObsId::Hermitian, mat1, {Qs[0], Qs[1]}); + ObsIdType h2 = sim->Observable(ObsId::Hermitian, mat2, {Qs[0]}); + ObsIdType tph1 = sim->TensorObservable({h1}); + ObsIdType tph2 = sim->TensorObservable({h2}); + + CHECK(sim->Var(tph1) == Approx(.0).margin(1e-5)); + CHECK(sim->Var(tph2) == Approx(1.0).margin(1e-5)); +} + +TEST_CASE("Var(TensorProd(HermitianObs[])) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 2; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + + std::vector> mat1(4, {1.0, 0}); + std::vector> mat2{ + {1.0, 0.0}, {0.0, 0.0}, {-1.0, 0.0}, {0.0, 0.0}}; + + ObsIdType h1 = sim->Observable(ObsId::Hermitian, mat1, {Qs[1]}); + ObsIdType h2 = sim->Observable(ObsId::Hermitian, mat2, {Qs[0]}); + ObsIdType tp = sim->TensorObservable({h1, h2}); + + CHECK(sim->Var(tp) == Approx(2.0).margin(1e-5)); +} + +TEST_CASE("Var(TensorProd(Obs[])) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + + std::vector> mat2{ + {1.0, 0.0}, {2.0, 0.0}, {-1.0, 0.0}, {3.0, 0.0}}; + + ObsIdType h = sim->Observable(ObsId::Hermitian, mat2, {Qs[0]}); + ObsIdType tp = sim->TensorObservable({px, h, pz}); + + CHECK(sim->Var(tp) == Approx(4.0).margin(1e-5)); +} + +TEST_CASE("Var(Tensor(Hamiltonian(NamedObs[]), NamedObs)) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[1]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[0]}); + ObsIdType hxy = sim->HamiltonianObservable({0.4, 0.8}, {px, py}); + ObsIdType thz = sim->TensorObservable({hxy, pz}); + + CHECK(sim->Var(thz) == Approx(0.64).margin(1e-5)); +} + +TEST_CASE("Var(Tensor(NamedObs[])) shots test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + constexpr std::size_t num_shots = 5000; + sim->SetDeviceShots(num_shots); + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[1]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[0]}); + ObsIdType thz = sim->TensorObservable({px, py, pz}); + + CHECK(sim->Var(thz) == Approx(0.99998976).margin(5e-2)); +} + +TEST_CASE("Var(Tensor(NamedObs[])) shots test without gates " + "(influenced from a bug in Lightning)", + "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 3; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + constexpr std::size_t num_shots = 5000; + sim->SetDeviceShots(num_shots); + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[1]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[0]}); + ObsIdType thz = sim->TensorObservable({px, py, pz}); + + CHECK(sim->Var(thz) == Approx(0.99966144).margin(5e-2)); +} + +TEST_CASE("Var(Tensor(HermitianObs, Hamiltonian()) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 3; + std::vector Qs = sim->AllocateQubits(n); + + sim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + + std::vector> mat2{ + {1.0, 0.0}, {0.0, 0.0}, {-1.0, 0.0}, {0.0, 0.0}}; + + ObsIdType her = sim->Observable(ObsId::Hermitian, mat2, {Qs[0]}); + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[1]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[2]}); + ObsIdType hxy = sim->HamiltonianObservable({0.4, 0.8}, {px, py}); + ObsIdType ten = sim->TensorObservable({her, hxy}); + + CHECK(sim->Var(ten) == Approx(0.8).margin(1e-5)); +} + +TEST_CASE("Var(Tensor(HermitianObs, Hamiltonian()) shots test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 3; + std::vector Qs = sim->AllocateQubits(n); + + constexpr std::size_t num_shots = 5000; + sim->SetDeviceShots(num_shots); + + sim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[1]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[2]}); + ObsIdType hxy = sim->HamiltonianObservable({0.4, 0.8}, {px, py}); + + CHECK(sim->Var(hxy) == Approx(0.8).margin(5e-2)); +} + +TEST_CASE("Var(Hamiltonian(NamedObs[])) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[1]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + ObsIdType hxyz = sim->HamiltonianObservable({0.4, 0.8, 0.2}, {px, py, pz}); + + CHECK(sim->Var(hxyz) == Approx(0.64).margin(1e-5)); +} + +TEST_CASE("Var(Hamiltonian(TensorObs[])) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType py = sim->Observable(ObsId::PauliY, {}, {Qs[1]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + ObsIdType tpxy = sim->TensorObservable({px, py}); + ObsIdType tpxz = sim->TensorObservable({px, pz}); + ObsIdType hxyz = sim->HamiltonianObservable({0.2, 0.6}, {tpxy, tpxz}); + + CHECK(sim->Var(hxyz) == Approx(0.04).margin(1e-5)); +} + +TEST_CASE("Var(Hamiltonian(Hermitian[])) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + + std::vector> mat2{ + {1.0, 0.0}, {2.0, 0.0}, {-1.0, 0.0}, {3.0, 0.0}}; + ObsIdType h = sim->Observable(ObsId::Hermitian, mat2, {Qs[0]}); + ObsIdType hxhz = sim->HamiltonianObservable({0.2, 0.3, 0.6}, {px, h, pz}); + + CHECK(sim->Var(hxhz) == Approx(0.36).margin(1e-5)); +} + +TEST_CASE("Var(Hamiltonian({TensorProd, Hermitian}[])) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + ObsIdType tp = sim->TensorObservable({px, pz}); + + std::vector> mat2{ + {1.0, 0.0}, {2.0, 0.0}, {-1.0, 0.0}, {3.0, 0.0}}; + ObsIdType h = sim->Observable(ObsId::Hermitian, mat2, {Qs[0]}); + ObsIdType hhtp = sim->HamiltonianObservable({0.5, 0.3}, {h, tp}); + + CHECK(sim->Var(hhtp) == Approx(1.0).margin(1e-5)); +} + +TEST_CASE("Var(Hamiltonian({Hamiltonian, Hermitian}[])) test", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + ObsIdType hp = sim->HamiltonianObservable({0.2, 0.6}, {px, pz}); + + std::vector> mat2{ + {1.0, 0.0}, {2.0, 0.0}, {-1.0, 0.0}, {3.0, 0.0}}; + ObsIdType h = sim->Observable(ObsId::Hermitian, mat2, {Qs[0]}); + ObsIdType hhtp = sim->HamiltonianObservable({0.5, 0.3}, {h, hp}); + + CHECK(sim->Var(hhtp) == Approx(1.0).margin(1e-5)); +} + +TEST_CASE("Var(Hamiltonian({Hamiltonian(Hamiltonian), Hermitian}[])) test", + "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("PauliX", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + ObsIdType px = sim->Observable(ObsId::PauliX, {}, {Qs[2]}); + ObsIdType pz = sim->Observable(ObsId::PauliZ, {}, {Qs[1]}); + ObsIdType hp = sim->HamiltonianObservable({0.2, 0.6}, {px, pz}); + ObsIdType hhp = sim->HamiltonianObservable({1}, {hp}); + + std::vector> mat2{ + {1.0, 0.0}, {2.0, 0.0}, {-1.0, 0.0}, {3.0, 0.0}}; + ObsIdType h = sim->Observable(ObsId::Hermitian, mat2, {Qs[0]}); + ObsIdType hhtp = sim->HamiltonianObservable({0.5, 0.3}, {hhp, h}); + + CHECK(sim->Var(hhtp) == Approx(0.36).margin(1e-5)); +} + +TEST_CASE("State test with incorrect size", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs = sim->AllocateQubits(n); + + std::vector> state(1U << (n - 1)); + DataView, 1> view(state); + REQUIRE_THROWS_WITH( + sim->State(view), + Catch::Contains("Invalid size for the pre-allocated state vector")); +} + +TEST_CASE("State test with numWires=4", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs = sim->AllocateQubits(n); + + sim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + std::vector> state(1U << sim->GetNumQubits()); + DataView, 1> view(state); + sim->State(view); + + for (std::size_t i = 0; i < 16; i++) { + if (i == 4 || i == 6 || i == 12 || i == 14) { + CHECK(std::real(state[i]) == Approx(0.).margin(1e-5)); + CHECK(std::imag(state[i]) == Approx(0.5).margin(1e-5)); + } else { + CHECK(std::real(state[i]) == Approx(0.).margin(1e-5)); + CHECK(std::imag(state[i]) == Approx(0.).margin(1e-5)); + } + } +} + +TEST_CASE("PartialProbs test with incorrect numWires and numAlloc", + "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + std::vector probs_vec(1); + DataView probs_view(probs_vec); + + REQUIRE_THROWS_WITH( + sim->PartialProbs(probs_view, {Qs[0], Qs[1], Qs[2], Qs[3], Qs[0]}), + Catch::Contains("Invalid number of wires")); + + REQUIRE_THROWS_WITH( + sim->PartialProbs(probs_view, {Qs[0]}), + Catch::Contains( + "Invalid size for the pre-allocated partial-probabilities")); + + REQUIRE_THROWS_WITH( + sim->Probs(probs_view), + Catch::Contains("Invalid size for the pre-allocated probabilities")); + + sim->ReleaseQubit(Qs[0]); + + REQUIRE_THROWS_WITH(sim->PartialProbs(probs_view, {Qs[0]}), + Catch::Contains("Invalid given wires to measure")); +} + +TEST_CASE("Probs and PartialProbs tests with numWires=0-4", "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + std::vector probs0(1); + DataView view0(probs0); + sim->PartialProbs(view0, std::vector{}); + + std::vector probs1(2); + DataView view1(probs1); + sim->PartialProbs(view1, std::vector{Qs[2]}); + + std::vector probs2(4); + DataView view2(probs2); + sim->PartialProbs(view2, std::vector{Qs[0], Qs[3]}); + + std::vector probs3(16); + DataView view3(probs3); + sim->PartialProbs(view3, Qs); + + std::vector probs4(16); + DataView view4(probs4); + sim->Probs(view4); + + CHECK(probs0.size() == 1); + CHECK(probs0[0] == Approx(1.0)); + CHECK(probs1[0] == Approx(0.5).margin(1e-5)); + CHECK(probs1[1] == Approx(0.5).margin(1e-5)); + for (std::size_t i = 0; i < 4; i++) { + if (i == 0 || i == 2) { + CHECK(probs2[i] == Approx(0.5).margin(1e-5)); + } else { + CHECK(probs2[i] == Approx(0.).margin(1e-5)); + } + } + for (std::size_t i = 0; i < 16; i++) { + if (i == 4 || i == 6 || i == 12 || i == 14) { + CHECK(probs3[i] == Approx(0.25).margin(1e-5)); + CHECK(probs4[i] == Approx(0.25).margin(1e-5)); + } else { + CHECK(probs3[i] == Approx(0.).margin(1e-5)); + CHECK(probs4[i] == Approx(0.).margin(1e-5)); + } + } +} + +TEST_CASE("Probs and PartialProbs shots tests with numWires=0-4", + "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + constexpr std::size_t num_shots = 10000; + sim->SetDeviceShots(num_shots); + + sim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + sim->NamedOperation("PauliY", {}, {Qs[1]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[2]}, false); + sim->NamedOperation("PauliZ", {}, {Qs[3]}, false); + + std::vector probs0(1); + DataView view0(probs0); + sim->PartialProbs(view0, std::vector{}); + + std::vector probs1(2); + DataView view1(probs1); + sim->PartialProbs(view1, std::vector{Qs[2]}); + + std::vector probs2(4); + DataView view2(probs2); + sim->PartialProbs(view2, std::vector{Qs[0], Qs[3]}); + + std::vector probs3(16); + DataView view3(probs3); + sim->PartialProbs(view3, Qs); + + std::vector probs4(16); + DataView view4(probs4); + sim->Probs(view4); + + CHECK(probs0.size() == 1); + CHECK(probs0[0] == Approx(1.0).margin(5e-2)); + CHECK(probs1[0] == Approx(0.5).margin(5e-2)); + CHECK(probs1[1] == Approx(0.5).margin(5e-2)); + for (std::size_t i = 0; i < 4; i++) { + if (i == 0 || i == 2) { + CHECK(probs2[i] == Approx(0.5).margin(5e-2)); + } else { + CHECK(probs2[i] == Approx(0.).margin(5e-2)); + } + } + for (std::size_t i = 0; i < 16; i++) { + if (i == 4 || i == 6 || i == 12 || i == 14) { + CHECK(probs3[i] == Approx(0.25).margin(5e-2)); + CHECK(probs4[i] == Approx(0.25).margin(5e-2)); + } else { + CHECK(probs3[i] == Approx(0.).margin(5e-2)); + CHECK(probs4[i] == Approx(0.).margin(5e-2)); + } + } +} + +TEST_CASE("PartialSample test with incorrect numWires and numAlloc", + "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + std::vector samples_vec(1); + MemRefT samples{samples_vec.data(), + samples_vec.data(), + 0, + {samples_vec.size(), 1}, + {1, 1}}; + DataView view(samples.data_aligned, samples.offset, + samples.sizes, samples.strides); + + REQUIRE_THROWS_WITH( + sim->PartialSample(view, {Qs[0], Qs[1], Qs[2], Qs[3], Qs[0]}, 4), + Catch::Contains("Invalid number of wires")); + + REQUIRE_THROWS_WITH( + sim->PartialSample(view, {Qs[0], Qs[1]}, 2), + Catch::Contains("Invalid size for the pre-allocated partial-samples")); + + REQUIRE_THROWS_WITH( + sim->Sample(view, 2), + Catch::Contains("Invalid size for the pre-allocated samples")); + + sim->ReleaseQubit(Qs[0]); + + REQUIRE_THROWS_WITH(sim->PartialSample(view, {Qs[0]}, 4), + Catch::Contains("Invalid given wires to measure")); +} + +TEST_CASE("PartialCounts test with incorrect numWires and numAlloc", + "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + std::vector eigvals_vec(1); + DataView eigvals_view(eigvals_vec); + + std::vector counts_vec(1); + DataView counts_view(counts_vec); + + REQUIRE_THROWS_WITH(sim->PartialCounts(eigvals_view, counts_view, + {Qs[0], Qs[1], Qs[2], Qs[3], Qs[0]}, + 4), + Catch::Contains("Invalid number of wires")); + + REQUIRE_THROWS_WITH( + sim->PartialCounts(eigvals_view, counts_view, {Qs[0]}, 1), + Catch::Contains("Invalid size for the pre-allocated partial-counts")); + + REQUIRE_THROWS_WITH( + sim->Counts(eigvals_view, counts_view, 1), + Catch::Contains("Invalid size for the pre-allocated counts")); + + sim->ReleaseQubit(Qs[0]); + + REQUIRE_THROWS_WITH( + sim->PartialCounts(eigvals_view, counts_view, {Qs[0]}, 4), + Catch::Contains("Invalid given wires to measure")); +} + +TEST_CASE("Sample and PartialSample tests with numWires=0-4 shots=100", + "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("RX", {0.5}, {Qs[0]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[1]}, false); + sim->NamedOperation("CNOT", {}, {Qs[0], Qs[1]}, false); + + std::size_t shots = 100; + + std::vector samples1(shots * 1); + MemRefT buffer1{ + samples1.data(), samples1.data(), 0, {shots, 1}, {1, 1}}; + DataView view1(buffer1.data_aligned, buffer1.offset, + buffer1.sizes, buffer1.strides); + sim->PartialSample(view1, std::vector{Qs[2]}, shots); + + std::vector samples2(shots * 2); + MemRefT buffer2{ + samples2.data(), samples2.data(), 0, {shots, 2}, {1, 1}}; + DataView view2(buffer2.data_aligned, buffer2.offset, + buffer2.sizes, buffer2.strides); + sim->PartialSample(view2, std::vector{Qs[0], Qs[3]}, shots); + + std::vector samples3(shots * 4); + MemRefT buffer3{ + samples3.data(), samples3.data(), 0, {shots, 4}, {1, 1}}; + DataView view3(buffer3.data_aligned, buffer3.offset, + buffer3.sizes, buffer3.strides); + sim->PartialSample(view3, Qs, shots); + + std::vector samples4(shots * 4); + MemRefT buffer4{ + samples4.data(), samples4.data(), 0, {shots, 4}, {1, 1}}; + DataView view4(buffer4.data_aligned, buffer4.offset, + buffer4.sizes, buffer4.strides); + sim->Sample(view4, shots); + + for (std::size_t i = 0; i < shots * 1; i++) { + CHECK((samples1[i] == 0. || samples1[i] == 1.)); + } + for (std::size_t i = 0; i < shots * 2; i++) { + CHECK((samples2[i] == 0. || samples2[i] == 1.)); + } + for (std::size_t i = 0; i < shots * 4; i++) { + CHECK((samples3[i] == 0. || samples3[i] == 1.)); + } + for (std::size_t i = 0; i < shots * 4; i++) { + CHECK((samples4[i] == 0. || samples4[i] == 1.)); + } +} + +TEST_CASE("Sample and PartialSample tests with numWires=0-4 " + "shots=1000 mcmc=True num_burnin=200", + "[Measures]") { + std::unique_ptr sim = + std::make_unique("{mcmc : True, num_burnin : 200}"); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs; + Qs.reserve(n); + for (std::size_t i = 0; i < n; i++) { + Qs.push_back(sim->AllocateQubit()); + } + + sim->NamedOperation("RX", {0.5}, {Qs[0]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[1]}, false); + sim->NamedOperation("CNOT", {}, {Qs[0], Qs[1]}, false); + + std::size_t shots = 100; + + std::vector samples1(shots * 1); + MemRefT buffer1{ + samples1.data(), samples1.data(), 0, {shots, 1}, {1, 1}}; + DataView view1(buffer1.data_aligned, buffer1.offset, + buffer1.sizes, buffer1.strides); + sim->PartialSample(view1, std::vector{Qs[2]}, shots); + + std::vector samples2(shots * 2); + MemRefT buffer2{ + samples2.data(), samples2.data(), 0, {shots, 2}, {1, 1}}; + DataView view2(buffer2.data_aligned, buffer2.offset, + buffer2.sizes, buffer2.strides); + sim->PartialSample(view2, std::vector{Qs[0], Qs[3]}, shots); + + std::vector samples3(shots * 4); + MemRefT buffer3{ + samples3.data(), samples3.data(), 0, {shots, 4}, {1, 1}}; + DataView view3(buffer3.data_aligned, buffer3.offset, + buffer3.sizes, buffer3.strides); + sim->PartialSample(view3, Qs, shots); + + std::vector samples4(shots * 4); + MemRefT buffer4{ + samples4.data(), samples4.data(), 0, {shots, 4}, {1, 1}}; + DataView view4(buffer4.data_aligned, buffer4.offset, + buffer4.sizes, buffer4.strides); + sim->Sample(view4, shots); + + for (std::size_t i = 0; i < shots * 1; i++) { + CHECK((samples1[i] == 0. || samples1[i] == 1.)); + } + for (std::size_t i = 0; i < shots * 2; i++) { + CHECK((samples2[i] == 0. || samples2[i] == 1.)); + } + for (std::size_t i = 0; i < shots * 4; i++) { + CHECK((samples3[i] == 0. || samples3[i] == 1.)); + } + for (std::size_t i = 0; i < shots * 4; i++) { + CHECK((samples4[i] == 0. || samples4[i] == 1.)); + } +} + +TEST_CASE("Counts and PartialCounts tests with numWires=0-4 shots=100", + "[Measures]") { + std::unique_ptr sim = std::make_unique(); + + // state-vector with #qubits = n + constexpr std::size_t n = 4; + std::vector Qs = sim->AllocateQubits(n); + + sim->NamedOperation("RX", {0.5}, {Qs[0]}, false); + sim->NamedOperation("Hadamard", {}, {Qs[1]}, false); + sim->NamedOperation("CNOT", {}, {Qs[0], Qs[1]}, false); + + std::size_t shots = 100; + + std::vector eigvals0(1); + std::vector counts0(1); + DataView eview0(eigvals0); + DataView cview0(counts0); + sim->PartialCounts(eview0, cview0, std::vector{}, shots); + + std::vector eigvals1(2); + std::vector counts1(2); + DataView eview1(eigvals1); + DataView cview1(counts1); + sim->PartialCounts(eview1, cview1, std::vector{Qs[2]}, shots); + + std::vector eigvals2(4); + std::vector counts2(4); + DataView eview2(eigvals2); + DataView cview2(counts2); + sim->PartialCounts(eview2, cview2, std::vector{Qs[0], Qs[3]}, + shots); + + std::vector eigvals3(16); + std::vector counts3(16); + DataView eview3(eigvals3); + DataView cview3(counts3); + sim->PartialCounts(eview3, cview3, Qs, shots); + + std::vector eigvals4(16); + std::vector counts4(16); + DataView eview4(eigvals4); + DataView cview4(counts4); + sim->Counts(eview4, cview4, shots); + + CHECK(eigvals0.size() == 1); + CHECK(eigvals0[0] == 0.0); + CHECK(counts0.size() == 1); + CHECK(counts0[0] == static_cast(shots)); + CHECK((eigvals1[0] == 0. && eigvals1[1] == 1.)); + CHECK((eigvals2[0] == 0. && eigvals2[1] == 1. && eigvals2[2] == 2. && + eigvals2[3] == 3.)); + for (std::size_t i = 0; i < 16; i++) { + CHECK(eigvals3[i] == static_cast(i)); + CHECK(eigvals4[i] == static_cast(i)); + } + + CHECK(counts1[0] + counts1[1] == static_cast(shots)); + CHECK(counts2[0] + counts2[1] + counts2[2] + counts2[3] == + static_cast(shots)); + std::size_t sum3 = 0; + std::size_t sum4 = 0; + for (std::size_t i = 0; i < 16; i++) { + sum3 += counts3[i]; + sum4 += counts4[i]; + } + CHECK(sum3 == shots); + CHECK(sum4 == shots); +} + +TEST_CASE("Measurement with a seeded device", "[Measures]") { + std::array, 2> sims; + std::vector gens{std::mt19937{37}, std::mt19937{37}}; + + auto circuit = [](LQSimulator &sim, std::mt19937 &gen) { + sim.SetDevicePRNG(&gen); + std::vector Qs; + Qs.reserve(1); + Qs.push_back(sim.AllocateQubit()); + sim.NamedOperation("Hadamard", {}, {Qs[0]}, false); + auto *m = sim.Measure(Qs[0]); + return m; + }; + + for (std::size_t trial = 0; trial < 5; trial++) { + sims[0] = std::make_unique(); + sims[1] = std::make_unique(); + + auto *m0 = circuit(*(sims[0]), gens[0]); + auto *m1 = circuit(*(sims[1]), gens[1]); + + CHECK(*m0 == *m1); + } +} + +TEST_CASE("Sample with a seeded device", "[Measures]") { + std::size_t shots = 100; + std::array, 2> sims; + std::vector> sample_vec(2, + std::vector(shots * 4)); + + std::vector> buffers{ + MemRefT{ + sample_vec[0].data(), sample_vec[0].data(), 0, {shots, 1}, {1, 1}}, + MemRefT{ + sample_vec[1].data(), sample_vec[1].data(), 0, {shots, 1}, {1, 1}}, + }; + std::vector> views{ + DataView(buffers[0].data_aligned, buffers[0].offset, + buffers[0].sizes, buffers[0].strides), + DataView(buffers[1].data_aligned, buffers[1].offset, + buffers[1].sizes, buffers[1].strides)}; + + std::vector gens{std::mt19937{37}, std::mt19937{37}}; + + auto circuit = [shots](LQSimulator &sim, DataView &view, + std::mt19937 &gen) { + sim.SetDevicePRNG(&gen); + std::vector Qs; + Qs.reserve(1); + Qs.push_back(sim.AllocateQubit()); + sim.NamedOperation("Hadamard", {}, {Qs[0]}, false); + sim.NamedOperation("RX", {0.5}, {Qs[0]}, false); + sim.Sample(view, shots); + }; + + for (std::size_t trial = 0; trial < 5; trial++) { + sims[0] = std::make_unique(); + sims[1] = std::make_unique(); + + for (std::size_t sim_idx = 0; sim_idx < sims.size(); sim_idx++) { + circuit(*(sims[sim_idx]), views[sim_idx], gens[sim_idx]); + } + + for (std::size_t i = 0; i < sample_vec[0].size(); i++) { + CHECK((sample_vec[0][i] == sample_vec[1][i])); + } + } +} diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/Test_LightningSimulator.cpp b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/Test_LightningSimulator.cpp new file mode 100644 index 000000000..ad89ea5ad --- /dev/null +++ b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/Test_LightningSimulator.cpp @@ -0,0 +1,703 @@ +// Copyright 2018 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. +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include + +#include "LightningSimulator.hpp" +#include "QuantumDevice.hpp" +#include "TestHelpers.hpp" + +/// @cond DEV +namespace { +using namespace Catalyst::Runtime::Simulator; +using namespace Pennylane::Util; +using LQSimulator = LightningSimulator; +using QDevice = Catalyst::Runtime::QuantumDevice; + +GENERATE_DEVICE_FACTORY(LightningSimulator, + Catalyst::Runtime::Simulator::LightningSimulator); +} // namespace +/// @endcond + +/** + * @brief Tests the LightningSimulator class. + * + */ +TEST_CASE("LightningSimulator::constructor", "[constructibility]") { + SECTION("LightningSimulator") { + REQUIRE(std::is_constructible::value); + } + SECTION("LightningSimulator(string))") { + REQUIRE(std::is_constructible::value); + } +} + +TEST_CASE("Test the device factory method", "[constructibility]") { + std::unique_ptr LQsim(LightningSimulatorFactory("")); + REQUIRE(LQsim->GetNumQubits() == 0); +} + +TEST_CASE("LightningSimulator::unit_tests", "[unit tests]") { + SECTION("Managing Qubits") { + std::unique_ptr LQsim = std::make_unique(); + std::vector Qs = LQsim->AllocateQubits(0); + REQUIRE(LQsim->GetNumQubits() == 0); + LQsim->AllocateQubits(2); + REQUIRE(LQsim->GetNumQubits() == 2); + LQsim->AllocateQubits(2); + REQUIRE(LQsim->GetNumQubits() == 4); + LQsim->ReleaseQubit(0); + REQUIRE(LQsim->GetNumQubits() == + 3); // releasing one qubit changes the total number. + LQsim->ReleaseAllQubits(); + REQUIRE(LQsim->GetNumQubits() == + 0); // releasing all qubits resets the simulator. + } + SECTION("Tape recording") { + std::unique_ptr LQsim = std::make_unique(); + std::vector Qs = LQsim->AllocateQubits(1); + REQUIRE_NOTHROW(LQsim->StartTapeRecording()); + REQUIRE_THROWS_WITH( + LQsim->StartTapeRecording(), + Catch::Matchers::Contains("Cannot re-activate the cache manager")); + REQUIRE_NOTHROW(LQsim->StopTapeRecording()); + REQUIRE_THROWS_WITH( + LQsim->StopTapeRecording(), + Catch::Matchers::Contains( + "Cannot stop an already stopped cache manager")); + } +} + +TEST_CASE("LightningSimulator::GateSet", "[GateSet]") { + SECTION("Identity gate") { + std::unique_ptr LQsim = std::make_unique(); + + constexpr std::size_t n_qubits = 10; + std::vector Qs; + Qs.reserve(n_qubits); + for (std::size_t ind = 0; ind < n_qubits; ind++) { + Qs[ind] = LQsim->AllocateQubit(); + } + + for (std::size_t ind = 0; ind < n_qubits; ind += 2) { + LQsim->NamedOperation("Identity", {}, {Qs[ind]}, false); + } + + std::vector> state(1U << LQsim->GetNumQubits()); + DataView, 1> view(state); + LQsim->State(view); + + CHECK(state.at(0) == std::complex{1, 0}); + + std::complex sum{0, 0}; + for (std::size_t ind = 1; ind < state.size(); ind++) { + sum += state[ind]; + } + + CHECK(sum == std::complex{0, 0}); + } + + SECTION("PauliX gate") { + std::unique_ptr LQsim = std::make_unique(); + + constexpr std::size_t n_qubits = 3; + std::vector Qs; + Qs.reserve(n_qubits); + for (std::size_t ind = 0; ind < n_qubits; ind++) { + Qs[ind] = LQsim->AllocateQubit(); + } + + for (std::size_t ind = 0; ind < n_qubits; ind++) { + LQsim->NamedOperation("PauliX", {}, {Qs[ind]}, false); + } + for (std::size_t ind = n_qubits; ind > 0; ind--) { + LQsim->NamedOperation("PauliX", {}, {Qs[ind - 1]}, false); + } + + std::vector> state(1U << LQsim->GetNumQubits()); + DataView, 1> view(state); + LQsim->State(view); + + CHECK(state.at(0) == std::complex{1, 0}); + + std::complex sum{0, 0}; + for (std::size_t ind = 1; ind < state.size(); ind++) { + sum += state[ind]; + } + + CHECK(sum == std::complex{0, 0}); + } + + SECTION("PauliY gate") { + std::unique_ptr LQsim = std::make_unique(); + + constexpr std::size_t n_qubits = 2; + std::vector Qs; + Qs.reserve(n_qubits); + for (std::size_t ind = 0; ind < n_qubits; ind++) { + Qs[ind] = LQsim->AllocateQubit(); + } + + for (std::size_t ind = 0; ind < n_qubits; ind++) { + LQsim->NamedOperation("PauliY", {}, {Qs[ind]}, false); + } + + std::vector> state(1U << LQsim->GetNumQubits()); + DataView, 1> view(state); + LQsim->State(view); + + CHECK(state.at(0) == std::complex{0, 0}); + CHECK(state.at(1) == std::complex{0, 0}); + CHECK(state.at(2) == std::complex{0, 0}); + CHECK(state.at(3) == std::complex{-1, 0}); + } + + SECTION("PauliY and PauliZ gates") { + std::unique_ptr LQsim = std::make_unique(); + + constexpr std::size_t n_qubits = 2; + std::vector Qs; + Qs.reserve(n_qubits); + for (std::size_t ind = 0; ind < n_qubits; ind++) { + Qs[ind] = LQsim->AllocateQubit(); + } + + LQsim->NamedOperation("PauliY", {}, {Qs[0]}, false); + LQsim->NamedOperation("PauliZ", {}, {Qs[1]}, false); + + std::vector> state(1U << LQsim->GetNumQubits()); + DataView, 1> view(state); + LQsim->State(view); + + CHECK(state.at(0) == std::complex{0, 0}); + CHECK(state.at(1) == std::complex{0, 0}); + CHECK(state.at(2) == std::complex{0, 1}); + CHECK(state.at(3) == std::complex{0, 0}); + } + + SECTION("Hadamard gate") { + std::unique_ptr LQsim = std::make_unique(); + + constexpr std::size_t n_qubits = 2; + std::vector Qs; + Qs.reserve(n_qubits); + for (std::size_t ind = 0; ind < n_qubits; ind++) { + Qs[ind] = LQsim->AllocateQubit(); + } + + for (std::size_t ind = 0; ind < n_qubits; ind++) { + LQsim->NamedOperation("Hadamard", {}, {Qs[ind]}, false); + } + + std::vector> state(1U << LQsim->GetNumQubits()); + DataView, 1> view(state); + LQsim->State(view); + + CHECK(state[0] == + PLApproxComplex(std::complex{0.5, 0}).epsilon(1e-5)); + CHECK(state.at(1) == state.at(0)); + CHECK(state.at(2) == state.at(0)); + CHECK(state.at(3) == state.at(0)); + } + + SECTION("R(X, Y, Z) and PauliX gates") { + std::unique_ptr LQsim = std::make_unique(); + + constexpr std::size_t n_qubits = 4; + std::vector Qs = LQsim->AllocateQubits(n_qubits); + + LQsim->NamedOperation("PauliX", {}, {Qs[0]}, false); + + LQsim->NamedOperation("RX", {0.123}, {Qs[1]}, false); + LQsim->NamedOperation("RY", {0.456}, {Qs[2]}, false); + LQsim->NamedOperation("RZ", {0.789}, {Qs[3]}, false); + + std::vector> state(1U << LQsim->GetNumQubits()); + DataView, 1> view(state); + LQsim->State(view); + + // calculated by pennylane. + CHECK(state.at(0) == std::complex{0, 0}); + CHECK(state.at(1) == std::complex{0, 0}); + CHECK(state.at(2) == std::complex{0, 0}); + CHECK(state.at(3) == std::complex{0, 0}); + CHECK(state.at(4) == std::complex{0, 0}); + CHECK(state.at(5) == std::complex{0, 0}); + CHECK(state.at(6) == std::complex{0, 0}); + CHECK(state.at(7) == std::complex{0, 0}); + CHECK(state[8] == + PLApproxComplex( + std::complex{0.8975969498074641, -0.3736920921192206}) + .epsilon(1e-5)); + CHECK(state.at(9) == std::complex{0, 0}); + CHECK(state[10] == + PLApproxComplex(std::complex{0.20827363966052723, + -0.08670953277495183}) + .epsilon(1e-5)); + CHECK(state.at(11) == std::complex{0, 0}); + CHECK(state[12] == + PLApproxComplex(std::complex{-0.023011082205037697, + -0.055271914055973925}) + .epsilon(1e-5)); + CHECK(state.at(13) == std::complex{0, 0}); + CHECK(state[14] == + PLApproxComplex(std::complex{-0.005339369573836912, + -0.012825002038956146}) + .epsilon(1e-5)); + CHECK(state.at(15) == std::complex{0, 0}); + } + + SECTION("Hadamard, RX, PhaseShift with cache manager") { + std::unique_ptr LQsim = std::make_unique(); + + constexpr std::size_t n_qubits = 2; + std::vector Qs; + Qs.reserve(n_qubits); + + Qs[0] = LQsim->AllocateQubit(); + Qs[1] = LQsim->AllocateQubit(); + + LQsim->StartTapeRecording(); + LQsim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + LQsim->NamedOperation("RX", {0.123}, {Qs[1]}, false); + LQsim->NamedOperation("PhaseShift", {0.456}, {Qs[0]}, false); + LQsim->StopTapeRecording(); + + std::vector> state(1U << LQsim->GetNumQubits()); + DataView, 1> view(state); + LQsim->State(view); + + // calculated by pennylane. + CHECK(state[0] == PLApproxComplex(std::complex{0.7057699753, 0}) + .epsilon(1e-5)); + CHECK(state[1] == PLApproxComplex(std::complex{0, -0.04345966}) + .epsilon(1e-5)); + CHECK(state[2] == + PLApproxComplex(std::complex{0.63365519, 0.31079312}) + .epsilon(1e-5)); + CHECK(state[3] == + PLApproxComplex(std::complex{0.01913791, -0.039019}) + .epsilon(1e-5)); + + std::tuple, std::vector> + expected{3, 0, 2, {"Hadamard", "RX", "PhaseShift"}, {}}; + REQUIRE(LQsim->CacheManagerInfo() == expected); + } + + // ============= 2-qubit operations ============= + + SECTION("PauliX and CNOT") { + std::unique_ptr LQsim = std::make_unique(); + + constexpr std::size_t n_qubits = 2; + std::vector Qs; + Qs.reserve(n_qubits); + + for (std::size_t i = 0; i < n_qubits; i++) { + Qs[i] = LQsim->AllocateQubit(); + } + + LQsim->NamedOperation("PauliX", {}, {Qs[0]}, false); + LQsim->NamedOperation("CNOT", {}, {Qs[0], Qs[1]}, false); + + std::vector> state(1U << LQsim->GetNumQubits()); + DataView, 1> view(state); + LQsim->State(view); + + CHECK(state.at(0) == std::complex{0, 0}); + CHECK(state.at(1) == std::complex{0, 0}); + CHECK(state.at(2) == std::complex{0, 0}); + CHECK(state.at(3) == std::complex{1, 0}); + } + + SECTION("Hadamard and CR(X, Y, Z)") { + std::unique_ptr LQsim = std::make_unique(); + + constexpr std::size_t n_qubits = 4; + std::vector Qs = LQsim->AllocateQubits(n_qubits); + + LQsim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + LQsim->NamedOperation("CRX", {0.123}, {Qs[0], Qs[1]}, false); + LQsim->NamedOperation("CRY", {0.456}, {Qs[0], Qs[2]}, false); + LQsim->NamedOperation("CRZ", {0.789}, {Qs[0], Qs[3]}, false); + + std::vector> state(1U << LQsim->GetNumQubits()); + DataView, 1> view(state); + LQsim->State(view); + + // calculated by pennylane. + CHECK( + state[0] == + PLApproxComplex(std::complex{M_SQRT1_2, 0}).epsilon(1e-5)); + CHECK(state.at(1) == std::complex{0, 0}); + CHECK(state.at(2) == std::complex{0, 0}); + CHECK(state.at(3) == std::complex{0, 0}); + CHECK(state.at(4) == std::complex{0, 0}); + CHECK(state.at(5) == std::complex{0, 0}); + CHECK(state.at(6) == std::complex{0, 0}); + CHECK(state.at(7) == std::complex{0, 0}); + CHECK(state[8] == + PLApproxComplex( + std::complex{0.6346968899812189, -0.2642402124132889}) + .epsilon(1e-5)); + CHECK(state.at(9) == std::complex{0, 0}); + CHECK(state[10] == + PLApproxComplex(std::complex{0.14727170294636227, + -0.061312898618685635}) + .epsilon(1e-5)); + CHECK(state.at(11) == std::complex{0, 0}); + CHECK(state[12] == + PLApproxComplex(std::complex{-0.016271292269623247, + -0.03908314523813921}) + .epsilon(1e-5)); + CHECK(state.at(13) == std::complex{0, 0}); + CHECK(state[14] == + PLApproxComplex(std::complex{-0.0037755044329212074, + -0.009068645910477189}) + .epsilon(1e-5)); + CHECK(state.at(15) == std::complex{0, 0}); + } + + SECTION("Hadamard and CRot") { + std::unique_ptr LQsim = std::make_unique(); + + constexpr std::size_t n_qubits = 2; + std::vector Qs = LQsim->AllocateQubits(n_qubits); + + LQsim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + LQsim->NamedOperation("CRot", {M_PI, M_PI_2, 0.5}, {Qs[0], Qs[1]}, + false); + + std::vector> state(1U << LQsim->GetNumQubits()); + DataView, 1> view(state); + LQsim->State(view); + + CHECK( + state[0] == + PLApproxComplex(std::complex{M_SQRT1_2, 0}).epsilon(1e-5)); + + CHECK(state[1] == + PLApproxComplex(std::complex{0, 0}).epsilon(1e-5)); + + CHECK(state[2] == PLApproxComplex(std::complex{-0.1237019796, + -0.4844562109}) + .epsilon(1e-5)); + CHECK(state[3] == + PLApproxComplex(std::complex{0.1237019796, -0.4844562109}) + .epsilon(1e-5)); + } + + SECTION("Hadamard, PauliZ, IsingXY, SWAP") { + std::unique_ptr LQsim = std::make_unique(); + + constexpr std::size_t n_qubits = 2; + std::vector Qs = LQsim->AllocateQubits(n_qubits); + + LQsim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + LQsim->NamedOperation("PauliZ", {}, {Qs[0]}, false); + LQsim->NamedOperation("IsingXY", {0.2}, {Qs[1], Qs[0]}, false); + LQsim->NamedOperation("SWAP", {}, {Qs[0], Qs[1]}, false); + + std::vector> state(1U << LQsim->GetNumQubits()); + DataView, 1> view(state); + LQsim->State(view); + + CHECK( + state[0] == + PLApproxComplex(std::complex{M_SQRT1_2, 0}).epsilon(1e-5)); + CHECK(state[1] == PLApproxComplex(std::complex{-0.70357419, 0}) + .epsilon(1e-5)); + CHECK(state[2] == PLApproxComplex(std::complex{0, -0.07059289}) + .epsilon(1e-5)); + CHECK(state[3] == + PLApproxComplex(std::complex{0, 0}).epsilon(1e-5)); + } + + SECTION("Hadamard, PauliX and Toffoli") { + std::unique_ptr LQsim = std::make_unique(); + + constexpr std::size_t n_qubits = 3; + std::vector Qs = LQsim->AllocateQubits(n_qubits); + + LQsim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + LQsim->NamedOperation("PauliX", {}, {Qs[1]}, false); + LQsim->NamedOperation("Toffoli", {}, {Qs[0], Qs[1], Qs[2]}, false); + + std::vector> state(1U << LQsim->GetNumQubits()); + DataView, 1> view(state); + LQsim->State(view); + + CHECK(state.at(0) == std::complex{0, 0}); + CHECK(state.at(1) == std::complex{0, 0}); + CHECK( + state[2] == + PLApproxComplex(std::complex{M_SQRT1_2, 0}).epsilon(1e-5)); + CHECK(state.at(3) == std::complex{0, 0}); + CHECK(state.at(4) == std::complex{0, 0}); + CHECK(state.at(5) == std::complex{0, 0}); + CHECK(state.at(6) == std::complex{0, 0}); + CHECK( + state[7] == + PLApproxComplex(std::complex{M_SQRT1_2, 0}).epsilon(1e-5)); + } + + SECTION("RX, Hadamard and MultiRZ") { + std::unique_ptr LQsim = std::make_unique(); + + constexpr std::size_t n_qubits = 2; + std::vector Qs = LQsim->AllocateQubits(n_qubits); + + LQsim->NamedOperation("RX", {M_PI}, {Qs[1]}, false); + LQsim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + LQsim->NamedOperation("Hadamard", {}, {Qs[1]}, false); + LQsim->NamedOperation("MultiRZ", {M_PI}, {Qs[0], Qs[1]}, false); + LQsim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + LQsim->NamedOperation("Hadamard", {}, {Qs[1]}, false); + + std::vector> state(1U << LQsim->GetNumQubits()); + DataView, 1> view(state); + LQsim->State(view); + + CHECK(state[2] == + PLApproxComplex(std::complex{-1, 0}).epsilon(1e-5)); + } + + SECTION("Hadamard, CNOT and Matrix") { + std::unique_ptr LQsim = std::make_unique(); + + constexpr std::size_t n_qubits = 2; + std::vector Qs = LQsim->AllocateQubits(n_qubits); + + LQsim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + LQsim->NamedOperation("CNOT", {}, {Qs[0], Qs[1]}, false); + + const std::vector wires = {Qs[0]}; + std::vector> matrix{ + {-0.6709485262524046, -0.6304426335363695}, + {-0.14885403153998722, 0.3608498832392019}, + {-0.2376311670004963, 0.3096798175687841}, + {-0.8818365947322423, -0.26456390390903695}, + }; + LQsim->MatrixOperation(matrix, wires, false); + + std::vector> state(1U << LQsim->GetNumQubits()); + DataView, 1> view(state); + LQsim->State(view); + + CHECK(state[0] == + PLApproxComplex(std::complex{-0.474432, -0.44579}) + .epsilon(1e-5)); + CHECK(state[1] == + PLApproxComplex(std::complex{-0.105256, 0.255159}) + .epsilon(1e-5)); + CHECK(state[2] == + PLApproxComplex(std::complex{-0.168031, 0.218977}) + .epsilon(1e-5)); + CHECK(state[3] == + PLApproxComplex(std::complex{-0.623553, -0.187075}) + .epsilon(1e-5)); + } + + SECTION("Hadamard, CR(X, Y, Z) and Matrix") { + std::unique_ptr LQsim = std::make_unique(); + + constexpr std::size_t n_qubits = 4; + std::vector Qs = LQsim->AllocateQubits(n_qubits); + + LQsim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + LQsim->NamedOperation("CRX", {0.123}, {Qs[0], Qs[1]}, false); + LQsim->NamedOperation("CRY", {0.456}, {Qs[0], Qs[2]}, false); + LQsim->NamedOperation("CRZ", {0.789}, {Qs[0], Qs[3]}, false); + + const std::vector wires = {Qs[0], Qs[1], Qs[2]}; + std::vector> matrix{ + {-0.14601911598243822, -0.18655250647340088}, + {-0.03917826201290317, -0.031161687050443518}, + {0.11497626236175404, 0.38310733543366354}, + {-0.0929691815340695, 0.1219804125497268}, + {0.07306514883467692, 0.017445444816725875}, + {-0.27330866098918355, -0.6007032759764033}, + {0.4530754397715841, -0.08267189625512258}, + {0.32125201986075, -0.036845158875036116}, + {0.032317572838307884, 0.02292755555300329}, + {-0.18775945295623664, -0.060215004737844156}, + {-0.3093351335745536, -0.2061961962889725}, + {0.4216087567144761, 0.010534488410902099}, + {0.2769943541718527, -0.26016137877135465}, + {0.18727884147867532, 0.02830415812286322}, + {0.3367562196770689, -0.5250999173939218}, + {0.05770014289220745, 0.26595514845958573}, + {0.37885720163317027, 0.3110931426403546}, + {0.13436510737129648, -0.4083415934958021}, + {-0.5443665467635203, 0.2458343977310266}, + {-0.050346912365833024, 0.08709833123617361}, + {0.11505259829552131, 0.010155858056939438}, + {-0.2930849061531229, 0.019339259194141145}, + {0.011825409829453282, 0.011597907736881019}, + {-0.10565527258356637, -0.3113689446440079}, + {0.0273191284561944, -0.2479498526173881}, + {-0.5528072425836249, -0.06114469689935285}, + {-0.20560364740746587, -0.3800208994544297}, + {-0.008236143958221483, 0.3017421511504845}, + {0.04817188123334976, 0.08550951191632741}, + {-0.24081054643565586, -0.3412671345149831}, + {-0.38913538197001885, 0.09288402897806938}, + {-0.07937578245883717, 0.013979426755633685}, + {0.22246583652015395, -0.18276674810033927}, + {0.22376666162382491, 0.2995723155125488}, + {-0.1727191441070097, -0.03880522034607489}, + {0.075780203819001, 0.2818783673816625}, + {-0.6161322400651016, 0.26067347179217193}, + {-0.021161519614267765, -0.08430919051054794}, + {0.1676500381348944, -0.30645601624407504}, + {-0.28858251997285883, 0.018089595494883842}, + {-0.19590767481842053, -0.12844366632033652}, + {0.18707834504831794, -0.1363932722670649}, + {-0.07224221779769334, -0.11267803536286894}, + {-0.23897684826459387, -0.39609971967853685}, + {-0.0032110880452929555, -0.29294331305690136}, + {-0.3188741682462722, -0.17338979346647143}, + {0.08194395032821632, -0.002944814673179825}, + {-0.5695791830944521, 0.33299548924055095}, + {-0.4983660307441444, -0.4222358493977972}, + {0.05533914327048402, -0.42575842134560576}, + {-0.2187623521182678, -0.03087596187054778}, + {0.11278255885846857, 0.07075886163492914}, + {-0.3054684775292515, -0.1739796870866232}, + {0.14151567663565712, 0.20399935744127418}, + {0.06720165377364941, 0.07543463072363207}, + {0.08019665306716581, -0.3473013434358584}, + {-0.2600167605995786, -0.08795704036197827}, + {0.125680477777759, 0.266342700305046}, + {-0.1586772594600269, 0.187360909108502}, + {-0.4653314704208982, 0.4048609954619629}, + {0.39992560380733094, -0.10029244177901954}, + {0.2533527906886461, 0.05222114898540775}, + {-0.15840033949128557, -0.2727320427534386}, + {-0.21590866323269536, -0.1191163626522938}, + }; + LQsim->MatrixOperation(matrix, wires, false); + + std::vector> state(1U << LQsim->GetNumQubits()); + DataView, 1> view(state); + LQsim->State(view); + + CHECK(state[0] == + PLApproxComplex(std::complex{-0.141499, -0.230993}) + .epsilon(1e-5)); + CHECK(state[2] == + PLApproxComplex(std::complex{0.135423, -0.235563}) + .epsilon(1e-5)); + CHECK(state[4] == + PLApproxComplex(std::complex{0.299458, 0.218321}) + .epsilon(1e-5)); + CHECK(state[6] == + PLApproxComplex(std::complex{0.0264869, -0.154913}) + .epsilon(1e-5)); + CHECK(state[8] == + PLApproxComplex(std::complex{-0.186607, 0.188884}) + .epsilon(1e-5)); + CHECK(state[10] == + PLApproxComplex(std::complex{-0.271843, -0.281136}) + .epsilon(1e-5)); + CHECK(state[12] == + PLApproxComplex(std::complex{-0.560499, -0.310176}) + .epsilon(1e-5)); + CHECK(state[14] == + PLApproxComplex(std::complex{0.0756372, -0.226334}) + .epsilon(1e-5)); + } + + SECTION("Hadamard and IsingZZ and cache manager") { + std::unique_ptr LQsim = std::make_unique(); + + constexpr std::size_t n_qubits = 2; + std::vector Qs = LQsim->AllocateQubits(n_qubits); + + LQsim->StartTapeRecording(); + LQsim->NamedOperation("Hadamard", {}, {Qs[0]}, false); + LQsim->NamedOperation("Hadamard", {}, {Qs[1]}, false); + LQsim->NamedOperation("IsingZZ", {M_PI_4}, {Qs[0], Qs[1]}, false); + LQsim->StopTapeRecording(); + + std::vector> state(1U << LQsim->GetNumQubits()); + DataView, 1> view(state); + LQsim->State(view); + + std::complex c1{0.4619397663, -0.1913417162}; + std::complex c2{0.4619397663, 0.1913417162}; + + CHECK(state[0] == PLApproxComplex(c1).epsilon(1e-5)); + CHECK(state[1] == PLApproxComplex(c2).epsilon(1e-5)); + CHECK(state[2] == PLApproxComplex(c2).epsilon(1e-5)); + CHECK(state[3] == PLApproxComplex(c1).epsilon(1e-5)); + + std::tuple, std::vector> + expected{3, 0, 1, {"Hadamard", "Hadamard", "IsingZZ"}, {}}; + REQUIRE(LQsim->CacheManagerInfo() == expected); + } + + SECTION("Test setStateVector") { + std::unique_ptr LQsim = std::make_unique(); + constexpr std::size_t n_qubits = 2; + std::vector Qs = LQsim->AllocateQubits(n_qubits); + + std::vector> data = {{0.5, 0.5}, {0.0, 0.0}}; + DataView, 1> data_view(data); + std::vector wires = {1}; + LQsim->SetState(data_view, wires); + + std::vector> state(1U << LQsim->GetNumQubits()); + DataView, 1> view(state); + LQsim->State(view); + + std::complex c1{0.5, 0.5}; + std::complex c2{0.0, 0.0}; + CHECK(state[0] == PLApproxComplex(c1).epsilon(1e-5)); + CHECK(state[1] == PLApproxComplex(c2).epsilon(1e-5)); + CHECK(state[2] == PLApproxComplex(c2).epsilon(1e-5)); + CHECK(state[3] == PLApproxComplex(c2).epsilon(1e-5)); + } + + SECTION("Test setBasisState") { + std::unique_ptr LQsim = std::make_unique(); + constexpr std::size_t n_qubits = 1; + std::vector Qs = LQsim->AllocateQubits(n_qubits); + + std::vector data = {0}; + DataView data_view(data); + std::vector wires = {0}; + LQsim->SetBasisState(data_view, wires); + + std::vector> state(1U << LQsim->GetNumQubits()); + DataView, 1> view(state); + LQsim->State(view); + + std::complex c1{1.0, 0.0}; + std::complex c2{0.0, 0.0}; + CHECK(state[0] == PLApproxComplex(c1).epsilon(1e-5)); + CHECK(state[1] == PLApproxComplex(c2).epsilon(1e-5)); + } +} diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/Test_SVDynamicCPU_Allocation.cpp b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/Test_SVDynamicCPU_Allocation.cpp new file mode 100644 index 000000000..f04774779 --- /dev/null +++ b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/Test_SVDynamicCPU_Allocation.cpp @@ -0,0 +1,275 @@ +// Copyright 2022 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. + +#include +#include +#include +#include +#include +#include +#include +#include + +#include "StateVectorLQubitDynamic.hpp" +#include "Util.hpp" +#include "cpu_kernels/GateImplementationsLM.hpp" +#include + +#include "TestHelpers.hpp" + +using namespace Pennylane::LightningQubit; + +TEMPLATE_TEST_CASE("StateVectorLQubitDynamic::getSubsystemPurity /allocation", + "[StateVectorLQubitDynamic]", float, double) { + using PrecisionT = TestType; + + SECTION("Test getSubsystemPurity for a state-vector with RX-RY") { + constexpr size_t num_qubits = 3; + StateVectorLQubitDynamic sv1(num_qubits); + + sv1.applyOperations({"RX", "RY"}, {{0}, {1}}, {false, false}, + {{0.1}, {0.2}}); + + CHECK(sv1.getSubsystemPurity(0) == + approx(std::complex{1, 0})); + CHECK(sv1.getSubsystemPurity(1) == + approx(std::complex{1, 0})); + CHECK(sv1.getSubsystemPurity(2) == + approx(std::complex{1, 0})); + } + + SECTION("Test checkSubsystemPurity for a state-vector with RX-RY") { + constexpr size_t num_qubits = 3; + StateVectorLQubitDynamic sv1(num_qubits); + + sv1.applyOperations({"RX", "RY"}, {{0}, {1}}, {false, false}, + {{0.1}, {0.2}}); + CHECK((sv1.checkSubsystemPurity(0) && sv1.checkSubsystemPurity(1))); + CHECK(sv1.checkSubsystemPurity(2)); + } + + SECTION("Test getSubsystemPurity for a state-vector with CNOT-RY") { + constexpr size_t num_qubits = 3; + StateVectorLQubitDynamic sv1(num_qubits); + + sv1.applyOperations({"CNOT", "RY"}, {{0, 1}, {1}}, {false, false}, + {{}, {0.2}}); + + CHECK(sv1.getSubsystemPurity(0) == + approx(std::complex{1, 0})); + CHECK(sv1.getSubsystemPurity(1) == + approx(std::complex{1, 0})); + CHECK(sv1.checkSubsystemPurity(2)); + } + + SECTION("Test getSubsystemPurity for a custom state-vector") { + constexpr size_t num_qubits = 2; + StateVectorLQubitDynamic sv1(num_qubits); + + std::vector> data{ + {1 / 2, 0}, {1 / 2, 0}, {-1 / 2, 0}, {-1 / 2, 0}}; + + sv1.updateData(data); + + CHECK(sv1.getSubsystemPurity(0) != + approx(std::complex{1, 0})); + CHECK(sv1.getSubsystemPurity(1) != + approx(std::complex{1, 0})); + } +} + +TEMPLATE_TEST_CASE("StateVectorLQubitDynamic::allocateWire /allocation", + "[StateVectorLQubitDynamic]", float, double) { + using PrecisionT = TestType; + + SECTION("applyOperations with released wires") { + constexpr size_t num_qubits = 3; + StateVectorLQubitDynamic sv1(num_qubits); + size_t new_idx = sv1.allocateWire(); + sv1.releaseWire(new_idx); + + CHECK((sv1.getNumQubits() == num_qubits && new_idx == num_qubits)); + + REQUIRE_NOTHROW(sv1.applyOperations( + {"PauliX", "PauliY"}, {{new_idx + 1}, {1}}, {false, false})); + } + + SECTION("Test counting wires for a simple state-vector") { + constexpr size_t num_qubits = 10; + StateVectorLQubitDynamic sv1(num_qubits); + + sv1.releaseWire(1); + sv1.releaseWire(3); + sv1.releaseWire(4); + sv1.releaseWire(6); + sv1.allocateWire(); + sv1.allocateWire(); + sv1.allocateWire(); + sv1.allocateWire(); + + CHECK(sv1.getNumQubits() == 10); + } + + SECTION( + "Test the validity of the shranked state vector in the second half") { + constexpr size_t num_qubits = 4; + StateVectorLQubitDynamic sv1(num_qubits); + std::vector> data(1U << num_qubits, + {0.0, 0.0}); + data[(1U << num_qubits) - 1] = {1.0, 0.0}; + sv1.updateData(data); + + sv1.releaseWire(0); + + CHECK(sv1.getNumQubits() == 3); + } + + SECTION("Test allocation/deallocation of a customed state-vector") { + StateVectorLQubitDynamic sv1(0); + size_t idx_0 = sv1.allocateWire(); // 1, 0 + + std::vector> expected_data{{1.0, 0.0}, + {0.0, 0.0}}; + + CHECK(sv1.getDataVector() == approx(expected_data)); + + sv1.applyOperation("Hadamard", {idx_0}, false, {}); + expected_data[0] = std::complex(0.707107, 0); + expected_data[1] = std::complex(0.707107, 0); + + CHECK(sv1.getDataVector() == approx(expected_data)); + + sv1.allocateWire(); + + std::vector> expected_data_n2{ + {0.707107, 0}, {0.0, 0.0}, {0.707107, 0}, {0.0, 0.0}}; + + CHECK(sv1.getDataVector() == approx(expected_data_n2)); + + sv1.allocateWire(); + + std::vector> expected_data_n3{ + {0.707107, 0}, {0.0, 0.0}, {0.0, 0.0}, {0.0, 0.0}, + {0.707107, 0}, {0.0, 0.0}, {0.0, 0.0}, {0.0, 0.0}}; + + CHECK(sv1.getDataVector() == approx(expected_data_n3)); + + sv1.releaseWire(0); + + CHECK(sv1.getDataVector() == approx(expected_data_n2)); + + sv1.releaseWire(0); + + CHECK(sv1.getDataVector() == approx(expected_data)); + + sv1.releaseWire(0); + + CHECK(sv1.getDataVector()[0] == + approx(std::complex(1.0, 0.0))); + } + + SECTION("Test allocation/deallocation of wires for a state-vector with " + "num_qubits=0") { + StateVectorLQubitDynamic sv1(0); + + std::vector> expected_data{{1, 0}}; + CHECK(sv1.getDataVector() == approx(expected_data)); + + size_t idx_0 = sv1.allocateWire(); + + expected_data.push_back({0, 0}); + CHECK(sv1.getDataVector() == approx(expected_data)); + CHECK(idx_0 == 0); + + sv1.applyOperation("Hadamard", {idx_0}, false, {}); + + StateVectorLQubitDynamic sv2 = sv1; + + size_t new_idx = sv1.allocateWire(); + sv1.applyOperation("RX", {new_idx}, false, {0.3}); + + sv1.releaseWire(0); + CHECK(sv1.getDataVector() == approx(sv2.getDataVector())); + } + + SECTION("Test allocation/deallocation of wires for a state-vector with " + "num_qubits=1") { + constexpr size_t num_qubits = 1; + StateVectorLQubitDynamic sv1(num_qubits); + sv1.applyOperation("Hadamard", {0}, false, {}); + + StateVectorLQubitDynamic sv2 = sv1; + + size_t new_idx = sv1.allocateWire(); + sv1.applyOperation("RX", {new_idx}, false, {0.3}); + + sv1.releaseWire(0); + CHECK(sv1.getDataVector() == approx(sv2.getDataVector())); + } + + SECTION("Test allocation/deallocation of wires for a state-vector with " + "num_qubits=2") { + constexpr size_t num_qubits = 2; + StateVectorLQubitDynamic sv1(num_qubits); + sv1.applyOperations({"RX", "CNOT"}, {{0}, {0, 1}}, {false, false}, + {{0.4}, {}}); + + StateVectorLQubitDynamic sv2 = sv1; + + size_t new_idx = sv1.allocateWire(); + sv1.applyOperation("RX", {new_idx}, false, {0.3}); + + sv1.releaseWire(0); + + CHECK(sv1.getDataVector() == approx(sv2.getDataVector())); + } + + SECTION("Test allocation/deallocation of wires for a state-vector with " + "num_qubits=3") { + constexpr size_t num_qubits = 3; + StateVectorLQubitDynamic sv1(num_qubits); + + sv1.applyOperations({"RX", "SWAP"}, {{0}, {0, 2}}, {false, false}, + {{0.4}, {}}); + + StateVectorLQubitDynamic sv2{num_qubits - 1}; + sv2.applyOperations({"RX", "SWAP"}, {{0}, {0, 1}}, {false, false}, + {{0.4}, {}}); + + sv1.releaseWire(1); + CHECK(sv1.getDataVector() == approx(sv2.getDataVector())); + } + + SECTION("Test allocation/deallocation of wires for a state-vector with " + "num_qubits=4") { + constexpr size_t num_qubits = 4; + StateVectorLQubitDynamic sv1(num_qubits); + + sv1.applyOperations({"RX", "SWAP", "RY", "Hadamard", "RZ", "CNOT"}, + {{0}, {1, 2}, {1}, {3}, {2}, {1, 3}}, + {false, false, false, false, false, false}, + {{0.4}, {}, {0.6}, {}, {0.8}, {}}); + + std::vector> result{ + {0.651289, -0.27536}, + {0.651289, -0.27536}, + }; + + sv1.releaseWire(1); + sv1.releaseWire(1); + sv1.releaseWire(1); + + CHECK(sv1.getDataVector() == approx(result)); + } +} diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/Test_SVDynamicCPU_Core.cpp b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/Test_SVDynamicCPU_Core.cpp new file mode 100644 index 000000000..0dad17e2d --- /dev/null +++ b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/Test_SVDynamicCPU_Core.cpp @@ -0,0 +1,205 @@ +// Copyright 2022 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. + +#include +#include +#include +#include +#include +#include +#include +#include + +#include "StateVectorLQubitDynamic.hpp" +#include "Util.hpp" +#include "cpu_kernels/GateImplementationsLM.hpp" +#include + +#include "TestHelpers.hpp" + +using namespace Pennylane::LightningQubit; + +TEMPLATE_TEST_CASE("StateVectorLQubitDynamic::StateVectorLQubitDynamic", + "[StateVectorLQubitDynamic]", float, double) { + using PrecisionT = TestType; + + SECTION("StateVectorLQubitDynamic") { + CHECK(!std::is_constructible_v>); + } + SECTION("StateVectorLQubitDynamic") { + CHECK(!std::is_constructible_v>); + } + SECTION("StateVectorLQubitDynamic {size_t}") { + CHECK(std::is_constructible_v, + size_t>); + const size_t num_qubits = 4; + StateVectorLQubitDynamic sv(num_qubits); + + CHECK(sv.getNumQubits() == 4); + CHECK(sv.getLength() == 16); + CHECK(sv.getDataVector().size() == 16); + } + SECTION("StateVectorLQubitDynamic {const " + "StateVectorLQubitDynamic&}") { + CHECK(std::is_copy_constructible_v>); + } + SECTION("StateVectorLQubitDynamic " + "{StateVectorLQubitDynamic&&}") { + CHECK(std::is_move_constructible_v>); + } + SECTION("Aligned 256bit statevector") { + const auto memory_model = CPUMemoryModel::Aligned256; + StateVectorLQubitDynamic sv(4, Threading::SingleThread, + memory_model); + /* Even when we allocate 256 bit aligend memory, it is possible that the + * alignment happens to be 512 bit */ + CHECK(((getMemoryModel(sv.getDataVector().data()) == + CPUMemoryModel::Aligned256) || + (getMemoryModel(sv.getDataVector().data()) == + CPUMemoryModel::Aligned512))); + } + + SECTION("Aligned 512bit statevector") { + const auto memory_model = CPUMemoryModel::Aligned512; + StateVectorLQubitDynamic sv(4, Threading::SingleThread, + memory_model); + CHECK((getMemoryModel(sv.getDataVector().data()) == + CPUMemoryModel::Aligned512)); + } +} + +TEMPLATE_TEST_CASE("StateVectorLQubitDynamic::applyMatrix with std::vector", + "[StateVectorLQubitDynamic]", float, double) { + using PrecisionT = TestType; + SECTION("Test wrong matrix size") { + std::vector> m(7, 0.0); + const size_t num_qubits = 4; + StateVectorLQubitDynamic sv(num_qubits); + REQUIRE_THROWS_WITH( + sv.applyMatrix(m, {0, 1}), + Catch::Contains( + "The size of matrix does not match with the given")); + } + + SECTION("Test wrong number of wires") { + std::vector> m(8, 0.0); + const size_t num_qubits = 4; + StateVectorLQubitDynamic sv(num_qubits); + REQUIRE_THROWS_WITH( + sv.applyMatrix(m, {0}), + Catch::Contains( + "The size of matrix does not match with the given")); + } +} + +TEMPLATE_TEST_CASE("StateVectorLQubitDynamic::applyMatrix with a pointer", + "[StateVectorLQubitDynamic]", float, double) { + using PrecisionT = TestType; + SECTION("Test wrong matrix") { + std::vector> m(8, 0.0); + const size_t num_qubits = 4; + StateVectorLQubitDynamic sv(num_qubits); + REQUIRE_THROWS_WITH(sv.applyMatrix(m.data(), {}), + Catch::Contains("must be larger than 0")); + } + + SECTION("Test with different number of wires") { + std::default_random_engine re{1337}; + const size_t num_qubits = 5; + for (size_t num_wires = 1; num_wires < num_qubits; num_wires++) { + StateVectorLQubitDynamic sv1(num_qubits); + StateVectorLQubitDynamic sv2(num_qubits); + + std::vector wires(num_wires); + std::iota(wires.begin(), wires.end(), 0); + + const auto m = + Pennylane::Util::randomUnitary(re, num_wires); + sv1.applyMatrix(m, wires); + Gates::GateImplementationsLM::applyMultiQubitOp( + sv2.getData(), num_qubits, m.data(), wires, false); + CHECK(sv1.getDataVector() == + approx(sv2.getDataVector()).margin(PrecisionT{1e-5})); + } + } +} + +TEMPLATE_TEST_CASE("StateVectorLQubitDynamic::applyOperations", + "[StateVectorLQubitDynamic]", float, double) { + using PrecisionT = TestType; + + std::mt19937 re{1337}; + + SECTION("Test invalid arguments without params") { + const size_t num_qubits = 4; + StateVectorLQubitDynamic sv(num_qubits); + REQUIRE_THROWS_WITH( + sv.applyOperations({"PauliX", "PauliY"}, {{0}}, {false, false}), + Catch::Contains("must all be equal")); // invalid wires + REQUIRE_THROWS_WITH( + sv.applyOperations({"PauliX", "PauliY"}, {{0}, {1}}, {false}), + Catch::Contains("must all be equal")); // invalid inverse + } + + SECTION("applyOperations without params works as expected") { + const size_t num_qubits = 3; + StateVectorLQubitDynamic sv1(num_qubits); + + sv1.updateData(Pennylane::Util::createRandomStateVectorData( + re, num_qubits)); + StateVectorLQubitDynamic sv2 = sv1; + + sv1.applyOperations({"PauliX", "PauliY"}, {{0}, {1}}, {false, false}); + + sv2.applyOperation("PauliX", {0}, false); + sv2.applyOperation("PauliY", {1}, false); + + CHECK(sv1.getDataVector() == approx(sv2.getDataVector())); + } + + SECTION("Test invalid arguments with params") { + const size_t num_qubits = 4; + StateVectorLQubitDynamic sv(num_qubits); + REQUIRE_THROWS_WITH( + sv.applyOperations({"RX", "RY"}, {{0}}, {false, false}, + {{0.0}, {0.0}}), + Catch::Contains("must all be equal")); // invalid wires + REQUIRE_THROWS_WITH( + sv.applyOperations({"RX", "RY"}, {{0}, {1}}, {false}, + {{0.0}, {0.0}}), + Catch::Contains("must all be equal")); // invalid inverse + + REQUIRE_THROWS_WITH( + sv.applyOperations({"RX", "RY"}, {{0}, {1}}, {false, false}, + {{0.0}}), + Catch::Contains("must all be equal")); // invalid params + } + + SECTION("applyOperations with params works as expected") { + const size_t num_qubits = 3; + StateVectorLQubitDynamic sv1(num_qubits); + + sv1.updateData(Pennylane::Util::createRandomStateVectorData( + re, num_qubits)); + StateVectorLQubitDynamic sv2 = sv1; + + sv1.applyOperations({"RX", "RY"}, {{0}, {1}}, {false, false}, + {{0.1}, {0.2}}); + + sv2.applyOperation("RX", {0}, false, {0.1}); + sv2.applyOperation("RY", {1}, false, {0.2}); + + CHECK(sv1.getDataVector() == approx(sv2.getDataVector())); + } +} diff --git a/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/runner_lightning_qubit_catalyst.cpp b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/runner_lightning_qubit_catalyst.cpp new file mode 100644 index 000000000..4ed06df1f --- /dev/null +++ b/pennylane_lightning/core/src/simulators/lightning_qubit/catalyst/tests/runner_lightning_qubit_catalyst.cpp @@ -0,0 +1,2 @@ +#define CATCH_CONFIG_MAIN +#include diff --git a/pennylane_lightning/lightning_qubit/lightning_qubit.py b/pennylane_lightning/lightning_qubit/lightning_qubit.py index abf080978..1cc2df3ec 100644 --- a/pennylane_lightning/lightning_qubit/lightning_qubit.py +++ b/pennylane_lightning/lightning_qubit/lightning_qubit.py @@ -15,6 +15,8 @@ This module contains the LightningQubit class, a PennyLane simulator device that interfaces with C++ for fast linear algebra calculations. """ +import os +import sys from dataclasses import replace from functools import reduce from pathlib import Path @@ -560,3 +562,45 @@ def simulate( state.reset_state() final_state = state.get_final_state(circuit) return self.LightningMeasurements(final_state, **mcmc).measure_final_state(circuit) + + @staticmethod + def get_c_interface(): + """Returns a tuple consisting of the device name, and + the location to the shared object with the C/C++ device implementation. + """ + + # The shared object file extension varies depending on the underlying operating system + file_extension = "" + OS = sys.platform + if OS == "linux": + file_extension = ".so" + elif OS == "darwin": + file_extension = ".dylib" + else: + raise RuntimeError( + f"'LightningSimulator' shared library not available for '{OS}' platform" + ) # pragma: no cover + + lib_name = "liblightning_qubit_catalyst" + file_extension + package_root = Path(__file__).parent + + # The absolute path of the plugin shared object varies according to the installation mode. + + # Wheel mode: + # Fixed location at the root of the project + wheel_mode_location = package_root.parent / lib_name + if wheel_mode_location.is_file(): + return "LightningSimulator", wheel_mode_location.as_posix() + + # Editable mode: + # The build directory contains a folder which varies according to the platform: + # lib.--" + # To avoid mismatching the folder name, we search for the shared object instead. + # TODO: locate where the naming convention of the folder is decided and replicate it here. + editable_mode_path = package_root.parent.parent / "build" + for path, _, files in os.walk(editable_mode_path): + if lib_name in files: + lib_location = (Path(path) / lib_name).as_posix() + return "LightningSimulator", lib_location + + raise RuntimeError("'LightningSimulator' shared library not found") # pragma: no cover diff --git a/setup.py b/setup.py index e0bf04ecc..30ba88d6b 100644 --- a/setup.py +++ b/setup.py @@ -162,7 +162,7 @@ def build_extension(self, ext: CMakeExtension): destination = os.path.join(os.getcwd(), f"build_{backend}") shutil.copy(source, destination) - if backend in ("lightning_kokkos"): + if backend in ("lightning_kokkos", "lightning_qubit"): if platform.system() in ["Linux", "Darwin"]: shared_lib_ext = {"Linux": ".so", "Darwin": ".dylib"}[platform.system()] source = os.path.join(f"{extdir}", f"lib{backend}_catalyst{shared_lib_ext}") diff --git a/tests/test_device.py b/tests/test_device.py index 50b631862..302eab6b6 100644 --- a/tests/test_device.py +++ b/tests/test_device.py @@ -168,3 +168,49 @@ def test_supported_macos_platform_kokkos(): assert dev_name == "LightningKokkosSimulator" assert "liblightning_kokkos_catalyst.dylib" in shared_lib_name + + +@pytest.mark.skipif( + (device_name != "lightning.qubit" or sys.platform != "win32"), + reason="This test is for LQ under Windows only.", +) +def test_unsupported_windows_platform_qubit(): + """Test unsupported Windows platform for LQ.""" + + dev = qml.device(device_name, wires=0) + + with pytest.raises( + RuntimeError, + match="'LightningSimulator' shared library not available for 'win32' platform", + ): + dev.get_c_interface() + + +@pytest.mark.skipif( + (device_name != "lightning.qubit" or sys.platform != "linux"), + reason="This test is for LQ under Linux only.", +) +def test_supported_linux_platform_qubit(): + """Test supported Linux platform for LQ.""" + + dev = qml.device(device_name, wires=0) + + dev_name, shared_lib_name = dev.get_c_interface() + + assert dev_name == "LightningSimulator" + assert "liblightning_qubit_catalyst.so" in shared_lib_name + + +@pytest.mark.skipif( + (device_name != "lightning.qubit" or sys.platform != "darwin"), + reason="This test is for LQ under MacOS only.", +) +def test_supported_macos_platform_qubit(): + """Test supported MacOS platform for LQ.""" + + dev = qml.device(device_name, wires=0) + + dev_name, shared_lib_name = dev.get_c_interface() + + assert dev_name == "LightningSimulator" + assert "liblightning_qubit_catalyst.dylib" in shared_lib_name