diff --git a/deps/fbitset b/deps/fbitset index 765b2f7..0a21596 160000 --- a/deps/fbitset +++ b/deps/fbitset @@ -1 +1 @@ -Subproject commit 765b2f740ca71e732d047aa65e067c1a9018feee +Subproject commit 0a21596ccbcb165237b0070a3b93584839ce1a21 diff --git a/include/libparenth.hpp b/include/libparenth.hpp index f79bb2a..87eaaa9 100644 --- a/include/libparenth.hpp +++ b/include/libparenth.hpp @@ -97,6 +97,7 @@ class Parenther { , n_sums_{ n_sums } , factors_with_{} , dims_on_{} + , if_classical_{ true } { assert(std::is_sorted(dims_.cbegin(), dims_.cbegin() + n_sums_)); @@ -116,6 +117,13 @@ class Parenther { } assert(factor_idx == n_factors); + + for (Size i = 0; i < n_sums_; ++i) { + if (factors_with_[i].count() != 2) { + if_classical_ = false; + break; + } + } } // Except the above constructor, normally we would put the basic @@ -750,7 +758,20 @@ class Parenther { } // Unnatural partition. - if (if_for_opt && chunks.size() > 2 + // + // This skipping is an acceleration for classical contractions, + // where every summation is involved by exactly two factors. There + // a partition shattering the factors into more than two chunks + // that are all already memoized is also reachable from a smaller + // set of broken summations, so nothing is lost by skipping it. + // When a summation is involved by more than two factors that no + // longer holds: the skipping then loses the optimal + // parenthesization, and skipping every candidate would leave the + // subproblem without any evaluation at all, making the read of + // `evals.front()` below undefined. Hence the restriction to + // classical problems, and the guard keeping one evaluation. + if (if_for_opt && if_classical_ && !evals.empty() + && chunks.size() > 2 && std::all_of( chunks.cbegin(), chunks.cend(), [&mem](const Subset& i) { return mem.count(i.factors) != 0; @@ -872,6 +893,14 @@ class Parenther { */ std::vector dims_on_; + + /** If every summation is involved by exactly two factors. + * + * This is the classical tensor contraction problem. Some of the + * accelerations of the search are only valid under this assumption. + */ + + bool if_classical_; }; } diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index dd78f0b..5fb6814 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -12,6 +12,7 @@ FetchContent_MakeAvailable(Catch2) # The main test driver. add_executable(testmain matrixchain.cpp + nonclassical.cpp ) target_link_libraries(testmain PRIVATE Catch2::Catch2WithMain) diff --git a/test/nonclassical.cpp b/test/nonclassical.cpp new file mode 100644 index 0000000..d45d3ab --- /dev/null +++ b/test/nonclassical.cpp @@ -0,0 +1,99 @@ +/** Tests on contractions where a summation is involved by more than two + * factors. + * + * A classical tensor contraction has every summation involved by exactly two + * factors. The search has an acceleration that skips a partition shattering + * the factors into more than two already-memoized chunks, which is only valid + * under that assumption. These tests cover problems outside it. + */ + +#include +#include + +#include + +// Disable Catch2's range detection for fbitset by providing stream insertion +// operator. +namespace fbitset { +template +inline std::ostream& operator<<(std::ostream& os, const Fbitset& fs) +{ + os << "Fbitset<" << N << ">{count=" << fs.count() << "}"; + return os; +} +} + +#include + +using namespace libparenth; + +using Dim = size_t; +using P = Parenther; + +/** Reads the cost of the optimal evaluation of the whole problem. + */ + +static Dim top_cost(const P::Mem& res, size_t n_factors) +{ + P::Factor_subset all(n_factors, true); + auto it = res.find(all); + REQUIRE(it != res.end()); + REQUIRE(!it->second.evals.empty()); + return it->second.evals.front().cost; +} + +TEST_CASE("A summation over four factors can be parenthesized") +{ + // s = sum_i x[i] y[i] z[i] w[i], with no external index at all. Every + // bipartition of the four factors breaks the single summation, so the + // partitions shatter into four chunks. Before the acceleration was + // restricted to classical problems, every candidate was skipped here and + // the subproblem was left without any evaluation. + std::vector dims = { 100 }; + std::vector> factors = { { 0 }, { 0 }, { 0 }, { 0 } }; + + P parenther( + dims.cbegin(), dims.cend(), 1, factors.cbegin(), factors.cend()); + + SECTION("The greedy strategy terminates with an evaluation") + { + auto res = parenther.opt(Mode::GREEDY, false); + CHECK(top_cost(res, 4) > 0); + } + + SECTION("The optimal strategy terminates with an evaluation") + { + auto res = parenther.opt(Mode::NORMAL, false); + CHECK(top_cost(res, 4) > 0); + } + + SECTION("All strategies agree on the optimal cost") + { + P p_normal( + dims.cbegin(), dims.cend(), 1, factors.cbegin(), factors.cend()); + P p_exhaust( + dims.cbegin(), dims.cend(), 1, factors.cbegin(), factors.cend()); + + auto normal = top_cost(p_normal.opt(Mode::NORMAL, false), 4); + auto exhaust = top_cost(p_exhaust.opt(Mode::EXHAUST, false), 4); + CHECK(normal == exhaust); + } +} + +TEST_CASE("A shared summation does not hide the optimal parenthesization") +{ + // t[e] = sum_{i, j} a[i, j, e] b[i, j] c[i, j] d[j]. Both summations are + // involved by more than two factors, so the optimal parenthesization is + // only reachable through a partition that the acceleration used to skip. + std::vector dims = { 10, 20, 30 }; + std::vector> factors + = { { 0, 1, 2 }, { 0, 1 }, { 0, 1 }, { 1 } }; + + P p_normal(dims.cbegin(), dims.cend(), 2, factors.cbegin(), factors.cend()); + P p_exhaust( + dims.cbegin(), dims.cend(), 2, factors.cbegin(), factors.cend()); + + auto normal = top_cost(p_normal.opt(Mode::NORMAL, false), 4); + auto exhaust = top_cost(p_exhaust.opt(Mode::EXHAUST, false), 4); + CHECK(normal == exhaust); +}