From 143827f5f181c297cc3945db4bf7a34ed7e482de Mon Sep 17 00:00:00 2001 From: Sachin Pisal Date: Tue, 14 Jan 2025 15:12:20 -0800 Subject: [PATCH] * Adding helpers functionlity * Aggregating parameters * Extracting documentation * Extracting positional and keyword arguments * Generating all quantum states for given degrees and dimensions * Permuting a given matrix * Canonicalizing degrees * Adding test cases for the above functions * Formatting Signed-off-by: Sachin Pisal --- runtime/cudaq/dynamics/CMakeLists.txt | 4 +- runtime/cudaq/dynamics/helpers.cpp | 135 +++++++++++++ runtime/cudaq/helpers.h | 42 ++++ runtime/cudaq/runge_kutta_integrator.h | 41 ++-- runtime/cudaq/runge_kutta_time_stepper.h | 24 +-- unittests/CMakeLists.txt | 1 + unittests/dynamics/runge_kutta_test_helpers.h | 4 +- unittests/dynamics/test_helpers.cpp | 189 ++++++++++++++++++ .../dynamics/test_runge_kutta_integrator.cpp | 134 ++++++------- .../test_runge_kutta_time_stepper.cpp | 152 +++++++------- 10 files changed, 549 insertions(+), 177 deletions(-) create mode 100644 runtime/cudaq/dynamics/helpers.cpp create mode 100644 runtime/cudaq/helpers.h create mode 100644 unittests/dynamics/test_helpers.cpp diff --git a/runtime/cudaq/dynamics/CMakeLists.txt b/runtime/cudaq/dynamics/CMakeLists.txt index 9709cd9a71..d608aba2c3 100644 --- a/runtime/cudaq/dynamics/CMakeLists.txt +++ b/runtime/cudaq/dynamics/CMakeLists.txt @@ -1,5 +1,5 @@ # ============================================================================ # -# Copyright (c) 2022 - 2024 NVIDIA Corporation & Affiliates. # +# Copyright (c) 2022 - 2025 NVIDIA Corporation & Affiliates. # # All rights reserved. # # # # This source code and the accompanying materials are made available under # @@ -11,7 +11,7 @@ set (CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Wno-ctad-maybe-unsupported") set(INTERFACE_POSITION_INDEPENDENT_CODE ON) set(CUDAQ_OPS_SRC - scalar_operators.cpp elementary_operators.cpp product_operators.cpp operator_sum.cpp schedule.cpp definition.cpp + scalar_operators.cpp elementary_operators.cpp product_operators.cpp operator_sum.cpp schedule.cpp definition.cpp helpers.cpp ) add_library(${LIBRARY_NAME} SHARED ${CUDAQ_OPS_SRC}) diff --git a/runtime/cudaq/dynamics/helpers.cpp b/runtime/cudaq/dynamics/helpers.cpp new file mode 100644 index 0000000000..ad86a2c26b --- /dev/null +++ b/runtime/cudaq/dynamics/helpers.cpp @@ -0,0 +1,135 @@ +/******************************************************************************* + * Copyright (c) 2022 - 2025 NVIDIA Corporation & Affiliates. * + * All rights reserved. * + * * + * This source code and the accompanying materials are made available under * + * the terms of the Apache License 2.0 which accompanies this distribution. * + ******************************************************************************/ + +#include "cudaq/helpers.h" +#include +#include + +namespace cudaq { +// Aggregate parameters from multiple mappings. +std::map OperatorHelpers::aggregate_parameters(const std::vector> ¶meter_mappings) { + std::map parameter_descriptions; + + for (const auto &descriptions : parameter_mappings) { + for (const auto &[key, new_desc] : descriptions) { + if (!parameter_descriptions[key].empty() && !new_desc.empty()) { + parameter_descriptions[key] += "\n---\n" + new_desc; + } else { + parameter_descriptions[key] = new_desc; + } + } + } + + return parameter_descriptions; +} + +// Extract documentation for a specific parameter from docstring. +std::string OperatorHelpers::parameter_docs(const std::string ¶m_name, const std::string &docs) { + if (param_name.empty() || docs.empty()) { + return ""; + } + + try { + std::regex keyword_pattern(R"(^\s*(Arguments|Args):\s*$)", std::regex::multiline); + std::regex param_pattern(R"(^\s*)" + param_name + R"(\s*(\(.*\))?:\s*(.*)$)", std::regex::multiline); + + std::smatch match; + std::sregex_iterator it(docs.begin(), docs.end(), keyword_pattern); + std::sregex_iterator end; + + if (it != end) { + std::string params_section = docs.substr(it->position() + it->length()); + if(std::regex_search(params_section, match, param_pattern)) { + std::string param_docs = match.str(2); + return std::regex_replace(param_docs, std::regex(R"(\s+)"), " "); + } + } + } catch (...) { + return ""; + } + + return ""; +} + +// Extract positional arguments and keyword-only arguments. +std::pair, std::map> OperatorHelpers::args_from_kwargs(const std::map &kwargs, + const std::vector &required_args, const std::vector &kwonly_args) { + std::vector extracted_args; + std::map kwonly_dict; + + for (const auto &arg : required_args) { + if (kwargs.count(arg)) { + extracted_args.push_back(kwargs.at(arg)); + } else { + throw std::invalid_argument("Missing required argument: " + arg); + } + } + + for (const auto &arg : kwonly_args) { + if (kwargs.count(arg)) { + kwonly_dict[arg] = kwargs.at(arg); + } + } + + return {extracted_args, kwonly_dict}; +} + +// Generate all possible quantum states for given degrees and dimensions +std::vector OperatorHelpers::generate_all_states(const std::vector °rees, const std::map &dimensions) { + if (degrees.empty()) { + return {}; + } + + std::vector> states; + for (int state = 0; state < dimensions.at(degrees[0]); state++) { + states.push_back({std::to_string(state)}); + } + + for (size_t i = 1; i < degrees.size(); i++) { + std::vector> new_states; + for (const auto ¤t : states) { + for (int state = 0; state < dimensions.at(degrees[i]); state++) { + auto new_entry = current; + new_entry.push_back(std::to_string(state)); + new_states.push_back(new_entry); + } + } + states = new_states; + } + + std::vector result; + for (const auto &state : states) { + std::ostringstream joined; + for (const auto &s : state) { + joined << s; + } + result.push_back(joined.str()); + } + return result; +} + +// Permute a given eigen matrix +void OperatorHelpers::permute_matrix(Eigen::MatrixXcd &matrix, const std::vector &permutation) { + Eigen::MatrixXcd permuted_matrix(matrix.rows(), matrix.cols()); + + for (size_t i = 0; i < permutation.size(); i++) { + for (size_t j = 0; j < permutation.size(); j++) { + permuted_matrix(i, j) = matrix(permutation[i], permutation[j]); + } + } + + matrix = permuted_matrix; +} + +// Canonicalize degrees by sorting in descending order +std::vector OperatorHelpers::canonicalize_degrees(const std::vector °rees) { + std::vector sorted_degrees = degrees; + std::sort(sorted_degrees.rbegin(), sorted_degrees.rend()); + return sorted_degrees; +} +} \ No newline at end of file diff --git a/runtime/cudaq/helpers.h b/runtime/cudaq/helpers.h new file mode 100644 index 0000000000..488ff3bf0c --- /dev/null +++ b/runtime/cudaq/helpers.h @@ -0,0 +1,42 @@ +/****************************************************************-*- C++ -*-**** + * Copyright (c) 2022 - 2025 NVIDIA Corporation & Affiliates. * + * All rights reserved. * + * * + * This source code and the accompanying materials are made available under * + * the terms of the Apache License 2.0 which accompanies this distribution. * + ******************************************************************************/ + +#pragma once + +#include +#include +#include +#include +#include +#include +#include +#include + +namespace cudaq { +class OperatorHelpers { +public: + // Aggregate parameters from multiple mappings. + static std::map aggregate_parameters(const std::vector> ¶meter_mappings); + + // Extract documentation for a specific parameter from docstring. + static std::string parameter_docs(const std::string ¶m_name, const std::string &docs); + + // Extract positional arguments and keyword-only arguments. + static std::pair, std::map> args_from_kwargs(const std::map &kwargs, + const std::vector &required_args, const std::vector &kwonly_args); + + // Generate all possible quantum states for given degrees and dimensions. + static std::vector generate_all_states(const std::vector °rees, const std::map &dimensions); + + // Permute a given Eigen matrix. + static void permute_matrix(Eigen::MatrixXcd &matrix, const std::vector &permutation); + + // Canonicalize degrees by sorting in descending order. + static std::vector canonicalize_degrees(const std::vector °rees); +}; +} \ No newline at end of file diff --git a/runtime/cudaq/runge_kutta_integrator.h b/runtime/cudaq/runge_kutta_integrator.h index 41d3da1715..9914258386 100644 --- a/runtime/cudaq/runge_kutta_integrator.h +++ b/runtime/cudaq/runge_kutta_integrator.h @@ -16,31 +16,32 @@ namespace cudaq { template class RungeKuttaIntegrator : public BaseIntegrator { public: - using DerivativeFunction = std::function; + using DerivativeFunction = std::function; - explicit RungeKuttaIntegrator(DerivativeFunction f) : stepper(std::make_shared>(f)) {} + explicit RungeKuttaIntegrator(DerivativeFunction f) + : stepper(std::make_shared>(f)) {} - // Initializes the integrator - void post_init() override { - if (!this->stepper) { - throw std::runtime_error("Time stepper is not set"); - } + // Initializes the integrator + void post_init() override { + if (!this->stepper) { + throw std::runtime_error("Time stepper is not set"); } + } - // Advances the system's state from current time to `t` - void integrate(double target_t) override { - if (!this->schedule || !this->hamiltonian) { - throw std::runtime_error("System is not properly set!"); - } - - while (this->t < target_t) { - stepper->compute(this->state, this->t); - // Time step size - this->t += 0.01; - } + // Advances the system's state from current time to `t` + void integrate(double target_t) override { + if (!this->schedule || !this->hamiltonian) { + throw std::runtime_error("System is not properly set!"); } + while (this->t < target_t) { + stepper->compute(this->state, this->t); + // Time step size + this->t += 0.01; + } + } + private: - std::shared_ptr> stepper; + std::shared_ptr> stepper; }; -} \ No newline at end of file +} // namespace cudaq \ No newline at end of file diff --git a/runtime/cudaq/runge_kutta_time_stepper.h b/runtime/cudaq/runge_kutta_time_stepper.h index 7718251a07..1dcd1f69cc 100644 --- a/runtime/cudaq/runge_kutta_time_stepper.h +++ b/runtime/cudaq/runge_kutta_time_stepper.h @@ -13,21 +13,21 @@ namespace cudaq { template class RungeKuttaTimeStepper : public BaseTimeStepper { public: - using DerivativeFunction = std::function; + using DerivativeFunction = std::function; - RungeKuttaTimeStepper(DerivativeFunction f) : derivativeFunc(f) {} + RungeKuttaTimeStepper(DerivativeFunction f) : derivativeFunc(f) {} - void compute(TState &state, double t, double dt = 0.01) override { - // 4th order Runge-Kutta method - TState k1 = derivativeFunc(state, t); - TState k2 = derivativeFunc(state + (dt / 2.0) * k1, t + dt / 2.0); - TState k3 = derivativeFunc(state + (dt / 2.0) * k2, t + dt / 2.0); - TState k4 = derivativeFunc(state + dt * k3, t + dt); + void compute(TState &state, double t, double dt = 0.01) override { + // 4th order Runge-Kutta method + TState k1 = derivativeFunc(state, t); + TState k2 = derivativeFunc(state + (dt / 2.0) * k1, t + dt / 2.0); + TState k3 = derivativeFunc(state + (dt / 2.0) * k2, t + dt / 2.0); + TState k4 = derivativeFunc(state + dt * k3, t + dt); - state = state + (dt / 6.0) * (k1 + 2 * k2 + 2 * k3 + k4); - } + state = state + (dt / 6.0) * (k1 + 2 * k2 + 2 * k3 + k4); + } private: - DerivativeFunction derivativeFunc; + DerivativeFunction derivativeFunc; }; -} \ No newline at end of file +} // namespace cudaq \ No newline at end of file diff --git a/unittests/CMakeLists.txt b/unittests/CMakeLists.txt index 5b18824956..afc452fd19 100644 --- a/unittests/CMakeLists.txt +++ b/unittests/CMakeLists.txt @@ -51,6 +51,7 @@ set(CUDAQ_RUNTIME_TEST_SOURCES dynamics/product_operators_arithmetic.cpp dynamics/test_runge_kutta_time_stepper.cpp dynamics/test_runge_kutta_integrator.cpp + dynamics/test_helpers.cpp ) # Make it so we can get function symbols diff --git a/unittests/dynamics/runge_kutta_test_helpers.h b/unittests/dynamics/runge_kutta_test_helpers.h index 28ec5c7723..4f93ffa242 100644 --- a/unittests/dynamics/runge_kutta_test_helpers.h +++ b/unittests/dynamics/runge_kutta_test_helpers.h @@ -15,10 +15,10 @@ using TestState = double; // Simple derivative function: dx/dt = -x (exponential decay) inline TestState simple_derivative(const TestState &state, double t) { - return -state; + return -state; } // A complex function: dx/dt = sin(t) inline TestState sine_derivative(const TestState &state, double t) { - return std::sin(t); + return std::sin(t); } diff --git a/unittests/dynamics/test_helpers.cpp b/unittests/dynamics/test_helpers.cpp new file mode 100644 index 0000000000..9e471a29bd --- /dev/null +++ b/unittests/dynamics/test_helpers.cpp @@ -0,0 +1,189 @@ +/******************************************************************************* + * Copyright (c) 2022 - 2025 NVIDIA Corporation & Affiliates. * + * All rights reserved. * + * * + * This source code and the accompanying materials are made available under * + * the terms of the Apache License 2.0 which accompanies this distribution. * + ******************************************************************************/ + +#include "cudaq/helpers.h" +#include +#include + +using namespace cudaq; + +TEST(OperatorHelpersTest, AggregateParameters_MultipleMappings) { + std::vector> mappings = { + {{"alpha", "Parameter A"}, {"beta", "Parameter B"}}, {{"alpha", "Updated Parameter A"}, {"gamma", "New Parameter"}} + }; + + auto result = OperatorHelpers::aggregate_parameters(mappings); + + EXPECT_EQ(result["alpha"], "Parameter A\n---\nUpdated Parameter A"); + EXPECT_EQ(result["beta"], "Parameter B"); + EXPECT_EQ(result["gamma"], "New Parameter"); +} + +TEST(OperatorHelpersTest, AggregateParameters_EmptyMappings) { + std::vector> mappings; + auto result = OperatorHelpers::aggregate_parameters(mappings); + + EXPECT_TRUE(result.empty()); +} + +TEST(OperatorHelpersTest, ParameterDocs_ValidExtraction) { + std::string docstring = + "Description of function.\n" + "Arguments:\n" + " alpha (float): The first parameter.\n" + " beta (int): The second parameter."; + + auto result = OperatorHelpers::parameter_docs("alpha", docstring); + EXPECT_EQ(result, "The first parameter."); + + result = OperatorHelpers::parameter_docs("beta", docstring); + EXPECT_EQ(result, "The second parameter."); +} + +TEST(OperatorHelpersTest, ParameterDocs_InvalidParam) { + std::string docstring = + "Description of function.\n" + "Arguments:\n" + " alpha (float): The first parameter.\n" + " beta (int): The second parameter."; + + auto result = OperatorHelpers::parameter_docs("gamma", docstring); + EXPECT_EQ(result, ""); +} + +TEST(OperatorHelpersTest, ParameterDocs_EmptyDocString) { + std::string docstring = ""; + auto result = OperatorHelpers::parameter_docs("alpha", docstring); + EXPECT_EQ(result, ""); +} + +TEST(OperatorHelpersTest, GenerateAllStates_TwoQubits) { + std::vector degrees = {0, 1}; + std::map dimensions = {{0, 2}, {1, 2}}; + + auto states = OperatorHelpers::generate_all_states(degrees, dimensions); + std::vector expected_states = {"00", "01", "10", "11"}; + + EXPECT_EQ(states, expected_states); +} + +TEST(OperatorHelpersTest, GenerateAllStates_ThreeQubits) { + std::vector degrees = {0, 1, 2}; + std::map dimensions = {{0, 2}, {1, 2}, {2, 2}}; + + auto states = OperatorHelpers::generate_all_states(degrees, dimensions); + std::vector expected_states = {"000", "001", "010", "011", "100", "101", "110", "111"}; + + EXPECT_EQ(states, expected_states); +} + +TEST(OperatorHelpersTest, GenerateAllStates_EmptyDegrees) { + std::vector degrees; + std::map dimensions; + + auto states = OperatorHelpers::generate_all_states(degrees, dimensions); + EXPECT_TRUE(states.empty()); +} + +TEST(OperatorHelpersTest, GenerateAllStates_MissingDegreesInMap) { + std::vector degrees = {0, 1, 2}; + std::map dimensions = {{0, 2}, {1, 2}}; + + EXPECT_THROW(OperatorHelpers::generate_all_states(degrees, dimensions), std::out_of_range); +} + +TEST(OperatorHelpersTest, PermuteMatrix_SingleSwap) { + Eigen::MatrixXcd matrix(2, 2); + matrix << 1, 2, + 3, 4; + + // Swap rows and columns + std::vector permutation = {1, 0}; + + OperatorHelpers::permute_matrix(matrix, permutation); + + Eigen::MatrixXcd expected(2, 2); + expected << 4, 3, + 2, 1; + + EXPECT_EQ(matrix, expected); +} + +TEST(OperatorHelpersTest, PermuteMatrix_IdentityPermutation) { + Eigen::MatrixXcd matrix(3, 3); + matrix << 1, 2, 3, + 4, 5, 6, + 7, 8, 9; + + // No change + std::vector permutation = {0, 1, 2}; + + OperatorHelpers::permute_matrix(matrix, permutation); + + Eigen::MatrixXcd expected(3, 3); + expected << 1, 2, 3, + 4, 5, 6, + 7, 8, 9; + + EXPECT_EQ(matrix, expected); +} + +TEST(OperatorHelpersTest, CanonicalizeDegrees_SortedDescending) { + std::vector degrees = {3, 1, 2}; + auto sorted = OperatorHelpers::canonicalize_degrees(degrees); + EXPECT_EQ(sorted, (std::vector{3, 2, 1})); +} + +TEST(OperatorHelpersTest, CanonicalizeDegrees_AlreadySorted) { + std::vector degrees = {5, 4, 3, 2, 1}; + auto sorted = OperatorHelpers::canonicalize_degrees(degrees); + EXPECT_EQ(sorted, (std::vector{5, 4, 3, 2, 1})); +} + +TEST(OperatorHelpersTest, CanonicalizeDegrees_EmptyList) { + std::vector degrees; + auto sorted = OperatorHelpers::canonicalize_degrees(degrees); + EXPECT_TRUE(sorted.empty()); +} + +TEST(OperatorHelpersTest, ArgsFromKwargs_ValidArgs) { + std::map kwargs = { + {"alpha", "0.5"}, + {"beta", "1.0"}, + {"gamma", "2.0"} + }; + + std::vector required_args = {"alpha", "beta"}; + std::vector kwonly_args = {"gamma"}; + + auto [args, kwonly] = OperatorHelpers::args_from_kwargs(kwargs, required_args, kwonly_args); + + EXPECT_EQ(args.size(), 2); + EXPECT_EQ(args[0], "0.5"); + EXPECT_EQ(args[1], "1.0"); + + EXPECT_EQ(kwonly.size(), 1); + EXPECT_EQ(kwonly["gamma"], "2.0"); +} + + +TEST(OperatorHelpersTest, ArgsFromKwargs_MissingRequiredArgs) { + std::map kwargs = { + {"beta", "1.0"}, + {"gamma", "2.0"} + }; + + std::vector required_args = {"alpha", "beta"}; + std::vector kwonly_args = {"gamma"}; + + EXPECT_THROW(OperatorHelpers::args_from_kwargs(kwargs, required_args, kwonly_args), std::invalid_argument); +} + + + + diff --git a/unittests/dynamics/test_runge_kutta_integrator.cpp b/unittests/dynamics/test_runge_kutta_integrator.cpp index 431ea1457a..c75a7c8d6d 100644 --- a/unittests/dynamics/test_runge_kutta_integrator.cpp +++ b/unittests/dynamics/test_runge_kutta_integrator.cpp @@ -1,88 +1,87 @@ // /******************************************************************************* -// * Copyright (c) 2022 - 2025 NVIDIA Corporation & Affiliates. * -// * All rights reserved. * -// * * -// * This source code and the accompanying materials are made available under * -// * the terms of the Apache License 2.0 which accompanies this distribution. * +// * Copyright (c) 2022 - 2025 NVIDIA Corporation & Affiliates. * +// * All rights reserved. * +// * * +// * This source code and the accompanying materials are made available under * +// * the terms of the Apache License 2.0 which accompanies this distribution. * // ******************************************************************************/ -#include "runge_kutta_test_helpers.h" #include "cudaq/runge_kutta_integrator.h" +#include "runge_kutta_test_helpers.h" +#include #include #include -#include using namespace cudaq; // Test fixture class class RungeKuttaIntegratorTest : public ::testing::Test { protected: - RungeKuttaIntegrator *integrator; - std::shared_ptr schedule; - std::shared_ptr hamiltonian; + RungeKuttaIntegrator *integrator; + std::shared_ptr schedule; + std::shared_ptr hamiltonian; - void SetUp() override { - integrator = new RungeKuttaIntegrator(simple_derivative); - // Initial state and time - integrator->set_state(1.0, 0.0); + void SetUp() override { + integrator = new RungeKuttaIntegrator(simple_derivative); + // Initial state and time + integrator->set_state(1.0, 0.0); - // A simple step sequence for the schedule - std::vector> steps = {0.1, 0.2, 0.3, 0.4, 0.5}; + // A simple step sequence for the schedule + std::vector> steps = {0.1, 0.2, 0.3, 0.4, 0.5}; - // Dummy parameters - std::vector parameters = {"param1"}; + // Dummy parameters + std::vector parameters = {"param1"}; - // A simple parameter function - auto value_function = [](const std::string ¶m, const std::complex &step) { - return step; - }; + // A simple parameter function + auto value_function = [](const std::string ¶m, + const std::complex &step) { return step; }; - // A valid schedule instance - schedule = std::make_shared(steps, parameters, value_function); + // A valid schedule instance + schedule = std::make_shared(steps, parameters, value_function); - // A simple hamiltonian as an operator_sum - hamiltonian = std::make_shared(); - *hamiltonian += 0.5 * elementary_operator::identity(0); - *hamiltonian += 0.5 * elementary_operator::number(0); + // A simple hamiltonian as an operator_sum + hamiltonian = std::make_shared(); + *hamiltonian += 0.5 * elementary_operator::identity(0); + *hamiltonian += 0.5 * elementary_operator::number(0); - // System with valid components - integrator->set_system({{0, 2}}, schedule, hamiltonian); - } + // System with valid components + integrator->set_system({{0, 2}}, schedule, hamiltonian); + } - void TearDown() override { - delete integrator; - } + void TearDown() override { delete integrator; } }; // Basic integration TEST_F(RungeKuttaIntegratorTest, BasicIntegration) { - integrator->integrate(1.0); + integrator->integrate(1.0); - // Expected result: x(t) = e^(-t) - double expected = std::exp(-1.0); + // Expected result: x(t) = e^(-t) + double expected = std::exp(-1.0); - EXPECT_NEAR(integrator->get_state().second, expected, 1e-3) << "Basic Runge-Kutta integration failed!"; + EXPECT_NEAR(integrator->get_state().second, expected, 1e-3) + << "Basic Runge-Kutta integration failed!"; } // Time evolution TEST_F(RungeKuttaIntegratorTest, TimeEvolution) { - integrator->integrate(2.0); + integrator->integrate(2.0); - double expected = 2.0; + double expected = 2.0; - EXPECT_NEAR(integrator->get_state().first, expected, 1e-5) << "Integrator did not correctly update time!"; + EXPECT_NEAR(integrator->get_state().first, expected, 1e-5) + << "Integrator did not correctly update time!"; } // Large step size TEST_F(RungeKuttaIntegratorTest, LargeStepSize) { - integrator->integrate(5.0); + integrator->integrate(5.0); - double expected = std::exp(-5.0); + double expected = std::exp(-5.0); - EXPECT_NEAR(integrator->get_state().second, expected, 1e-1) << "Runge-Kutta integration failed for large step size!!"; + EXPECT_NEAR(integrator->get_state().second, expected, 1e-1) + << "Runge-Kutta integration failed for large step size!!"; } - // // Integrating Sine function // TEST_F(RungeKuttaIntegratorTest, SineFunction) { // integrator = new RungeKuttaIntegrator(sine_derivative); @@ -93,40 +92,41 @@ TEST_F(RungeKuttaIntegratorTest, LargeStepSize) { // double expected = std::cos(M_PI / 2); -// EXPECT_NEAR(integrator->get_state().second, expected, 1e-3) << "Runge-Kutta integration for sine function failed!"; +// EXPECT_NEAR(integrator->get_state().second, expected, 1e-3) << +// "Runge-Kutta integration for sine function failed!"; // } - // Small step size TEST_F(RungeKuttaIntegratorTest, SmallStepIntegration) { - integrator->set_state(1.0, 0.0); - integrator->set_system({{0, 2}}, schedule, hamiltonian); + integrator->set_state(1.0, 0.0); + integrator->set_system({{0, 2}}, schedule, hamiltonian); - double step_size = 0.001; - while (integrator->get_state().first < 1.0) { - integrator->integrate(integrator->get_state().first + step_size); - } + double step_size = 0.001; + while (integrator->get_state().first < 1.0) { + integrator->integrate(integrator->get_state().first + step_size); + } - double expected = std::exp(-1.0); + double expected = std::exp(-1.0); - EXPECT_NEAR(integrator->get_state().second, expected, 5e-4) << "Runge-Kutta integration for small step size failed!"; + EXPECT_NEAR(integrator->get_state().second, expected, 5e-4) + << "Runge-Kutta integration for small step size failed!"; } - // Large step size TEST_F(RungeKuttaIntegratorTest, LargeStepIntegration) { - integrator->set_state(1.0, 0.0); - integrator->set_system({{0, 2}}, schedule, hamiltonian); + integrator->set_state(1.0, 0.0); + integrator->set_system({{0, 2}}, schedule, hamiltonian); - double step_size = 0.5; - double t = 0.0; - double target_t = 1.0; - while (t < target_t) { - integrator->integrate(std::min(t + step_size, target_t)); - t += step_size; - } + double step_size = 0.5; + double t = 0.0; + double target_t = 1.0; + while (t < target_t) { + integrator->integrate(std::min(t + step_size, target_t)); + t += step_size; + } - double expected = std::exp(-1.0); + double expected = std::exp(-1.0); - EXPECT_NEAR(integrator->get_state().second, expected, 1e-2) << "Runge-Kutta integration for large step size failed!"; + EXPECT_NEAR(integrator->get_state().second, expected, 1e-2) + << "Runge-Kutta integration for large step size failed!"; } diff --git a/unittests/dynamics/test_runge_kutta_time_stepper.cpp b/unittests/dynamics/test_runge_kutta_time_stepper.cpp index bc32bb587c..4c4c7b7588 100644 --- a/unittests/dynamics/test_runge_kutta_time_stepper.cpp +++ b/unittests/dynamics/test_runge_kutta_time_stepper.cpp @@ -6,75 +6,80 @@ * the terms of the Apache License 2.0 which accompanies this distribution. * ******************************************************************************/ -#include "runge_kutta_test_helpers.h" #include "cudaq/runge_kutta_time_stepper.h" +#include "runge_kutta_test_helpers.h" +#include #include #include -#include // Test fixture class class RungeKuttaTimeStepperTest : public ::testing::Test { protected: - std::shared_ptr> stepper; + std::shared_ptr> stepper; - void SetUp() override { - stepper = std::make_shared>(simple_derivative); - } + void SetUp() override { + stepper = std::make_shared>( + simple_derivative); + } }; // Single step integration TEST_F(RungeKuttaTimeStepperTest, SingleStep) { - // Initial values - double state = 1.0; - double t = 0.0; - double dt = 0.1; + // Initial values + double state = 1.0; + double t = 0.0; + double dt = 0.1; - stepper->compute(state, t, dt); + stepper->compute(state, t, dt); - // Expected result using analytical solution: x(t) = e^(-t) - double expected = std::exp(-dt); + // Expected result using analytical solution: x(t) = e^(-t) + double expected = std::exp(-dt); - EXPECT_NEAR(state, expected, 1e-3) << "Single step Runge-Kutta integration failed!"; + EXPECT_NEAR(state, expected, 1e-3) + << "Single step Runge-Kutta integration failed!"; } // Multiple step integration TEST_F(RungeKuttaTimeStepperTest, MultipleSteps) { - // Initial values - double state = 1.0; - double t = 0.0; - double dt = 0.1; - int steps = 10; + // Initial values + double state = 1.0; + double t = 0.0; + double dt = 0.1; + int steps = 10; - for (int i = 0; i < steps; i++) { - stepper->compute(state, t, dt); - } + for (int i = 0; i < steps; i++) { + stepper->compute(state, t, dt); + } - // Expected result: x(t) = e^(-t) - double expected = std::exp(-1.0); + // Expected result: x(t) = e^(-t) + double expected = std::exp(-1.0); - EXPECT_NEAR(state, expected, 1e-2) << "Multiple step Runge-Kutta integration failed!"; + EXPECT_NEAR(state, expected, 1e-2) + << "Multiple step Runge-Kutta integration failed!"; } // Convergence to Analytical Solution TEST_F(RungeKuttaTimeStepperTest, Convergence) { - // Initial values - double state = 1.0; - double t = 0.0; - double dt = 0.01; - int steps = 100; + // Initial values + double state = 1.0; + double t = 0.0; + double dt = 0.01; + int steps = 100; - for (int i = 0; i < steps; i++) { - stepper->compute(state, t, dt); - } + for (int i = 0; i < steps; i++) { + stepper->compute(state, t, dt); + } - double expected = std::exp(-1.0); + double expected = std::exp(-1.0); - EXPECT_NEAR(state, expected, 1e-3) << "Runge-Kutta integration does not converge!"; + EXPECT_NEAR(state, expected, 1e-3) + << "Runge-Kutta integration does not converge!"; } // // Integrating Sine function // TEST_F(RungeKuttaTimeStepperTest, SineFunction) { -// auto sine_stepper = std::make_shared>(sine_derivative); +// auto sine_stepper = +// std::make_shared>(sine_derivative); // // Initial values // double state = 0.0; @@ -89,60 +94,59 @@ TEST_F(RungeKuttaTimeStepperTest, Convergence) { // // Expected integral of sin(t) over [0, 1] is 1 - cos(1) // double expected = 1 - std::cos(1); -// EXPECT_NEAR(state, expected, 1e-2) << "Runge-Kutta integration for sine function failed!"; +// EXPECT_NEAR(state, expected, 1e-2) << "Runge-Kutta integration for sine +// function failed!"; // } // Handling small steps sizes TEST_F(RungeKuttaTimeStepperTest, SmallStepSize) { - // Initial values - double state = 1.0; - double t = 0.0; - double dt = 1e-5; - int steps = 100000; + // Initial values + double state = 1.0; + double t = 0.0; + double dt = 1e-5; + int steps = 100000; - for (int i = 0; i < steps; i++) { - stepper->compute(state, t, dt); - } + for (int i = 0; i < steps; i++) { + stepper->compute(state, t, dt); + } - double expected = std::exp(-1.0); + double expected = std::exp(-1.0); - EXPECT_NEAR(state, expected, 1e-3) << "Runge-Kutta fails with small step sizes!"; + EXPECT_NEAR(state, expected, 1e-3) + << "Runge-Kutta fails with small step sizes!"; } // Handling large steps sizes TEST_F(RungeKuttaTimeStepperTest, LargeStepSize) { - // Initial values - double state = 1.0; - double t = 0.0; - double dt = 1.0; + // Initial values + double state = 1.0; + double t = 0.0; + double dt = 1.0; - stepper->compute(state, t, dt); + stepper->compute(state, t, dt); - double expected = std::exp(-1.0); + double expected = std::exp(-1.0); - EXPECT_NEAR(state, expected, 1e-1) << "Runge-Kutta is unstable with large step sizes!"; + EXPECT_NEAR(state, expected, 1e-1) + << "Runge-Kutta is unstable with large step sizes!"; } // Constant derivative (dx/dt = 0) TEST_F(RungeKuttaTimeStepperTest, ConstantFunction) { - auto constant_stepper = std::make_shared>( - [](const TestState &state, double t) { - return 0.0; - } - ); - - // Initial values - double state = 5.0; - double t = 0.0; - double dt = 0.1; - int steps = 10; - - for (int i = 0; i < steps; i++) { - constant_stepper->compute(state, t, dt); - } - - EXPECT_NEAR(state, 5.0, 1e-6) << "Runge-Kutta should not change a constant function!"; + auto constant_stepper = + std::make_shared>( + [](const TestState &state, double t) { return 0.0; }); + + // Initial values + double state = 5.0; + double t = 0.0; + double dt = 0.1; + int steps = 10; + + for (int i = 0; i < steps; i++) { + constant_stepper->compute(state, t, dt); + } + + EXPECT_NEAR(state, 5.0, 1e-6) + << "Runge-Kutta should not change a constant function!"; } - - -