llvmbot wrote:
<!--LLVM PR SUMMARY COMMENT--> @llvm/pr-subscribers-mlir Author: None (Abhinav271828) <details> <summary>Changes</summary> We implement a function to compute the generating function for a signed unimodular cone with a parametric vertex. These functions, summed over tangential cones, provide the generating function corresponding to a parametric polyhedron. --- Full diff: https://github.com/llvm/llvm-project/pull/77199.diff 8 Files Affected: - (added) mlir/include/mlir/Analysis/Presburger/Barvinok.h (+90) - (renamed) mlir/include/mlir/Analysis/Presburger/GeneratingFunction.h (+4-1) - (added) mlir/lib/Analysis/Presburger/Barvinok.cpp (+134) - (modified) mlir/lib/Analysis/Presburger/CMakeLists.txt (+1) - (added) mlir/unittests/Analysis/Presburger/BarvinokTest.cpp (+48) - (modified) mlir/unittests/Analysis/Presburger/CMakeLists.txt (+2) - (added) mlir/unittests/Analysis/Presburger/GeneratingFunctionTest.cpp (+39) - (modified) mlir/unittests/Analysis/Presburger/Utils.h (+35-1) ``````````diff diff --git a/mlir/include/mlir/Analysis/Presburger/Barvinok.h b/mlir/include/mlir/Analysis/Presburger/Barvinok.h new file mode 100644 index 00000000000000..93b29e2d718e59 --- /dev/null +++ b/mlir/include/mlir/Analysis/Presburger/Barvinok.h @@ -0,0 +1,90 @@ +//===- Barvinok.h - Barvinok's Algorithm ------------------------*- C++ -*-===// +// +// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions. +// See https://llvm.org/LICENSE.txt for license information. +// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception +// +//===----------------------------------------------------------------------===// +// +// Implementation of Barvinok's algorithm and related utility functions. +// Currently a work in progress. +// These include functions to manipulate cones (define a cone object, get its +// dual, and find its index). +// +// The implementation is based on: +// 1. Barvinok, Alexander, and James E. Pommersheim. "An algorithmic theory of +// lattice points in polyhedra." New perspectives in algebraic combinatorics +// 38 (1999): 91-147. +// 2. Verdoolaege, Sven, et al. "Counting integer points in parametric +// polytopes using Barvinok's rational functions." Algorithmica 48 (2007): +// 37-66. +// +//===----------------------------------------------------------------------===// + +#ifndef MLIR_ANALYSIS_PRESBURGER_BARVINOK_H +#define MLIR_ANALYSIS_PRESBURGER_BARVINOK_H + +#include "mlir/Analysis/Presburger/GeneratingFunction.h" +#include "mlir/Analysis/Presburger/IntegerRelation.h" +#include "mlir/Analysis/Presburger/Matrix.h" +#include <optional> + +namespace mlir { +namespace presburger { +namespace detail { + +/// A polyhedron in H-representation is a set of inequalities +/// in d variables with integer coefficients. +using PolyhedronH = IntegerRelation; + +/// A polyhedron in V-representation is a set of rays and points, i.e., +/// vectors, stored as rows of a matrix. +using PolyhedronV = IntMatrix; + +/// A cone in either representation is a special case of +/// a polyhedron in that representation. +using ConeH = PolyhedronH; +using ConeV = PolyhedronV; + +inline ConeH defineHRep(int numVars) { + // We don't distinguish between domain and range variables, so + // we set the number of domain variables as 0 and the number of + // range variables as the number of actual variables. + // There are no symbols (we don't work with parametric cones) and no local + // (existentially quantified) variables. + // Once the cone is defined, we use `addInequality()` to set inequalities. + return ConeH(PresburgerSpace::getSetSpace(/*numDims=*/numVars, + /*numSymbols=*/0, + /*numLocals=*/0)); +} + +/// Get the index of a cone, i.e., the volume of the parallelepiped +/// spanned by its generators, which is equal to the number of integer +/// points in its fundamental parallelepiped. +/// If the index is 1, the cone is unimodular. +/// Barvinok, A., and J. E. Pommersheim. "An algorithmic theory of lattice +/// points in polyhedra." p. 107 If it has more rays than the dimension, return +/// 0. +MPInt getIndex(ConeV cone); + +/// Given a cone in H-representation, return its dual. The dual cone is in +/// V-representation. +/// This assumes that the input is pointed at the origin; it assert-fails +/// otherwise. +ConeV getDual(ConeH cone); + +/// Given a cone in V-representation, return its dual. The dual cone is in +/// H-representation. +/// The returned cone is pointed at the origin. +ConeH getDual(ConeV cone); + +/// Compute the generating function for a unimodular cone. +/// It assert-fails if the input cone is not unimodular. +GeneratingFunction unimodularConeGeneratingFunction(ParamPoint vertex, int sign, + ConeH cone); + +} // namespace detail +} // namespace presburger +} // namespace mlir + +#endif // MLIR_ANALYSIS_PRESBURGER_BARVINOK_H diff --git a/mlir/lib/Analysis/Presburger/GeneratingFunction.h b/mlir/include/mlir/Analysis/Presburger/GeneratingFunction.h similarity index 96% rename from mlir/lib/Analysis/Presburger/GeneratingFunction.h rename to mlir/include/mlir/Analysis/Presburger/GeneratingFunction.h index f7deba921ea51e..4dd692c251563b 100644 --- a/mlir/lib/Analysis/Presburger/GeneratingFunction.h +++ b/mlir/include/mlir/Analysis/Presburger/GeneratingFunction.h @@ -19,6 +19,7 @@ namespace mlir { namespace presburger { +namespace detail { // A parametric point is a vector, each of whose elements // is an affine function of n parameters. Each row @@ -83,7 +84,8 @@ class GeneratingFunction { std::vector<std::vector<Point>> sumDenominators = denominators; sumDenominators.insert(sumDenominators.end(), gf.denominators.begin(), gf.denominators.end()); - return GeneratingFunction(0, sumSigns, sumNumerators, sumDenominators); + return GeneratingFunction(numParam, sumSigns, sumNumerators, + sumDenominators); } llvm::raw_ostream &print(llvm::raw_ostream &os) const { @@ -128,6 +130,7 @@ class GeneratingFunction { std::vector<std::vector<Point>> denominators; }; +} // namespace detail } // namespace presburger } // namespace mlir diff --git a/mlir/lib/Analysis/Presburger/Barvinok.cpp b/mlir/lib/Analysis/Presburger/Barvinok.cpp new file mode 100644 index 00000000000000..31fafbda1129bf --- /dev/null +++ b/mlir/lib/Analysis/Presburger/Barvinok.cpp @@ -0,0 +1,134 @@ +//===- Barvinok.cpp - Barvinok's Algorithm ----------------------*- C++ -*-===// +// +// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions. +// See https://llvm.org/LICENSE.txt for license information. +// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception +// +//===----------------------------------------------------------------------===// + +#include "mlir/Analysis/Presburger/Barvinok.h" +#include "mlir/Analysis/Presburger/GeneratingFunction.h" + +using namespace mlir; +using namespace presburger; +using namespace mlir::presburger::detail; + +/// Assuming that the input cone is pointed at the origin, +/// converts it to its dual in V-representation. +/// Essentially we just remove the all-zeroes constant column. +ConeV mlir::presburger::detail::getDual(ConeH cone) { + unsigned numIneq = cone.getNumInequalities(); + unsigned numVar = cone.getNumCols() - 1; + ConeV dual(numIneq, numVar, 0, 0); + // Assuming that an inequality of the form + // a1*x1 + ... + an*xn + b ≥ 0 + // is represented as a row [a1, ..., an, b] + // and that b = 0. + + for (unsigned i = 0; i < numIneq; ++i) { + assert(cone.atIneq(i, numVar) == 0 && + "H-representation of cone is not centred at the origin!"); + for (unsigned j = 0; j < numVar; ++j) { + dual.at(i, j) = cone.atIneq(i, j); + } + } + + // Now dual is of the form [ [a1, ..., an] , ... ] + // which is the V-representation of the dual. + return dual; +} + +/// Converts a cone in V-representation to the H-representation +/// of its dual, pointed at the origin (not at the original vertex). +/// Essentially adds a column consisting only of zeroes to the end. +ConeH mlir::presburger::detail::getDual(ConeV cone) { + unsigned rows = cone.getNumRows(); + unsigned columns = cone.getNumColumns(); + ConeH dual = defineHRep(columns); + // Add a new column (for constants) at the end. + // This will be initialized to zero. + cone.insertColumn(columns); + + for (unsigned i = 0; i < rows; ++i) + dual.addInequality(cone.getRow(i)); + + // Now dual is of the form [ [a1, ..., an, 0] , ... ] + // which is the H-representation of the dual. + return dual; +} + +/// Find the index of a cone in V-representation. +MPInt mlir::presburger::detail::getIndex(ConeV cone) { + if (cone.getNumRows() > cone.getNumColumns()) + return MPInt(0); + + return cone.determinant(); +} + +/// Compute the generating function for a unimodular cone. +GeneratingFunction mlir::presburger::detail::unimodularConeGeneratingFunction( + ParamPoint vertex, int sign, ConeH cone) { + // `cone` is assumed to be unimodular. + assert(getIndex(getDual(cone)) == 1 && "input cone is not unimodular!"); + + unsigned numVar = cone.getNumVars(); + unsigned numIneq = cone.getNumInequalities(); + + // Thus its ray matrix, U, is the inverse of the + // transpose of its inequality matrix, `cone`. + FracMatrix transp(numVar, numIneq); + for (unsigned i = 0; i < numVar; i++) + for (unsigned j = 0; j < numIneq; j++) + transp(j, i) = Fraction(cone.atIneq(i, j), 1); + + FracMatrix generators(numVar, numIneq); + transp.determinant(&generators); // This is the U-matrix. + + // The denominators of the generating function + // are given by the generators of the cone, i.e., + // the rows of the matrix U. + std::vector<Point> denominator(numIneq); + ArrayRef<Fraction> row; + for (unsigned i = 0; i < numVar; i++) { + row = generators.getRow(i); + denominator[i] = Point(row); + } + + // The vertex is v : [d, n+1]. + // We need to find affine functions of parameters λi(p) + // such that v = Σ λi(p)*ui. + // The λi are given by the columns of Λ = v^T @ U^{-1} = v^T @ transp. + // Then the numerator will be Σ -floor(-λi(p))*u_i. + // Thus we store the numerator as the affine function -Λ, + // since the generators are already stored in the denominator. + // Note that the outer -1 will have to be accounted for, as it is not stored. + // See end for an example. + + unsigned numColumns = vertex.getNumColumns(); + unsigned numRows = vertex.getNumRows(); + ParamPoint numerator(numColumns, numRows); + SmallVector<Fraction> ithCol(numRows); + for (unsigned i = 0; i < numColumns; i++) { + for (unsigned j = 0; j < numRows; j++) + ithCol[j] = vertex(j, i); + numerator.setRow(i, transp.preMultiplyWithRow(ithCol)); + numerator.negateRow(i); + } + + return GeneratingFunction(numColumns, SmallVector<int>(1, sign), + std::vector({numerator}), + std::vector({denominator})); + + // Suppose the vertex is given by the matrix [ 2 2 0], with 2 params + // [-1 -1/2 1] + // and the cone has H-representation [0 -1] => U-matrix [ 2 -1] + // [-1 -2] [-1 0] + // Therefore Λ will be given by [ 1 0 ] and the negation of this will be + // stored as the numerator. + // [ 1/2 -1 ] + // [ -1 -2 ] + + // Algebraically, the numerator exponent is + // [ -2 ⌊ -N - M/2 +1 ⌋ + 1 ⌊ 0 +M +2 ⌋ ] -> first COLUMN of U is [2, -1] + // [ 1 ⌊ -N - M/2 +1 ⌋ + 0 ⌊ 0 +M +2 ⌋ ] -> second COLUMN of U is [-1, 0] +} diff --git a/mlir/lib/Analysis/Presburger/CMakeLists.txt b/mlir/lib/Analysis/Presburger/CMakeLists.txt index e77e1623dae175..83d0514c9e7d17 100644 --- a/mlir/lib/Analysis/Presburger/CMakeLists.txt +++ b/mlir/lib/Analysis/Presburger/CMakeLists.txt @@ -1,4 +1,5 @@ add_mlir_library(MLIRPresburger + Barvinok.cpp IntegerRelation.cpp LinearTransform.cpp Matrix.cpp diff --git a/mlir/unittests/Analysis/Presburger/BarvinokTest.cpp b/mlir/unittests/Analysis/Presburger/BarvinokTest.cpp new file mode 100644 index 00000000000000..b88baa6c6b48a4 --- /dev/null +++ b/mlir/unittests/Analysis/Presburger/BarvinokTest.cpp @@ -0,0 +1,48 @@ +#include "mlir/Analysis/Presburger/Barvinok.h" +#include "./Utils.h" +#include <gmock/gmock.h> +#include <gtest/gtest.h> + +using namespace mlir; +using namespace presburger; +using namespace mlir::presburger::detail; + +/// The following are 3 randomly generated vectors with 4 +/// entries each and define a cone's H-representation +/// using these numbers. We check that the dual contains +/// the same numbers. +/// We do the same in the reverse case. +TEST(BarvinokTest, getDual) { + ConeH cone1 = defineHRep(4); + cone1.addInequality({1, 2, 3, 4, 0}); + cone1.addInequality({3, 4, 2, 5, 0}); + cone1.addInequality({6, 2, 6, 1, 0}); + + ConeV dual1 = getDual(cone1); + + EXPECT_EQ_INT_MATRIX( + dual1, makeIntMatrix(3, 4, {{1, 2, 3, 4}, {3, 4, 2, 5}, {6, 2, 6, 1}})); + + ConeV cone2 = makeIntMatrix(3, 4, {{3, 6, 1, 5}, {3, 1, 7, 2}, {9, 3, 2, 7}}); + + ConeH dual2 = getDual(cone2); + + ConeH expected = defineHRep(4); + expected.addInequality({3, 6, 1, 5, 0}); + expected.addInequality({3, 1, 7, 2, 0}); + expected.addInequality({9, 3, 2, 7, 0}); + + EXPECT_TRUE(dual2.isEqual(expected)); +} + +/// We randomly generate a nxn matrix to use as a cone +/// with n inequalities in n variables and check for +/// the determinant being equal to the index. +TEST(BarvinokTest, getIndex) { + ConeV cone = makeIntMatrix(3, 3, {{4, 2, 1}, {5, 2, 7}, {4, 1, 6}}); + EXPECT_EQ(getIndex(cone), cone.determinant()); + + cone = makeIntMatrix( + 4, 4, {{4, 2, 5, 1}, {4, 1, 3, 6}, {8, 2, 5, 6}, {5, 2, 5, 7}}); + EXPECT_EQ(getIndex(cone), cone.determinant()); +} diff --git a/mlir/unittests/Analysis/Presburger/CMakeLists.txt b/mlir/unittests/Analysis/Presburger/CMakeLists.txt index e37133354e53ca..c98668f63fa5dc 100644 --- a/mlir/unittests/Analysis/Presburger/CMakeLists.txt +++ b/mlir/unittests/Analysis/Presburger/CMakeLists.txt @@ -1,5 +1,7 @@ add_mlir_unittest(MLIRPresburgerTests + BarvinokTest.cpp FractionTest.cpp + GeneratingFunctionTest.cpp IntegerPolyhedronTest.cpp IntegerRelationTest.cpp LinearTransformTest.cpp diff --git a/mlir/unittests/Analysis/Presburger/GeneratingFunctionTest.cpp b/mlir/unittests/Analysis/Presburger/GeneratingFunctionTest.cpp new file mode 100644 index 00000000000000..5df1a5a0f64c05 --- /dev/null +++ b/mlir/unittests/Analysis/Presburger/GeneratingFunctionTest.cpp @@ -0,0 +1,39 @@ +//===- MatrixTest.cpp - Tests for QuasiPolynomial -------------------------===// +// +// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions. +// See https://llvm.org/LICENSE.txt for license information. +// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception +// +//===----------------------------------------------------------------------===// + +#include "mlir/Analysis/Presburger/GeneratingFunction.h" +#include "./Utils.h" +#include <gmock/gmock.h> +#include <gtest/gtest.h> + +using namespace mlir; +using namespace presburger; +using namespace mlir::presburger::detail; + +TEST(GeneratingFunctionTest, sum) { + GeneratingFunction gf1(2, {1, -1}, + {makeFracMatrix(2, 3, {{1, 2, 5}, {7, 2, 6}}), + makeFracMatrix(2, 3, {{5, 2, 5}, {3, 7, 2}})}, + {{{3, 6}, {7, 2}}, {{2, 8}, {6, 3}}}); + GeneratingFunction gf2(2, {1, 1}, + {makeFracMatrix(2, 3, {{6, 2, 1}, {4, 2, 6}}), + makeFracMatrix(2, 3, {{3, 2, 6}, {9, 2, 5}})}, + {{{3, 7}, {5, 1}}, {{5, 2}, {6, 2}}}); + + GeneratingFunction sum = gf1 + gf2; + EXPECT_EQ_REPR_GENERATINGFUNCTION( + sum, GeneratingFunction(2, {1, -1, 1, 1}, + {makeFracMatrix(2, 3, {{1, 2, 5}, {7, 2, 6}}), + makeFracMatrix(2, 3, {{5, 2, 5}, {3, 7, 2}}), + makeFracMatrix(2, 3, {{6, 2, 1}, {4, 2, 6}}), + makeFracMatrix(2, 3, {{3, 2, 6}, {9, 2, 5}})}, + {{{3, 6}, {7, 2}}, + {{2, 8}, {6, 3}}, + {{3, 7}, {5, 1}}, + {{5, 2}, {6, 2}}})); +} diff --git a/mlir/unittests/Analysis/Presburger/Utils.h b/mlir/unittests/Analysis/Presburger/Utils.h index 2a9966c7ce2ea5..6b00898a7e2749 100644 --- a/mlir/unittests/Analysis/Presburger/Utils.h +++ b/mlir/unittests/Analysis/Presburger/Utils.h @@ -13,6 +13,7 @@ #ifndef MLIR_UNITTESTS_ANALYSIS_PRESBURGER_UTILS_H #define MLIR_UNITTESTS_ANALYSIS_PRESBURGER_UTILS_H +#include "mlir/Analysis/Presburger/GeneratingFunction.h" #include "mlir/Analysis/Presburger/IntegerRelation.h" #include "mlir/Analysis/Presburger/Matrix.h" #include "mlir/Analysis/Presburger/PWMAFunction.h" @@ -72,9 +73,42 @@ inline void EXPECT_EQ_FRAC_MATRIX(FracMatrix a, FracMatrix b) { EXPECT_EQ(a(row, col), b(row, col)); } +// Check the coefficients (in order) of two generating functions. +// Note that this is not a true equality check. +inline void EXPECT_EQ_REPR_GENERATINGFUNCTION(detail::GeneratingFunction a, + detail::GeneratingFunction b) { + EXPECT_EQ(a.getNumParams(), b.getNumParams()); + + SmallVector<int> aSigns = a.getSigns(); + SmallVector<int> bSigns = b.getSigns(); + EXPECT_EQ(aSigns.size(), bSigns.size()); + for (unsigned i = 0, e = aSigns.size(); i < e; i++) + EXPECT_EQ(aSigns[i], bSigns[i]); + + std::vector<detail::ParamPoint> aNums = a.getNumerators(); + std::vector<detail::ParamPoint> bNums = b.getNumerators(); + EXPECT_EQ(aNums.size(), bNums.size()); + for (unsigned i = 0, e = aNums.size(); i < e; i++) + EXPECT_EQ_FRAC_MATRIX(aNums[i], bNums[i]); + + std::vector<std::vector<detail::Point>> aDens = a.getDenominators(); + std::vector<std::vector<detail::Point>> bDens = b.getDenominators(); + EXPECT_EQ(aDens.size(), bDens.size()); + for (unsigned i = 0, e = aDens.size(); i < e; i++) { + EXPECT_EQ(aDens[i].size(), bDens[i].size()); + for (unsigned j = 0, f = aDens[i].size(); j < f; j++) { + EXPECT_EQ(aDens[i][j].size(), bDens[i][j].size()); + for (unsigned k = 0, g = aDens[i][j].size(); k < g; k++) { + EXPECT_EQ(aDens[i][j][k], bDens[i][j][k]); + } + } + } +} + // Check the coefficients (in order) of two quasipolynomials. // Note that this is not a true equality check. -inline void EXPECT_EQ_REPR_QUASIPOLYNOMIAL(QuasiPolynomial a, QuasiPolynomial b) { +inline void EXPECT_EQ_REPR_QUASIPOLYNOMIAL(QuasiPolynomial a, + QuasiPolynomial b) { EXPECT_EQ(a.getNumInputs(), b.getNumInputs()); SmallVector<Fraction> aCoeffs = a.getCoefficients(), `````````` </details> https://github.com/llvm/llvm-project/pull/77199 _______________________________________________ cfe-commits mailing list cfe-commits@lists.llvm.org https://lists.llvm.org/cgi-bin/mailman/listinfo/cfe-commits
