From 585e9ad7f72dbccd64cc94742c14aa2949124a4a Mon Sep 17 00:00:00 2001 From: Jean-Paul Pelteret Date: Mon, 27 Dec 2021 10:13:40 +0100 Subject: [PATCH] Add SD::BatchOptimizer::copy_from() function --- .../changes/minor/20211127Jean-PaulPelteret | 4 + .../differentiation/sd/symengine_optimizer.h | 13 ++ .../differentiation/sd/symengine_optimizer.cc | 17 ++ tests/symengine/batch_optimizer_08_1.cc | 42 +++++ tests/symengine/batch_optimizer_08_1.output | 74 +++++++++ tests/symengine/batch_optimizer_08_2.cc | 41 +++++ tests/symengine/batch_optimizer_08_2.output | 74 +++++++++ tests/symengine/batch_optimizer_08_3.cc | 42 +++++ ...mizer_08_3.with_symengine_with_llvm.output | 38 +++++ .../sd_common_tests/batch_optimizer_08.h | 153 ++++++++++++++++++ 10 files changed, 498 insertions(+) create mode 100644 doc/news/changes/minor/20211127Jean-PaulPelteret create mode 100644 tests/symengine/batch_optimizer_08_1.cc create mode 100644 tests/symengine/batch_optimizer_08_1.output create mode 100644 tests/symengine/batch_optimizer_08_2.cc create mode 100644 tests/symengine/batch_optimizer_08_2.output create mode 100644 tests/symengine/batch_optimizer_08_3.cc create mode 100644 tests/symengine/batch_optimizer_08_3.with_symengine_with_llvm.output create mode 100644 tests/symengine/sd_common_tests/batch_optimizer_08.h diff --git a/doc/news/changes/minor/20211127Jean-PaulPelteret b/doc/news/changes/minor/20211127Jean-PaulPelteret new file mode 100644 index 0000000000..0870a8d973 --- /dev/null +++ b/doc/news/changes/minor/20211127Jean-PaulPelteret @@ -0,0 +1,4 @@ +New: The SD::BatchOptimizer::copy_from() function can be used to duplicate the +data from one batch optimizer instance to another. +
+(Jean-Paul Pelteret, 2021/11/27) diff --git a/include/deal.II/differentiation/sd/symengine_optimizer.h b/include/deal.II/differentiation/sd/symengine_optimizer.h index 89e203c647..0eef522569 100644 --- a/include/deal.II/differentiation/sd/symengine_optimizer.h +++ b/include/deal.II/differentiation/sd/symengine_optimizer.h @@ -1498,6 +1498,19 @@ namespace Differentiation */ ~BatchOptimizer() = default; + /** + * Duplicate the data stored in an @p other BatchOptimizer instance. + * + * @note The optimized data and results from previous substitutions + * executed by the @p other optimizer instance are not copied over. + * It is therefore necessary to call optimize() before it is possible to + * substitute() values and evaluate() data. One may, however, still + * extract() values using @p this optimizer instance if those results are + * stored elsewhere. + */ + void + copy_from(const BatchOptimizer &other); + /** * Print some information on state of the internal data * structures stored in the class. diff --git a/source/differentiation/sd/symengine_optimizer.cc b/source/differentiation/sd/symengine_optimizer.cc index db0c0840bb..cc2a84243f 100644 --- a/source/differentiation/sd/symengine_optimizer.cc +++ b/source/differentiation/sd/symengine_optimizer.cc @@ -67,6 +67,23 @@ namespace Differentiation + template + void + BatchOptimizer::copy_from( + const BatchOptimizer &other) + { + method = other.method; + flags = other.flags; + independent_variables_symbols = other.independent_variables_symbols; + dependent_variables_functions = other.dependent_variables_functions; + dependent_variables_output.clear(); + map_dep_expr_vec_entry = other.map_dep_expr_vec_entry; + ready_for_value_extraction = false; + has_been_serialized = false; + } + + + template void BatchOptimizer::set_optimization_method( diff --git a/tests/symengine/batch_optimizer_08_1.cc b/tests/symengine/batch_optimizer_08_1.cc new file mode 100644 index 0000000000..09dbae317c --- /dev/null +++ b/tests/symengine/batch_optimizer_08_1.cc @@ -0,0 +1,42 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +// Check that the BatchOptimizer::copy_from() function works as expected. +// +// Here we use only dictionary substitution, and invoke no symbolic +// optimizations. + +#include "../tests.h" + +#include "sd_common_tests/batch_optimizer_08.h" + +int +main() +{ + initlog(); + deallog << std::setprecision(10); + + // Show the difference between a SymEngine "value" and + // an evaluated, floating point number + // deallog << std::setprecision(3); + + const enum SD::OptimizerType opt_method = SD::OptimizerType::dictionary; + const enum SD::OptimizationFlags opt_flags = + SD::OptimizationFlags::optimize_all; + + const int dim = 2; + run_tests(); +} diff --git a/tests/symengine/batch_optimizer_08_1.output b/tests/symengine/batch_optimizer_08_1.output new file mode 100644 index 0000000000..7f0578c4fc --- /dev/null +++ b/tests/symengine/batch_optimizer_08_1.output @@ -0,0 +1,74 @@ + +DEAL:Copy:Float::f: x*y +DEAL:Copy:Float::g: x/y +DEAL:Copy:Float::f(x,y): 20.00000000 +DEAL:Copy:Float::g(x,y): 5.000000000 +DEAL:Copy:Float:: +DEAL:Copy:Float::Before copy +DEAL:Copy:Float::f: x1*y1 +DEAL:Copy:Float::g: x1/y1 +DEAL:Copy:Float::f(x,y): 60.00000000 +DEAL:Copy:Float::g(x,y): 6.666666508 +DEAL:Copy:Float:: +DEAL:Copy:Float::After copy +DEAL:Copy:Float::f: x*y +DEAL:Copy:Float::g: x/y +DEAL:Copy:Float::f(x,y): 20.00000000 +DEAL:Copy:Float::g(x,y): 5.000000000 +DEAL:Copy:Float:: +DEAL:Copy:Float::OK +DEAL:Copy:Double::f: x*y +DEAL:Copy:Double::g: x/y +DEAL:Copy:Double::f(x,y): 20.00000000 +DEAL:Copy:Double::g(x,y): 5.000000000 +DEAL:Copy:Double:: +DEAL:Copy:Double::Before copy +DEAL:Copy:Double::f: x1*y1 +DEAL:Copy:Double::g: x1/y1 +DEAL:Copy:Double::f(x,y): 60.00000000 +DEAL:Copy:Double::g(x,y): 6.666666667 +DEAL:Copy:Double:: +DEAL:Copy:Double::After copy +DEAL:Copy:Double::f: x*y +DEAL:Copy:Double::g: x/y +DEAL:Copy:Double::f(x,y): 20.00000000 +DEAL:Copy:Double::g(x,y): 5.000000000 +DEAL:Copy:Double:: +DEAL:Copy:Double::OK +DEAL:Copy:Complex float::f: x*y +DEAL:Copy:Complex float::g: x/y +DEAL:Copy:Complex float::f(x,y): (20.00000000,0.000000000) +DEAL:Copy:Complex float::g(x,y): (5.000000000,0.000000000) +DEAL:Copy:Complex float:: +DEAL:Copy:Complex float::Before copy +DEAL:Copy:Complex float::f: x1*y1 +DEAL:Copy:Complex float::g: x1/y1 +DEAL:Copy:Complex float::f(x,y): (60.00000000,0.000000000) +DEAL:Copy:Complex float::g(x,y): (6.666666508,0.000000000) +DEAL:Copy:Complex float:: +DEAL:Copy:Complex float::After copy +DEAL:Copy:Complex float::f: x*y +DEAL:Copy:Complex float::g: x/y +DEAL:Copy:Complex float::f(x,y): (20.00000000,0.000000000) +DEAL:Copy:Complex float::g(x,y): (5.000000000,0.000000000) +DEAL:Copy:Complex float:: +DEAL:Copy:Complex float::OK +DEAL:Copy:Complex double::f: x*y +DEAL:Copy:Complex double::g: x/y +DEAL:Copy:Complex double::f(x,y): (20.00000000,0.000000000) +DEAL:Copy:Complex double::g(x,y): (5.000000000,0.000000000) +DEAL:Copy:Complex double:: +DEAL:Copy:Complex double::Before copy +DEAL:Copy:Complex double::f: x1*y1 +DEAL:Copy:Complex double::g: x1/y1 +DEAL:Copy:Complex double::f(x,y): (60.00000000,0.000000000) +DEAL:Copy:Complex double::g(x,y): (6.666666667,0.000000000) +DEAL:Copy:Complex double:: +DEAL:Copy:Complex double::After copy +DEAL:Copy:Complex double::f: x*y +DEAL:Copy:Complex double::g: x/y +DEAL:Copy:Complex double::f(x,y): (20.00000000,0.000000000) +DEAL:Copy:Complex double::g(x,y): (5.000000000,0.000000000) +DEAL:Copy:Complex double:: +DEAL:Copy:Complex double::OK +DEAL::OK diff --git a/tests/symengine/batch_optimizer_08_2.cc b/tests/symengine/batch_optimizer_08_2.cc new file mode 100644 index 0000000000..3fd2bc2ddc --- /dev/null +++ b/tests/symengine/batch_optimizer_08_2.cc @@ -0,0 +1,41 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +// Check that the BatchOptimizer::copy_from() function works as expected. +// +// Here we use only lambda substitution, and invoke no symbolic optimizations. + +#include "../tests.h" + +#include "sd_common_tests/batch_optimizer_08.h" + +int +main() +{ + initlog(); + deallog << std::setprecision(10); + + // Show the difference between a SymEngine "value" and + // an evaluated, floating point number + // deallog << std::setprecision(3); + + const enum SD::OptimizerType opt_method = SD::OptimizerType::lambda; + const enum SD::OptimizationFlags opt_flags = + SD::OptimizationFlags::optimize_default; + + const int dim = 2; + run_tests(); +} diff --git a/tests/symengine/batch_optimizer_08_2.output b/tests/symengine/batch_optimizer_08_2.output new file mode 100644 index 0000000000..7f0578c4fc --- /dev/null +++ b/tests/symengine/batch_optimizer_08_2.output @@ -0,0 +1,74 @@ + +DEAL:Copy:Float::f: x*y +DEAL:Copy:Float::g: x/y +DEAL:Copy:Float::f(x,y): 20.00000000 +DEAL:Copy:Float::g(x,y): 5.000000000 +DEAL:Copy:Float:: +DEAL:Copy:Float::Before copy +DEAL:Copy:Float::f: x1*y1 +DEAL:Copy:Float::g: x1/y1 +DEAL:Copy:Float::f(x,y): 60.00000000 +DEAL:Copy:Float::g(x,y): 6.666666508 +DEAL:Copy:Float:: +DEAL:Copy:Float::After copy +DEAL:Copy:Float::f: x*y +DEAL:Copy:Float::g: x/y +DEAL:Copy:Float::f(x,y): 20.00000000 +DEAL:Copy:Float::g(x,y): 5.000000000 +DEAL:Copy:Float:: +DEAL:Copy:Float::OK +DEAL:Copy:Double::f: x*y +DEAL:Copy:Double::g: x/y +DEAL:Copy:Double::f(x,y): 20.00000000 +DEAL:Copy:Double::g(x,y): 5.000000000 +DEAL:Copy:Double:: +DEAL:Copy:Double::Before copy +DEAL:Copy:Double::f: x1*y1 +DEAL:Copy:Double::g: x1/y1 +DEAL:Copy:Double::f(x,y): 60.00000000 +DEAL:Copy:Double::g(x,y): 6.666666667 +DEAL:Copy:Double:: +DEAL:Copy:Double::After copy +DEAL:Copy:Double::f: x*y +DEAL:Copy:Double::g: x/y +DEAL:Copy:Double::f(x,y): 20.00000000 +DEAL:Copy:Double::g(x,y): 5.000000000 +DEAL:Copy:Double:: +DEAL:Copy:Double::OK +DEAL:Copy:Complex float::f: x*y +DEAL:Copy:Complex float::g: x/y +DEAL:Copy:Complex float::f(x,y): (20.00000000,0.000000000) +DEAL:Copy:Complex float::g(x,y): (5.000000000,0.000000000) +DEAL:Copy:Complex float:: +DEAL:Copy:Complex float::Before copy +DEAL:Copy:Complex float::f: x1*y1 +DEAL:Copy:Complex float::g: x1/y1 +DEAL:Copy:Complex float::f(x,y): (60.00000000,0.000000000) +DEAL:Copy:Complex float::g(x,y): (6.666666508,0.000000000) +DEAL:Copy:Complex float:: +DEAL:Copy:Complex float::After copy +DEAL:Copy:Complex float::f: x*y +DEAL:Copy:Complex float::g: x/y +DEAL:Copy:Complex float::f(x,y): (20.00000000,0.000000000) +DEAL:Copy:Complex float::g(x,y): (5.000000000,0.000000000) +DEAL:Copy:Complex float:: +DEAL:Copy:Complex float::OK +DEAL:Copy:Complex double::f: x*y +DEAL:Copy:Complex double::g: x/y +DEAL:Copy:Complex double::f(x,y): (20.00000000,0.000000000) +DEAL:Copy:Complex double::g(x,y): (5.000000000,0.000000000) +DEAL:Copy:Complex double:: +DEAL:Copy:Complex double::Before copy +DEAL:Copy:Complex double::f: x1*y1 +DEAL:Copy:Complex double::g: x1/y1 +DEAL:Copy:Complex double::f(x,y): (60.00000000,0.000000000) +DEAL:Copy:Complex double::g(x,y): (6.666666667,0.000000000) +DEAL:Copy:Complex double:: +DEAL:Copy:Complex double::After copy +DEAL:Copy:Complex double::f: x*y +DEAL:Copy:Complex double::g: x/y +DEAL:Copy:Complex double::f(x,y): (20.00000000,0.000000000) +DEAL:Copy:Complex double::g(x,y): (5.000000000,0.000000000) +DEAL:Copy:Complex double:: +DEAL:Copy:Complex double::OK +DEAL::OK diff --git a/tests/symengine/batch_optimizer_08_3.cc b/tests/symengine/batch_optimizer_08_3.cc new file mode 100644 index 0000000000..6fe8dc844f --- /dev/null +++ b/tests/symengine/batch_optimizer_08_3.cc @@ -0,0 +1,42 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +// Check that the BatchOptimizer::copy_from() function works as expected. +// +// Here we invoke the LLVM optimizer before symbolic evaluation takes place, +// and we invoke no additional symbolic optimizations as well. + +#include "../tests.h" + +#include "sd_common_tests/batch_optimizer_08.h" + +int +main() +{ + initlog(); + deallog << std::setprecision(10); + + // Show the difference between a SymEngine "value" and + // an evaluated, floating point number + // deallog << std::setprecision(3); + + const enum SD::OptimizerType opt_method = SD::OptimizerType::llvm; + const enum SD::OptimizationFlags opt_flags = + SD::OptimizationFlags::optimize_default; + + const int dim = 2; + run_tests(); +} diff --git a/tests/symengine/batch_optimizer_08_3.with_symengine_with_llvm.output b/tests/symengine/batch_optimizer_08_3.with_symengine_with_llvm.output new file mode 100644 index 0000000000..83eabc376f --- /dev/null +++ b/tests/symengine/batch_optimizer_08_3.with_symengine_with_llvm.output @@ -0,0 +1,38 @@ + +DEAL:Copy:Float::f: x*y +DEAL:Copy:Float::g: x/y +DEAL:Copy:Float::f(x,y): 20.00000000 +DEAL:Copy:Float::g(x,y): 5.000000000 +DEAL:Copy:Float:: +DEAL:Copy:Float::Before copy +DEAL:Copy:Float::f: x1*y1 +DEAL:Copy:Float::g: x1/y1 +DEAL:Copy:Float::f(x,y): 60.00000000 +DEAL:Copy:Float::g(x,y): 6.666666508 +DEAL:Copy:Float:: +DEAL:Copy:Float::After copy +DEAL:Copy:Float::f: x*y +DEAL:Copy:Float::g: x/y +DEAL:Copy:Float::f(x,y): 20.00000000 +DEAL:Copy:Float::g(x,y): 5.000000000 +DEAL:Copy:Float:: +DEAL:Copy:Float::OK +DEAL:Copy:Double::f: x*y +DEAL:Copy:Double::g: x/y +DEAL:Copy:Double::f(x,y): 20.00000000 +DEAL:Copy:Double::g(x,y): 5.000000000 +DEAL:Copy:Double:: +DEAL:Copy:Double::Before copy +DEAL:Copy:Double::f: x1*y1 +DEAL:Copy:Double::g: x1/y1 +DEAL:Copy:Double::f(x,y): 60.00000000 +DEAL:Copy:Double::g(x,y): 6.666666667 +DEAL:Copy:Double:: +DEAL:Copy:Double::After copy +DEAL:Copy:Double::f: x*y +DEAL:Copy:Double::g: x/y +DEAL:Copy:Double::f(x,y): 20.00000000 +DEAL:Copy:Double::g(x,y): 5.000000000 +DEAL:Copy:Double:: +DEAL:Copy:Double::OK +DEAL::OK diff --git a/tests/symengine/sd_common_tests/batch_optimizer_08.h b/tests/symengine/sd_common_tests/batch_optimizer_08.h new file mode 100644 index 0000000000..fa0f7df173 --- /dev/null +++ b/tests/symengine/sd_common_tests/batch_optimizer_08.h @@ -0,0 +1,153 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 - 2021 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +// Check that the BatchOptimizer::copy_from() function works as expected. + +#include + +#include "../../tests.h" + +using namespace dealii; +namespace SD = Differentiation::SD; + +template +void +test_copy() +{ + namespace SD = Differentiation::SD; + using SD_number_t = SD::Expression; + + // Define + const SD_number_t x("x"), y("y"); + const SD_number_t f = x * y; + const SD_number_t g = x / y; + + // Substitution map + SD::types::substitution_map sub_vals = + SD::make_substitution_map(std::make_pair(x, NumberType(10.0)), + std::make_pair(y, NumberType(2.0))); + + // Optimize + SD::BatchOptimizer optimizer(opt_method, opt_flags); + optimizer.register_symbols(sub_vals); // Independent symbols + optimizer.register_functions(f, g); + optimizer.optimize(); + + // Substitute + optimizer.substitute(sub_vals); + + // Evaluate + const NumberType val_f = optimizer.evaluate(f); + const NumberType val_g = optimizer.evaluate(g); + deallog << "f: " << f << std::endl + << "g: " << g << std::endl + << "f(x,y): " << val_f << std::endl + << "g(x,y): " << val_g << std::endl + << std::endl; + + // Define something new + { + deallog << "Before copy" << std::endl; + + const SD_number_t x1("x1"), y1("y1"); + const SD_number_t f1 = x1 * y1; + const SD_number_t g1 = x1 / y1; + + SD::add_to_substitution_map(sub_vals, + std::make_pair(x1, NumberType(20.0)), + std::make_pair(y1, NumberType(3.0))); + + // Optimizer + SD::BatchOptimizer new_optimizer; + new_optimizer.register_symbols(sub_vals); // Independent symbols + new_optimizer.register_functions(f1, g1); + new_optimizer.optimize(); + + // Substitute + new_optimizer.substitute(sub_vals); + + // Evaluate + const NumberType val_f1 = new_optimizer.evaluate(f1); + const NumberType val_g1 = new_optimizer.evaluate(g1); + deallog << "f: " << f1 << std::endl + << "g: " << g1 << std::endl + << "f(x,y): " << val_f1 << std::endl + << "g(x,y): " << val_g1 << std::endl + << std::endl; + + // Perform copy + deallog << "After copy" << std::endl; + new_optimizer.copy_from(optimizer); + Assert(new_optimizer.optimized() == false, + ExcMessage("Expected new optimizer to be unoptimized.")); + + // Extract + const NumberType val_f = new_optimizer.extract(f, optimizer.evaluate()); + const NumberType val_g = new_optimizer.extract(g, optimizer.evaluate()); + deallog << "f: " << f << std::endl + << "g: " << g << std::endl + << "f(x,y): " << val_f << std::endl + << "g(x,y): " << val_g << std::endl + << std::endl; + } + + deallog << "OK" << std::endl; +} + + +template +void +run_tests() +{ + // Show the difference between a SymEngine "value" and + // an evaluated, floating point number + // deallog << std::setprecision(3); + + deallog.push("Copy"); + { + // deallog.push("Integer"); + // test_tensor(n_runs, timer); + // deallog.pop(); + + deallog.push("Float"); + test_copy(); + deallog.pop(); + + deallog.push("Double"); + test_copy(); + deallog.pop(); + + // The LLVM optimizer does not currently support complex numbers. + if (opt_method != SD::OptimizerType::llvm) + { + deallog.push("Complex float"); + test_copy, opt_method, opt_flags>(); + deallog.pop(); + + deallog.push("Complex double"); + test_copy, opt_method, opt_flags>(); + deallog.pop(); + } + } + deallog.pop(); + + deallog << "OK" << std::endl; +} -- 2.39.5