[llvm-branch-commits] [mlir] 6ebeba8 - Support emptiness checks for unbounded FlatAffineConstraints.
Alex Zinenko via llvm-branch-commits
llvm-branch-commits at lists.llvm.org
Thu Jan 14 10:38:34 PST 2021
Author: Arjun P
Date: 2021-01-14T19:33:37+01:00
New Revision: 6ebeba88f51959d763a8f274cdfecea46d51d28c
URL: https://github.com/llvm/llvm-project/commit/6ebeba88f51959d763a8f274cdfecea46d51d28c
DIFF: https://github.com/llvm/llvm-project/commit/6ebeba88f51959d763a8f274cdfecea46d51d28c.diff
LOG: Support emptiness checks for unbounded FlatAffineConstraints.
With this, we have complete support for emptiness checks. This also paves the way for future support to check if two FlatAffineConstraints are equal.
Reviewed By: ftynse
Differential Revision: https://reviews.llvm.org/D94272
Added:
mlir/include/mlir/Analysis/LinearTransform.h
mlir/lib/Analysis/LinearTransform.cpp
mlir/unittests/Analysis/LinearTransformTest.cpp
Modified:
mlir/include/mlir/Analysis/AffineStructures.h
mlir/include/mlir/Analysis/Presburger/Fraction.h
mlir/include/mlir/Analysis/Presburger/Matrix.h
mlir/include/mlir/Analysis/Presburger/Simplex.h
mlir/lib/Analysis/AffineStructures.cpp
mlir/lib/Analysis/CMakeLists.txt
mlir/lib/Analysis/Presburger/CMakeLists.txt
mlir/lib/Analysis/Presburger/Matrix.cpp
mlir/lib/Analysis/Presburger/Simplex.cpp
mlir/unittests/Analysis/AffineStructuresTest.cpp
mlir/unittests/Analysis/CMakeLists.txt
Removed:
################################################################################
diff --git a/mlir/include/mlir/Analysis/AffineStructures.h b/mlir/include/mlir/Analysis/AffineStructures.h
index 25071db100e3..fa80db7d4b63 100644
--- a/mlir/include/mlir/Analysis/AffineStructures.h
+++ b/mlir/include/mlir/Analysis/AffineStructures.h
@@ -13,6 +13,7 @@
#ifndef MLIR_ANALYSIS_AFFINE_STRUCTURES_H
#define MLIR_ANALYSIS_AFFINE_STRUCTURES_H
+#include "mlir/Analysis/Presburger/Matrix.h"
#include "mlir/IR/AffineExpr.h"
#include "mlir/IR/OpDefinition.h"
#include "mlir/Support/LogicalResult.h"
@@ -153,6 +154,12 @@ class FlatAffineConstraints {
/// false if a solution exists or all tests were inconclusive.
bool isIntegerEmpty() const;
+ // Returns a matrix where each row is a vector along which the polytope is
+ // bounded. The span of the returned vectors is guaranteed to contain all
+ // such vectors. The returned vectors are NOT guaranteed to be linearly
+ // independent. This function should not be called on empty sets.
+ Matrix getBoundedDirections() const;
+
/// Find a sample point satisfying the constraints. This uses a branch and
/// bound algorithm with generalized basis reduction, which always works if
/// the set is bounded. This should not be called for unbounded sets.
diff --git a/mlir/include/mlir/Analysis/LinearTransform.h b/mlir/include/mlir/Analysis/LinearTransform.h
new file mode 100644
index 000000000000..0850f5a00609
--- /dev/null
+++ b/mlir/include/mlir/Analysis/LinearTransform.h
@@ -0,0 +1,48 @@
+//===- LinearTransform.h - MLIR LinearTransform Class -----------*- 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
+//
+//===----------------------------------------------------------------------===//
+//
+// Support for linear transforms and applying them to FlatAffineConstraints.
+//
+//===----------------------------------------------------------------------===//
+
+#ifndef MLIR_ANALYSIS_LINEARTRANSFORM_H
+#define MLIR_ANALYSIS_LINEARTRANSFORM_H
+
+#include "mlir/Analysis/AffineStructures.h"
+#include "mlir/Analysis/Presburger/Matrix.h"
+#include "llvm/ADT/SmallVector.h"
+
+namespace mlir {
+
+class LinearTransform {
+public:
+ explicit LinearTransform(Matrix &&oMatrix);
+ explicit LinearTransform(const Matrix &oMatrix);
+
+ // Returns a linear transform T such that MT is M in column echelon form.
+ // Also returns the number of non-zero columns in MT.
+ //
+ // Specifically, T is such that in every column the first non-zero row is
+ // strictly below that of the previous column, and all columns which have only
+ // zeros are at the end.
+ static std::pair<unsigned, LinearTransform>
+ makeTransformToColumnEchelon(Matrix m);
+
+ // Returns a FlatAffineConstraints having a constraint vector vT for every
+ // constraint vector v in fac, where T is this transform.
+ FlatAffineConstraints applyTo(const FlatAffineConstraints &fac);
+
+ // Post-multiply the given vector v with this transform, say T, returning vT.
+ SmallVector<int64_t, 8> applyTo(ArrayRef<int64_t> v);
+
+private:
+ Matrix matrix;
+};
+
+} // namespace mlir
+#endif // MLIR_ANALYSIS_LINEARTRANSFORM_H
diff --git a/mlir/include/mlir/Analysis/Presburger/Fraction.h b/mlir/include/mlir/Analysis/Presburger/Fraction.h
index 09996c486ef3..61b0915e559e 100644
--- a/mlir/include/mlir/Analysis/Presburger/Fraction.h
+++ b/mlir/include/mlir/Analysis/Presburger/Fraction.h
@@ -64,6 +64,8 @@ inline bool operator<=(Fraction x, Fraction y) { return compare(x, y) <= 0; }
inline bool operator==(Fraction x, Fraction y) { return compare(x, y) == 0; }
+inline bool operator!=(Fraction x, Fraction y) { return compare(x, y) != 0; }
+
inline bool operator>(Fraction x, Fraction y) { return compare(x, y) > 0; }
inline bool operator>=(Fraction x, Fraction y) { return compare(x, y) >= 0; }
diff --git a/mlir/include/mlir/Analysis/Presburger/Matrix.h b/mlir/include/mlir/Analysis/Presburger/Matrix.h
index 7bc29f81a834..8ed40bb9c026 100644
--- a/mlir/include/mlir/Analysis/Presburger/Matrix.h
+++ b/mlir/include/mlir/Analysis/Presburger/Matrix.h
@@ -58,6 +58,12 @@ class Matrix {
/// Add `scale` multiples of the source row to the target row.
void addToRow(unsigned sourceRow, unsigned targetRow, int64_t scale);
+ /// Add `scale` multiples of the source column to the target column.
+ void addToColumn(unsigned sourceColumn, unsigned targetColumn, int64_t scale);
+
+ /// Negate the specified column.
+ void negateColumn(unsigned column);
+
/// Resize the matrix to the specified dimensions. If a dimension is smaller,
/// the values are truncated; if it is bigger, the new values are default
/// initialized.
diff --git a/mlir/include/mlir/Analysis/Presburger/Simplex.h b/mlir/include/mlir/Analysis/Presburger/Simplex.h
index 05d241e60958..370035cbc7ba 100644
--- a/mlir/include/mlir/Analysis/Presburger/Simplex.h
+++ b/mlir/include/mlir/Analysis/Presburger/Simplex.h
@@ -17,10 +17,12 @@
#include "mlir/Analysis/AffineStructures.h"
#include "mlir/Analysis/Presburger/Fraction.h"
#include "mlir/Analysis/Presburger/Matrix.h"
+#include "mlir/IR/Location.h"
#include "mlir/Support/LogicalResult.h"
#include "llvm/ADT/ArrayRef.h"
#include "llvm/ADT/Optional.h"
#include "llvm/ADT/SmallVector.h"
+#include "llvm/Support/StringSaver.h"
#include "llvm/Support/raw_ostream.h"
namespace mlir {
@@ -84,7 +86,7 @@ class GBRSimplex;
///
/// The unknowns in row position are represented in terms of the basis unknowns.
/// If the basis unknowns are u_1, u_2, ... u_m, and a row in the tableau is
-/// d, c, a_1, a_2, ... a_m, this representats the unknown for that row as
+/// d, c, a_1, a_2, ... a_m, this represents the unknown for that row as
/// (c + a_1*u_1 + a_2*u_2 + ... + a_m*u_m)/d. In our running example, if the
/// basis is the initial basis of x, y, then the constraint 1 + 2x + 3y >= 0
/// would be represented by the row [1, 1, 2, 3].
@@ -173,20 +175,25 @@ class Simplex {
void intersectFlatAffineConstraints(const FlatAffineConstraints &fac);
/// Compute the maximum or minimum value of the given row, depending on
- /// direction. The specified row is never pivoted.
+ /// direction. The specified row is never pivoted. On return, the row may
+ /// have a negative sample value if the direction is down.
///
- /// Returns a (num, den) pair denoting the optimum, or None if no
- /// optimum exists, i.e., if the expression is unbounded in this direction.
+ /// Returns a Fraction denoting the optimum, or a null value if no optimum
+ /// exists, i.e., if the expression is unbounded in this direction.
Optional<Fraction> computeRowOptimum(Direction direction, unsigned row);
/// Compute the maximum or minimum value of the given expression, depending on
- /// direction.
+ /// direction. Should not be called when the Simplex is empty.
///
- /// Returns a (num, den) pair denoting the optimum, or a null value if no
- /// optimum exists, i.e., if the expression is unbounded in this direction.
+ /// Returns a Fraction denoting the optimum, or a null value if no optimum
+ /// exists, i.e., if the expression is unbounded in this direction.
Optional<Fraction> computeOptimum(Direction direction,
ArrayRef<int64_t> coeffs);
+ /// Returns whether the perpendicular of the specified constraint is a
+ /// is a direction along which the polytope is bounded.
+ bool isBoundedAlongConstraint(unsigned constraintIndex);
+
/// Returns whether the specified constraint has been marked as redundant.
/// Constraints are numbered from 0 starting at the first added inequality.
/// Equalities are added as a pair of inequalities and so correspond to two
@@ -299,6 +306,15 @@ class Simplex {
/// sample value, false otherwise.
LogicalResult restoreRow(Unknown &u);
+ /// Compute the maximum or minimum of the specified Unknown, depending on
+ /// direction. The specified unknown may be pivoted. If the unknown is
+ /// restricted, it will have a non-negative sample value on return.
+ /// Should not be called if the Simplex is empty.
+ ///
+ /// Returns a Fraction denoting the optimum, or a null value if no optimum
+ /// exists, i.e., if the expression is unbounded in this direction.
+ Optional<Fraction> computeOptimum(Direction direction, Unknown &u);
+
/// Mark the specified unknown redundant. This operation is added to the undo
/// log and will be undone by rollbacks. The specified unknown must be in row
/// orientation.
diff --git a/mlir/lib/Analysis/AffineStructures.cpp b/mlir/lib/Analysis/AffineStructures.cpp
index 51141e6f6184..12c90fbcfc54 100644
--- a/mlir/lib/Analysis/AffineStructures.cpp
+++ b/mlir/lib/Analysis/AffineStructures.cpp
@@ -11,6 +11,7 @@
//===----------------------------------------------------------------------===//
#include "mlir/Analysis/AffineStructures.h"
+#include "mlir/Analysis/LinearTransform.h"
#include "mlir/Analysis/Presburger/Simplex.h"
#include "mlir/Dialect/Affine/IR/AffineOps.h"
#include "mlir/Dialect/Affine/IR/AffineValueMap.h"
@@ -20,6 +21,7 @@
#include "mlir/Support/LLVM.h"
#include "mlir/Support/MathExtras.h"
#include "llvm/ADT/SmallPtrSet.h"
+#include "llvm/ADT/SmallVector.h"
#include "llvm/Support/Debug.h"
#include "llvm/Support/raw_ostream.h"
@@ -1034,21 +1036,152 @@ bool FlatAffineConstraints::isEmptyByGCDTest() const {
return false;
}
-// First, try the GCD test heuristic.
+// Returns a matrix where each row is a vector along which the polytope is
+// bounded. The span of the returned vectors is guaranteed to contain all
+// such vectors. The returned vectors are NOT guaranteed to be linearly
+// independent. This function should not be called on empty sets.
//
-// If that doesn't find the set empty, check if the set is unbounded. If it is,
-// we cannot use the GBR algorithm and we conservatively return false.
-//
-// If the set is bounded, we use the complete emptiness check for this case
-// provided by Simplex::findIntegerSample(), which gives a definitive answer.
+// It is sufficient to check the perpendiculars of the constraints, as the set
+// of perpendiculars which are bounded must span all bounded directions.
+Matrix FlatAffineConstraints::getBoundedDirections() const {
+ // Note that it is necessary to add the equalities too (which the constructor
+ // does) even though we don't need to check if they are bounded; whether an
+ // inequality is bounded or not depends on what other constraints, including
+ // equalities, are present.
+ Simplex simplex(*this);
+
+ assert(!simplex.isEmpty() && "It is not meaningful to ask whether a "
+ "direction is bounded in an empty set.");
+
+ SmallVector<unsigned, 8> boundedIneqs;
+ // The constructor adds the inequalities to the simplex first, so this
+ // processes all the inequalities.
+ for (unsigned i = 0, e = getNumInequalities(); i < e; ++i) {
+ if (simplex.isBoundedAlongConstraint(i))
+ boundedIneqs.push_back(i);
+ }
+
+ // The direction vector is given by the coefficients and does not include the
+ // constant term, so the matrix has one fewer column.
+ unsigned dirsNumCols = getNumCols() - 1;
+ Matrix dirs(boundedIneqs.size() + getNumEqualities(), dirsNumCols);
+
+ // Copy the bounded inequalities.
+ unsigned row = 0;
+ for (unsigned i : boundedIneqs) {
+ for (unsigned col = 0; col < dirsNumCols; ++col)
+ dirs(row, col) = atIneq(i, col);
+ ++row;
+ }
+
+ // Copy the equalities. All the equalities' perpendiculars are bounded.
+ for (unsigned i = 0, e = getNumEqualities(); i < e; ++i) {
+ for (unsigned col = 0; col < dirsNumCols; ++col)
+ dirs(row, col) = atEq(i, col);
+ ++row;
+ }
+
+ return dirs;
+}
+
+bool eqInvolvesSuffixDims(const FlatAffineConstraints &fac, unsigned eqIndex,
+ unsigned numDims) {
+ for (unsigned e = fac.getNumDimIds(), j = e - numDims; j < e; ++j)
+ if (fac.atEq(eqIndex, j) != 0)
+ return true;
+ return false;
+}
+bool ineqInvolvesSuffixDims(const FlatAffineConstraints &fac,
+ unsigned ineqIndex, unsigned numDims) {
+ for (unsigned e = fac.getNumDimIds(), j = e - numDims; j < e; ++j)
+ if (fac.atIneq(ineqIndex, j) != 0)
+ return true;
+ return false;
+}
+
+void removeConstraintsInvolvingSuffixDims(FlatAffineConstraints &fac,
+ unsigned unboundedDims) {
+ // We iterate backwards so that whether we remove constraint i - 1 or not, the
+ // next constraint to be tested is always i - 2.
+ for (unsigned i = fac.getNumEqualities(); i > 0; i--)
+ if (eqInvolvesSuffixDims(fac, i - 1, unboundedDims))
+ fac.removeEquality(i - 1);
+ for (unsigned i = fac.getNumInequalities(); i > 0; i--)
+ if (ineqInvolvesSuffixDims(fac, i - 1, unboundedDims))
+ fac.removeInequality(i - 1);
+}
+
+/// Let this set be S. If S is bounded then we directly call into the GBR
+/// sampling algorithm. Otherwise, there are some unbounded directions, i.e.,
+/// vectors v such that S extends to infininty along v or -v. In this case we
+/// use an algorithm described in the integer set library (isl) manual and used
+/// by the isl_set_sample function in that library. The algorithm is:
+///
+/// 1) Apply a unimodular transform T to S to obtain S*T, such that all
+/// dimensions in which S*T is bounded lie in the linear span of a prefix of the
+/// dimensions.
+///
+/// 2) Construct a set transformedSet by removing all constraints that involve
+/// the unbounded dimensions and also deleting the unbounded dimensions. Note
+/// that this is a bounded set.
+///
+/// 3) Check if transformedSet is empty using the GBR sampling algorithm.
+///
+/// 4) return S is empty iff transformedSet is empty.
+///
+/// Since T is unimodular, a vector v is a solution to S*T iff T*v is a
+/// solution to S. The following is a sketch of a proof that S*T is empty
+/// iff transformedSet is empty:
+///
+/// If transformedSet is empty, then S*T is certainly empty since transformedSet
+/// was obtained by removing constraints and deleting dimensions from S*T.
+///
+/// If transformedSet contains a sample, consider the set C obtained by
+/// substituting the sample for the bounded dimensions of S*T. All the
+/// constraints of S*T that did not involve unbounded dimensions are
+/// satisfied by this substitution.
+///
+/// In step 1, all dimensions in the linear span of the dimensions outside the
+/// prefix are unbounded in S*T. Substituting values for the bounded dimensions
+/// cannot makes these dimensions bounded, and these are the only remaining
+/// dimensions in C, so C is unbounded along every vector. C is hence a
+/// full-dimensional cone and therefore always contains an integer point, which
+/// we can then substitute to get a full solution to S*T.
bool FlatAffineConstraints::isIntegerEmpty() const {
+ // First, try the GCD test heuristic.
if (isEmptyByGCDTest())
return true;
Simplex simplex(*this);
- if (simplex.isUnbounded())
- return false;
- return !simplex.findIntegerSample().hasValue();
+ if (simplex.isEmpty())
+ return true;
+
+ // For a bounded set, we directly call into the GBR sampling algorithm.
+ if (!simplex.isUnbounded())
+ return !simplex.findIntegerSample().hasValue();
+
+ // The set is unbounded. We cannot directly use the GBR algorithm.
+ //
+ // m is a matrix containing, in each row, a vector in which S is
+ // bounded, such that the linear span of all these dimensions contains all
+ // bounded dimensions in S.
+ Matrix m = getBoundedDirections();
+ // In column echelon form, each row of m occupies only the first rank(m)
+ // columns and has zeros on the other columns. The transform T that brings S
+ // to column echelon form is unimodular as well, so this is a suitable
+ // transform to use in step 1 of the algorithm.
+ std::pair<unsigned, LinearTransform> result =
+ LinearTransform::makeTransformToColumnEchelon(std::move(m));
+ FlatAffineConstraints transformedSet = result.second.applyTo(*this);
+
+ unsigned numBoundedDims = result.first;
+ unsigned numUnboundedDims = getNumIds() - numBoundedDims;
+ removeConstraintsInvolvingSuffixDims(transformedSet, numUnboundedDims);
+
+ // Remove all the unbounded dimensions.
+ transformedSet.removeIdRange(numBoundedDims, transformedSet.getNumIds());
+
+ return !Simplex(transformedSet).findIntegerSample().hasValue();
}
Optional<SmallVector<int64_t, 8>>
diff --git a/mlir/lib/Analysis/CMakeLists.txt b/mlir/lib/Analysis/CMakeLists.txt
index 3247ef1f56b0..585ba2aa8baf 100644
--- a/mlir/lib/Analysis/CMakeLists.txt
+++ b/mlir/lib/Analysis/CMakeLists.txt
@@ -3,6 +3,7 @@ set(LLVM_OPTIONAL_SOURCES
AffineStructures.cpp
BufferAliasAnalysis.cpp
CallGraph.cpp
+ LinearTransform.cpp
Liveness.cpp
LoopAnalysis.cpp
NestedMatcher.cpp
@@ -36,6 +37,7 @@ add_mlir_library(MLIRAnalysis
add_mlir_library(MLIRLoopAnalysis
AffineAnalysis.cpp
AffineStructures.cpp
+ LinearTransform.cpp
LoopAnalysis.cpp
NestedMatcher.cpp
PresburgerSet.cpp
diff --git a/mlir/lib/Analysis/LinearTransform.cpp b/mlir/lib/Analysis/LinearTransform.cpp
new file mode 100644
index 000000000000..7176cb01231f
--- /dev/null
+++ b/mlir/lib/Analysis/LinearTransform.cpp
@@ -0,0 +1,156 @@
+//===- LinearTransform.cpp - MLIR LinearTransform Class -------------------===//
+//
+// 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/LinearTransform.h"
+#include "mlir/Analysis/AffineStructures.h"
+
+namespace mlir {
+
+LinearTransform::LinearTransform(Matrix &&oMatrix) : matrix(oMatrix) {}
+LinearTransform::LinearTransform(const Matrix &oMatrix) : matrix(oMatrix) {}
+
+// Set M(row, targetCol) to its remainder on division by M(row, sourceCol)
+// by subtracting from column targetCol an appropriate integer multiple of
+// sourceCol. This brings M(row, targetCol) to the range [0, M(row, sourceCol)).
+// Apply the same column operation to otherMatrix, with the same integer
+// multiple.
+static void modEntryColumnOperation(Matrix &m, unsigned row, unsigned sourceCol,
+ unsigned targetCol, Matrix &otherMatrix) {
+ assert(m(row, sourceCol) != 0 && "Cannot divide by zero!");
+ assert((m(row, sourceCol) > 0 && m(row, targetCol) > 0) &&
+ "Operands must be positive!");
+ int64_t ratio = m(row, targetCol) / m(row, sourceCol);
+ m.addToColumn(sourceCol, targetCol, -ratio);
+ otherMatrix.addToColumn(sourceCol, targetCol, -ratio);
+}
+
+std::pair<unsigned, LinearTransform>
+LinearTransform::makeTransformToColumnEchelon(Matrix m) {
+ // We start with an identity result matrix and perform operations on m
+ // until m is in column echelon form. We apply the same sequence of operations
+ // on resultMatrix to obtain a transform that takes m to column echelon
+ // form.
+ Matrix resultMatrix = Matrix::identity(m.getNumColumns());
+
+ unsigned echelonCol = 0;
+ // Invariant: in all rows above row, all columns from echelonCol onwards
+ // are all zero elements. In an iteration, if the curent row has any non-zero
+ // elements echelonCol onwards, we bring one to echelonCol and use it to
+ // make all elements echelonCol + 1 onwards zero.
+ for (unsigned row = 0; row < m.getNumRows(); ++row) {
+ // Search row for a non-empty entry, starting at echelonCol.
+ unsigned nonZeroCol = echelonCol;
+ for (unsigned e = m.getNumColumns(); nonZeroCol < e; ++nonZeroCol) {
+ if (m(row, nonZeroCol) == 0)
+ continue;
+ break;
+ }
+
+ // Continue to the next row with the same echelonCol if this row is all
+ // zeros from echelonCol onwards.
+ if (nonZeroCol == m.getNumColumns())
+ continue;
+
+ // Bring the non-zero column to echelonCol. This doesn't affect rows
+ // above since they are all zero at these columns.
+ if (nonZeroCol != echelonCol) {
+ m.swapColumns(nonZeroCol, echelonCol);
+ resultMatrix.swapColumns(nonZeroCol, echelonCol);
+ }
+
+ // Make m(row, echelonCol) non-negative.
+ if (m(row, echelonCol) < 0) {
+ m.negateColumn(echelonCol);
+ resultMatrix.negateColumn(echelonCol);
+ }
+
+ // Make all the entries in row after echelonCol zero.
+ for (unsigned i = echelonCol + 1, e = m.getNumColumns(); i < e; ++i) {
+ // We make m(row, i) non-negative, and then apply the Euclidean GCD
+ // algorithm to (row, i) and (row, echelonCol). At the end, one of them
+ // has value equal to the gcd of the two entries, and the other is zero.
+
+ if (m(row, i) < 0) {
+ m.negateColumn(i);
+ resultMatrix.negateColumn(i);
+ }
+
+ unsigned targetCol = i, sourceCol = echelonCol;
+ // At every step, we set m(row, targetCol) %= m(row, sourceCol), and
+ // swap the indices sourceCol and targetCol. (not the columns themselves)
+ // This modulo is implemented as a subtraction
+ // m(row, targetCol) -= quotient * m(row, sourceCol),
+ // where quotient = floor(m(row, targetCol) / m(row, sourceCol)),
+ // which brings m(row, targetCol) to the range [0, m(row, sourceCol)).
+ //
+ // We are only allowed column operations; we perform the above
+ // for every row, i.e., the above subtraction is done as a column
+ // operation. This does not affect any rows above us since they are
+ // guaranteed to be zero at these columns.
+ while (m(row, targetCol) != 0 && m(row, sourceCol) != 0) {
+ modEntryColumnOperation(m, row, sourceCol, targetCol, resultMatrix);
+ std::swap(targetCol, sourceCol);
+ }
+
+ // One of (row, echelonCol) and (row, i) is zero and the other is the gcd.
+ // Make it so that (row, echelonCol) holds the non-zero value.
+ if (m(row, echelonCol) == 0) {
+ m.swapColumns(i, echelonCol);
+ resultMatrix.swapColumns(i, echelonCol);
+ }
+ }
+
+ ++echelonCol;
+ }
+
+ return {echelonCol, LinearTransform(std::move(resultMatrix))};
+}
+
+SmallVector<int64_t, 8> LinearTransform::applyTo(ArrayRef<int64_t> v) {
+ assert(v.size() == matrix.getNumRows() &&
+ "vector dimension should be matrix output dimension");
+
+ SmallVector<int64_t, 8> result;
+ result.reserve(v.size());
+ for (unsigned col = 0, e = matrix.getNumColumns(); col < e; ++col) {
+ int64_t elem = 0;
+ for (unsigned i = 0, e = matrix.getNumRows(); i < e; ++i)
+ elem += v[i] * matrix(i, col);
+ result.push_back(elem);
+ }
+ return result;
+}
+
+FlatAffineConstraints
+LinearTransform::applyTo(const FlatAffineConstraints &fac) {
+ FlatAffineConstraints result(fac.getNumDimIds());
+
+ for (unsigned i = 0, e = fac.getNumEqualities(); i < e; ++i) {
+ ArrayRef<int64_t> eq = fac.getEquality(i);
+
+ int64_t c = eq.back();
+
+ SmallVector<int64_t, 8> newEq = applyTo(eq.drop_back());
+ newEq.push_back(c);
+ result.addEquality(newEq);
+ }
+
+ for (unsigned i = 0, e = fac.getNumInequalities(); i < e; ++i) {
+ ArrayRef<int64_t> ineq = fac.getInequality(i);
+
+ int64_t c = ineq.back();
+
+ SmallVector<int64_t, 8> newIneq = applyTo(ineq.drop_back());
+ newIneq.push_back(c);
+ result.addInequality(newIneq);
+ }
+
+ return result;
+}
+
+} // namespace mlir
diff --git a/mlir/lib/Analysis/Presburger/CMakeLists.txt b/mlir/lib/Analysis/Presburger/CMakeLists.txt
index 49cdd5ac1431..2561013696d9 100644
--- a/mlir/lib/Analysis/Presburger/CMakeLists.txt
+++ b/mlir/lib/Analysis/Presburger/CMakeLists.txt
@@ -1,4 +1,4 @@
add_mlir_library(MLIRPresburger
Simplex.cpp
Matrix.cpp
- )
\ No newline at end of file
+ )
diff --git a/mlir/lib/Analysis/Presburger/Matrix.cpp b/mlir/lib/Analysis/Presburger/Matrix.cpp
index 213f1111e2a3..4a5a53921548 100644
--- a/mlir/lib/Analysis/Presburger/Matrix.cpp
+++ b/mlir/lib/Analysis/Presburger/Matrix.cpp
@@ -79,6 +79,19 @@ void Matrix::addToRow(unsigned sourceRow, unsigned targetRow, int64_t scale) {
return;
}
+void Matrix::addToColumn(unsigned sourceColumn, unsigned targetColumn,
+ int64_t scale) {
+ if (scale == 0)
+ return;
+ for (unsigned row = 0, e = getNumRows(); row < e; ++row)
+ at(row, targetColumn) += scale * at(row, sourceColumn);
+}
+
+void Matrix::negateColumn(unsigned column) {
+ for (unsigned row = 0, e = getNumRows(); row < e; ++row)
+ at(row, column) = -at(row, column);
+}
+
void Matrix::print(raw_ostream &os) const {
for (unsigned row = 0; row < nRows; ++row) {
for (unsigned column = 0; column < nColumns; ++column)
diff --git a/mlir/lib/Analysis/Presburger/Simplex.cpp b/mlir/lib/Analysis/Presburger/Simplex.cpp
index 47e199baba2a..2cfe5929e21d 100644
--- a/mlir/lib/Analysis/Presburger/Simplex.cpp
+++ b/mlir/lib/Analysis/Presburger/Simplex.cpp
@@ -9,6 +9,7 @@
#include "mlir/Analysis/Presburger/Simplex.h"
#include "mlir/Analysis/Presburger/Matrix.h"
#include "mlir/Support/MathExtras.h"
+#include "llvm/ADT/Optional.h"
namespace mlir {
using Direction = Simplex::Direction;
@@ -482,7 +483,7 @@ Optional<Fraction> Simplex::computeRowOptimum(Direction direction,
/// or None if it is unbounded.
Optional<Fraction> Simplex::computeOptimum(Direction direction,
ArrayRef<int64_t> coeffs) {
- assert(!empty && "Tableau should not be empty");
+ assert(!empty && "Simplex should not be empty");
unsigned snapshot = getSnapshot();
unsigned conIndex = addRow(coeffs);
@@ -492,6 +493,34 @@ Optional<Fraction> Simplex::computeOptimum(Direction direction,
return optimum;
}
+Optional<Fraction> Simplex::computeOptimum(Direction direction, Unknown &u) {
+ assert(!empty && "Simplex should not be empty!");
+ if (u.orientation == Orientation::Column) {
+ unsigned column = u.pos;
+ Optional<unsigned> pivotRow = findPivotRow({}, direction, column);
+ // If no pivot is returned, the constraint is unbounded in the specified
+ // direction.
+ if (!pivotRow)
+ return {};
+ pivot(*pivotRow, column);
+ }
+
+ unsigned row = u.pos;
+ Optional<Fraction> optimum = computeRowOptimum(direction, row);
+ if (u.restricted && direction == Direction::Down &&
+ (!optimum || *optimum < Fraction(0, 1)))
+ restoreRow(u);
+ return optimum;
+}
+
+bool Simplex::isBoundedAlongConstraint(unsigned constraintIndex) {
+ assert(!empty && "It is not meaningful to ask whether a direction is bounded "
+ "in an empty set.");
+ // The constraint's perpendicular is already bounded below, since it is a
+ // constraint. If it is also bounded above, we can return true.
+ return computeOptimum(Direction::Up, con[constraintIndex]).hasValue();
+}
+
/// Redundant constraints are those that are in row orientation and lie in
/// rows 0 to nRedundant - 1.
bool Simplex::isMarkedRedundant(unsigned constraintIndex) const {
diff --git a/mlir/unittests/Analysis/AffineStructuresTest.cpp b/mlir/unittests/Analysis/AffineStructuresTest.cpp
index 6fcb1c489cfc..ac11c90ec15b 100644
--- a/mlir/unittests/Analysis/AffineStructuresTest.cpp
+++ b/mlir/unittests/Analysis/AffineStructuresTest.cpp
@@ -15,22 +15,36 @@
namespace mlir {
-/// If 'hasValue' is true, check that findIntegerSample returns a valid sample
+enum class TestFunction { Sample, Empty };
+
+/// If fn is TestFunction::Sample (default):
+/// If hasSample is true, check that findIntegerSample returns a valid sample
/// for the FlatAffineConstraints fac.
+/// If hasSample is false, check that findIntegerSample returns None.
///
-/// If hasValue is false, check that findIntegerSample does not return None.
-static void checkSample(bool hasValue, const FlatAffineConstraints &fac) {
- Optional<SmallVector<int64_t, 8>> maybeSample = fac.findIntegerSample();
- if (!hasValue) {
- EXPECT_FALSE(maybeSample.hasValue());
- if (maybeSample.hasValue()) {
- for (auto x : *maybeSample)
- llvm::errs() << x << ' ';
- llvm::errs() << '\n';
+/// If fn is TestFunction::Empty, check that isIntegerEmpty returns the
+/// opposite of hasSample.
+static void checkSample(bool hasSample, const FlatAffineConstraints &fac,
+ TestFunction fn = TestFunction::Sample) {
+ Optional<SmallVector<int64_t, 8>> maybeSample;
+ switch (fn) {
+ case TestFunction::Sample:
+ maybeSample = fac.findIntegerSample();
+ if (!hasSample) {
+ EXPECT_FALSE(maybeSample.hasValue());
+ if (maybeSample.hasValue()) {
+ for (auto x : *maybeSample)
+ llvm::errs() << x << ' ';
+ llvm::errs() << '\n';
+ }
+ } else {
+ ASSERT_TRUE(maybeSample.hasValue());
+ EXPECT_TRUE(fac.containsPoint(*maybeSample));
}
- } else {
- ASSERT_TRUE(maybeSample.hasValue());
- EXPECT_TRUE(fac.containsPoint(*maybeSample));
+ break;
+ case TestFunction::Empty:
+ EXPECT_EQ(!hasSample, fac.isIntegerEmpty());
+ break;
}
}
@@ -52,9 +66,11 @@ makeFACFromConstraints(unsigned dims, ArrayRef<SmallVector<int64_t, 4>> ineqs,
/// orderings may cause the algorithm to proceed
diff erently. At least some of
///.these permutations should make it past the heuristics and test the
/// implementation of the GBR algorithm itself.
-static void checkPermutationsSample(bool hasValue, unsigned nDim,
+/// Use TestFunction fn to test.
+static void checkPermutationsSample(bool hasSample, unsigned nDim,
ArrayRef<SmallVector<int64_t, 4>> ineqs,
- ArrayRef<SmallVector<int64_t, 4>> eqs) {
+ ArrayRef<SmallVector<int64_t, 4>> eqs,
+ TestFunction fn = TestFunction::Sample) {
SmallVector<unsigned, 4> perm(nDim);
std::iota(perm.begin(), perm.end(), 0);
auto permute = [&perm](ArrayRef<int64_t> coeffs) {
@@ -71,8 +87,8 @@ static void checkPermutationsSample(bool hasValue, unsigned nDim,
for (const auto &eq : eqs)
permutedEqs.push_back(permute(eq));
- checkSample(hasValue,
- makeFACFromConstraints(nDim, permutedIneqs, permutedEqs));
+ checkSample(hasSample,
+ makeFACFromConstraints(nDim, permutedIneqs, permutedEqs), fn);
} while (std::next_permutation(perm.begin(), perm.end()));
}
@@ -206,19 +222,158 @@ TEST(FlatAffineConstraintsTest, IsIntegerEmptyTest) {
EXPECT_FALSE(
makeFACFromConstraints(1, {{5, -1}, {-5, 9}}, {}).isIntegerEmpty());
- // An unbounded set, which isIntegerEmpty should detect as unbounded and
- // return without calling findIntegerSample.
+ // Unbounded sets.
+ EXPECT_TRUE(makeFACFromConstraints(3,
+ {
+ {2, 0, 0, -1}, // 2x >= 1
+ {-2, 0, 0, 1}, // 2x <= 1
+ {0, 2, 0, -1}, // 2y >= 1
+ {0, -2, 0, 1}, // 2y <= 1
+ {0, 0, 2, -1}, // 2z >= 1
+ },
+ {})
+ .isIntegerEmpty());
+
EXPECT_FALSE(makeFACFromConstraints(3,
{
- {2, 0, 0, -1},
- {-2, 0, 0, 1},
- {0, 2, 0, -1},
- {0, -2, 0, 1},
- {0, 0, 2, -1},
+ {2, 0, 0, -1}, // 2x >= 1
+ {-3, 0, 0, 3}, // 3x <= 3
+ {0, 0, 5, -6}, // 5z >= 6
+ {0, 0, -7, 17}, // 7z <= 17
+ {0, 3, 0, -2}, // 3y >= 2
},
{})
.isIntegerEmpty());
+ // 2D cone with apex at (10000, 10000) and
+ // edges passing through (1/3, 0) and (2/3, 0).
+ EXPECT_FALSE(
+ makeFACFromConstraints(
+ 2, {{300000, -299999, -100000}, {-300000, 299998, 200000}}, {})
+ .isIntegerEmpty());
+
+ // Cartesian product of a tetrahedron and a 2D cone.
+ // The tetrahedron has vertices at
+ // (1/3, 0, 0), (2/3, 0, 0), (2/3, 0, 10000), and (10000, 10000, 10000).
+ // The first three points form a triangular base on the xz plane with the
+ // apex at the fourth point, which is the only integer point.
+ // The cone has apex at (10000, 10000) and
+ // edges passing through (1/3, 0) and (2/3, 0).
+ checkPermutationsSample(
+ true /* not empty */, 5,
+ {
+ // Tetrahedron contraints:
+ {0, 1, 0, 0, 0, 0}, // y >= 0
+ {0, -1, 1, 0, 0, 0}, // z >= y
+ // -300000x + 299998y + 100000 + z <= 0.
+ {300000, -299998, -1, 0, 0, -100000},
+ // -150000x + 149999y + 100000 >= 0.
+ {-150000, 149999, 0, 0, 0, 100000},
+
+ // Triangle constraints:
+ // 300000p - 299999q >= 100000
+ {0, 0, 0, 300000, -299999, -100000},
+ // -300000p + 299998q + 200000 >= 0
+ {0, 0, 0, -300000, 299998, 200000},
+ },
+ {}, TestFunction::Empty);
+
+ // Cartesian product of same tetrahedron as above and {(p, q) : 1/3 <= p <=
+ // 2/3}. Since the second set is empty, the whole set is too.
+ checkPermutationsSample(
+ false /* empty */, 5,
+ {
+ // Tetrahedron contraints:
+ {0, 1, 0, 0, 0, 0}, // y >= 0
+ {0, -1, 1, 0, 0, 0}, // z >= y
+ // -300000x + 299998y + 100000 + z <= 0.
+ {300000, -299998, -1, 0, 0, -100000},
+ // -150000x + 149999y + 100000 >= 0.
+ {-150000, 149999, 0, 0, 0, 100000},
+
+ // Second set constraints:
+ // 3p >= 1
+ {0, 0, 0, 3, 0, -1},
+ // 3p <= 2
+ {0, 0, 0, -3, 0, 2},
+ },
+ {}, TestFunction::Empty);
+
+ // Cartesian product of same tetrahedron as above and
+ // {(p, q, r) : 1 <= p <= 2 and p = 3q + 3r}.
+ // Since the second set is empty, the whole set is too.
+ checkPermutationsSample(
+ false /* empty */, 5,
+ {
+ // Tetrahedron contraints:
+ {0, 1, 0, 0, 0, 0, 0}, // y >= 0
+ {0, -1, 1, 0, 0, 0, 0}, // z >= y
+ // -300000x + 299998y + 100000 + z <= 0.
+ {300000, -299998, -1, 0, 0, 0, -100000},
+ // -150000x + 149999y + 100000 >= 0.
+ {-150000, 149999, 0, 0, 0, 0, 100000},
+
+ // Second set constraints:
+ // p >= 1
+ {0, 0, 0, 1, 0, 0, -1},
+ // p <= 2
+ {0, 0, 0, -1, 0, 0, 2},
+ },
+ {
+ {0, 0, 0, 1, -3, -3, 0}, // p = 3q + 3r
+ },
+ TestFunction::Empty);
+
+ // Cartesian product of a tetrahedron and a 2D cone.
+ // The tetrahedron is empty and has vertices at
+ // (1/3, 0, 0), (2/3, 0, 0), (2/3, 0, 100), and (100, 100 - 1/3, 100).
+ // The cone has apex at (10000, 10000) and
+ // edges passing through (1/3, 0) and (2/3, 0).
+ // Since the tetrahedron is empty, the Cartesian product is too.
+ checkPermutationsSample(false /* empty */, 5,
+ {
+ // Tetrahedron contraints:
+ {0, 1, 0, 0, 0, 0},
+ {0, -300, 299, 0, 0, 0},
+ {300 * 299, -89400, -299, 0, 0, -100 * 299},
+ {-897, 894, 0, 0, 0, 598},
+
+ // Triangle constraints:
+ // 300000p - 299999q >= 100000
+ {0, 0, 0, 300000, -299999, -100000},
+ // -300000p + 299998q + 200000 >= 0
+ {0, 0, 0, -300000, 299998, 200000},
+ },
+ {}, TestFunction::Empty);
+
+ // Cartesian product of same tetrahedron as above and
+ // {(p, q) : 1/3 <= p <= 2/3}.
+ checkPermutationsSample(false /* empty */, 5,
+ {
+ // Tetrahedron contraints:
+ {0, 1, 0, 0, 0, 0},
+ {0, -300, 299, 0, 0, 0},
+ {300 * 299, -89400, -299, 0, 0, -100 * 299},
+ {-897, 894, 0, 0, 0, 598},
+
+ // Second set constraints:
+ // 3p >= 1
+ {0, 0, 0, 3, 0, -1},
+ // 3p <= 2
+ {0, 0, 0, -3, 0, 2},
+ },
+ {}, TestFunction::Empty);
+
+ EXPECT_FALSE(makeFACFromConstraints(3,
+ {
+ {2, 0, 0, -1}, // 2x >= 1
+ },
+ {{
+ {1, -1, 0, -1}, // y = x - 1
+ {0, 1, -1, 0}, // z = y
+ }})
+ .isIntegerEmpty());
+
// FlatAffineConstraints::isEmpty() does not detect the following sets to be
// empty.
diff --git a/mlir/unittests/Analysis/CMakeLists.txt b/mlir/unittests/Analysis/CMakeLists.txt
index 6317aeb8df89..0df0af866d66 100644
--- a/mlir/unittests/Analysis/CMakeLists.txt
+++ b/mlir/unittests/Analysis/CMakeLists.txt
@@ -1,5 +1,6 @@
add_mlir_unittest(MLIRAnalysisTests
AffineStructuresTest.cpp
+ LinearTransformTest.cpp
PresburgerSetTest.cpp
)
diff --git a/mlir/unittests/Analysis/LinearTransformTest.cpp b/mlir/unittests/Analysis/LinearTransformTest.cpp
new file mode 100644
index 000000000000..598c84920d5d
--- /dev/null
+++ b/mlir/unittests/Analysis/LinearTransformTest.cpp
@@ -0,0 +1,87 @@
+//===- LinearTransformTest.cpp - Tests for LinearTransform ----------------===//
+//
+// 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/LinearTransform.h"
+#include <gmock/gmock.h>
+#include <gtest/gtest.h>
+
+namespace mlir {
+
+void testColumnEchelonForm(const Matrix &m, unsigned expectedRank) {
+ unsigned lastAllowedNonZeroCol = 0;
+ std::pair<unsigned, LinearTransform> result =
+ LinearTransform::makeTransformToColumnEchelon(m);
+ unsigned rank = result.first;
+ EXPECT_EQ(rank, expectedRank);
+ LinearTransform transform = result.second;
+ // In column echelon form, each row's last non-zero value can be at most one
+ // column to the right of the last non-zero column among the previous rows.
+ for (unsigned row = 0, nRows = m.getNumRows(); row < nRows; ++row) {
+ SmallVector<int64_t, 8> rowVec = transform.applyTo(m.getRow(row));
+ for (unsigned col = lastAllowedNonZeroCol + 1, nCols = m.getNumColumns();
+ col < nCols; ++col) {
+ EXPECT_EQ(rowVec[col], 0);
+ if (rowVec[col] != 0) {
+ llvm::errs() << "Failed at input matrix:\n";
+ m.dump();
+ }
+ }
+ if (rowVec[lastAllowedNonZeroCol] != 0)
+ lastAllowedNonZeroCol++;
+ }
+ // The final value of lastAllowedNonZeroCol is the index of the first
+ // all-zeros column, so it must be equal to the rank.
+ EXPECT_EQ(lastAllowedNonZeroCol, rank);
+}
+
+TEST(LinearTransformTest, transformToColumnEchelonTest) {
+ // m1, m2, m3 are rank 1 matrices -- the first and second rows are identical.
+ Matrix m1(2, 2);
+ m1(0, 0) = 4;
+ m1(0, 1) = -7;
+ m1(1, 0) = 4;
+ m1(1, 1) = -7;
+ testColumnEchelonForm(m1, 1u);
+
+ Matrix m2(2, 2);
+ m2(0, 0) = -4;
+ m2(0, 1) = 7;
+ m2(1, 0) = 4;
+ m2(1, 1) = -7;
+ testColumnEchelonForm(m2, 1u);
+
+ Matrix m3(2, 2);
+ m3(0, 0) = -4;
+ m3(0, 1) = -7;
+ m3(1, 0) = -4;
+ m3(1, 1) = -7;
+ testColumnEchelonForm(m3, 1u);
+
+ // m4, m5, m6 are rank 2 matrices -- the first and second rows are
diff erent.
+ Matrix m4(2, 2);
+ m4(0, 0) = 4;
+ m4(0, 1) = -7;
+ m4(1, 0) = -4;
+ m4(1, 1) = -7;
+ testColumnEchelonForm(m4, 2u);
+
+ Matrix m5(2, 2);
+ m5(0, 0) = -4;
+ m5(0, 1) = 7;
+ m5(1, 0) = 4;
+ m5(1, 1) = 7;
+ testColumnEchelonForm(m5, 2u);
+
+ Matrix m6(2, 2);
+ m6(0, 0) = -4;
+ m6(0, 1) = -7;
+ m6(1, 0) = 4;
+ m6(1, 1) = -7;
+ testColumnEchelonForm(m5, 2u);
+}
+} // namespace mlir
More information about the llvm-branch-commits
mailing list