Google OR-Tools v9.11
a fast and portable software suite for combinatorial optimization
|
#include <sparse.h>
Public Member Functions | |
TriangularMatrix () | |
TriangularMatrix (const TriangularMatrix &)=delete | |
This type is neither copyable nor movable. | |
TriangularMatrix & | operator= (const TriangularMatrix &)=delete |
void | PopulateFromTranspose (const TriangularMatrix &input) |
void | Swap (TriangularMatrix *other) |
bool | IsEmpty () const |
RowIndex | num_rows () const |
ColIndex | num_cols () const |
EntryIndex | num_entries () const |
void | Reset (RowIndex num_rows, ColIndex col_capacity) |
void | PopulateFromTriangularSparseMatrix (const SparseMatrix &input) |
void | AddTriangularColumn (const ColumnView &column, RowIndex diagonal_row) |
void | AddTriangularColumnWithGivenDiagonalEntry (const SparseColumn &column, RowIndex diagonal_row, Fractional diagonal_value) |
void | AddDiagonalOnlyColumn (Fractional diagonal_value) |
void | AddAndNormalizeTriangularColumn (const SparseColumn &column, RowIndex diagonal_row, Fractional diagonal_coefficient) |
void | ApplyRowPermutationToNonDiagonalEntries (const RowPermutation &row_perm) |
Applies the given row permutation to all entries except the diagonal ones. | |
void | CopyColumnToSparseColumn (ColIndex col, SparseColumn *output) const |
Copy a triangular column with its diagonal entry to the given SparseColumn. | |
void | CopyToSparseMatrix (SparseMatrix *output) const |
Copy a triangular matrix to the given SparseMatrix. | |
ColIndex | GetFirstNonIdentityColumn () const |
Fractional | GetDiagonalCoefficient (ColIndex col) const |
Returns the diagonal coefficient of the given column. | |
bool | ColumnIsDiagonalOnly (ColIndex col) const |
Returns true iff the column contains no non-diagonal entries. | |
void | LowerSolve (DenseColumn *rhs) const |
void | UpperSolve (DenseColumn *rhs) const |
Solves the system U.x = rhs for an upper triangular matrix. | |
void | TransposeUpperSolve (DenseColumn *rhs) const |
void | LowerSolveStartingAt (ColIndex start, DenseColumn *rhs) const |
This assumes that the rhs is all zero before the given position. | |
void | TransposeLowerSolve (DenseColumn *rhs) const |
void | HyperSparseSolve (DenseColumn *rhs, RowIndexVector *non_zero_rows) const |
void | HyperSparseSolveWithReversedNonZeros (DenseColumn *rhs, RowIndexVector *non_zero_rows) const |
void | TransposeHyperSparseSolve (DenseColumn *rhs, RowIndexVector *non_zero_rows) const |
void | TransposeHyperSparseSolveWithReversedNonZeros (DenseColumn *rhs, RowIndexVector *non_zero_rows) const |
void | ComputeRowsToConsiderWithDfs (RowIndexVector *non_zero_rows) const |
void | ComputeRowsToConsiderInSortedOrder (RowIndexVector *non_zero_rows) const |
void | PermutedLowerSolve (const SparseColumn &rhs, const RowPermutation &row_perm, const RowMapping &partial_inverse_row_perm, SparseColumn *lower, SparseColumn *upper) const |
void | PermutedLowerSparseSolve (const ColumnView &rhs, const RowPermutation &row_perm, SparseColumn *lower, SparseColumn *upper) |
int64_t | NumFpOperationsInLastPermutedLowerSparseSolve () const |
This is used to compute the deterministic time of a matrix factorization. | |
bool | IsLowerTriangular () const |
bool | IsUpperTriangular () const |
void | PermutedComputeRowsToConsider (const ColumnView &rhs, const RowPermutation &row_perm, RowIndexVector *lower_column_rows, RowIndexVector *upper_column_rows) |
Fractional | ComputeInverseInfinityNormUpperBound () const |
Fractional | ComputeInverseInfinityNorm () const |
Specialization of a CompactSparseMatrix used for triangular matrices. To be able to solve triangular systems as efficiently as possible, the diagonal entries are stored in a separate vector and not in the underlying CompactSparseMatrix.
Advanced usage: this class also support matrices that can be permuted into a triangular matrix and some functions work directly on such matrices.
|
inline |
|
delete |
This type is neither copyable nor movable.
void operations_research::glop::TriangularMatrix::AddAndNormalizeTriangularColumn | ( | const SparseColumn & | column, |
RowIndex | diagonal_row, | ||
Fractional | diagonal_coefficient ) |
Adds the given sparse column divided by diagonal_coefficient. The diagonal_row is assumed to be present and its value should be the same as the one given in diagonal_coefficient. Note that this function tests for zero coefficients in the input column and removes them.
void operations_research::glop::TriangularMatrix::AddDiagonalOnlyColumn | ( | Fractional | diagonal_value | ) |
void operations_research::glop::TriangularMatrix::AddTriangularColumn | ( | const ColumnView & | column, |
RowIndex | diagonal_row ) |
Functions to create a triangular matrix incrementally, column by column. A client needs to call Reset(num_rows) first, and then each column must be added by calling one of the 3 functions below.
void operations_research::glop::TriangularMatrix::AddTriangularColumnWithGivenDiagonalEntry | ( | const SparseColumn & | column, |
RowIndex | diagonal_row, | ||
Fractional | diagonal_value ) |
void operations_research::glop::TriangularMatrix::ApplyRowPermutationToNonDiagonalEntries | ( | const RowPermutation & | row_perm | ) |
|
inline |
Fractional operations_research::glop::TriangularMatrix::ComputeInverseInfinityNorm | ( | ) | const |
Fractional operations_research::glop::TriangularMatrix::ComputeInverseInfinityNormUpperBound | ( | ) | const |
The upper bound is computed using one of the algorithm presented in "A Survey of Condition Number Estimation for Triangular Matrices" https:epubs.siam.org/doi/pdf/10.1137/1029112/
A known upper bound for the infinity norm of T^{-1} is the infinity norm of y where T'*y = x with:
Identity matrix
void operations_research::glop::TriangularMatrix::ComputeRowsToConsiderInSortedOrder | ( | RowIndexVector * | non_zero_rows | ) | const |
Same as TriangularComputeRowsToConsider() but always returns the non-zeros sorted by rows. It is up to the client to call the direct or reverse hyper-sparse solve function depending if the matrix is upper or lower triangular.
void operations_research::glop::TriangularMatrix::ComputeRowsToConsiderWithDfs | ( | RowIndexVector * | non_zero_rows | ) | const |
Given the positions of the non-zeros of a vector, computes the non-zero positions of the vector after a solve by this triangular matrix. The order of the returned non-zero positions will be in the REVERSE elimination order. If the function detects that there are too many non-zeros, then it aborts early and non_zero_rows is cleared.
We don't start the DFS if the initial number of non-zeros is under the sparsity_threshold. During the DFS, we abort it if the number of floating points operations get larger than the num_ops_threshold.
In both cases, we make sure to clear non_zero_rows so that the solving part will use the non-hypersparse version of the code.
Initialize using the non-zero positions of the input.
Topological sort based on Depth-First-Search. Same remarks as the version implemented in PermutedComputeRowsToConsider().
If the depth-first search from the current node is finished, we store the node. This will store the node in reverse topological order.
If the node is already stored, skip.
Go one level forward in the depth-first search, and store the 'adjacent' node on nodes_to_explore_ for further processing.
We reverse the sign of nodes_to_explore_.back() to detect when the DFS will be back on this node.
Abort if the number of operations is not negligible compared to the number of rows. Note that this test also prevents the code from cycling in case the matrix is actually not triangular.
Clear stored_.
If we aborted, clear the result.
void operations_research::glop::TriangularMatrix::CopyColumnToSparseColumn | ( | ColIndex | col, |
SparseColumn * | output ) const |
Copy a triangular column with its diagonal entry to the given SparseColumn.
void operations_research::glop::TriangularMatrix::CopyToSparseMatrix | ( | SparseMatrix * | output | ) | const |
Copy a triangular matrix to the given SparseMatrix.
|
inline |
|
inline |
void operations_research::glop::TriangularMatrix::HyperSparseSolve | ( | DenseColumn * | rhs, |
RowIndexVector * | non_zero_rows ) const |
Hyper-sparse version of the triangular solve functions. The passed non_zero_rows should contain the positions of the symbolic non-zeros of the result in the order in which they need to be accessed (or in the reverse order for the Reverse*() versions).
The non-zero vector is mutable so that the symbolic non-zeros that are actually zero because of numerical cancellations can be removed.
The non-zeros can be computed by one of these two methods:
For a given solve, here is the required order:
For a general discussion of hyper-sparsity in LP, see: J.A.J. Hall, K.I.M. McKinnon, "Exploiting hyper-sparsity in the revised simplex method", December 1999, MS 99-014. http://www.maths.ed.ac.uk/hall/MS-99/MS9914.pdf
void operations_research::glop::TriangularMatrix::HyperSparseSolveWithReversedNonZeros | ( | DenseColumn * | rhs, |
RowIndexVector * | non_zero_rows ) const |
|
inline |
bool operations_research::glop::TriangularMatrix::IsLowerTriangular | ( | ) | const |
To be used in DEBUG mode by the client code. This check that the matrix is lower- (resp. upper-) triangular without any permutation and that there is no zero on the diagonal. We can't do that on each Solve() that require so, otherwise it will be too slow in debug.
bool operations_research::glop::TriangularMatrix::IsUpperTriangular | ( | ) | const |
void operations_research::glop::TriangularMatrix::LowerSolve | ( | DenseColumn * | rhs | ) | const |
void operations_research::glop::TriangularMatrix::LowerSolveStartingAt | ( | ColIndex | start, |
DenseColumn * | rhs ) const |
|
inline |
|
inline |
|
inline |
|
inline |
|
delete |
void operations_research::glop::TriangularMatrix::PermutedComputeRowsToConsider | ( | const ColumnView & | rhs, |
const RowPermutation & | row_perm, | ||
RowIndexVector * | lower_column_rows, | ||
RowIndexVector * | upper_column_rows ) |
Visible for testing. This is used by PermutedLowerSparseSolve() to compute the non-zero indices of the result. The output is as follow:
This function is non-const because it prunes the underlying dependency graph of the lower matrix while doing a solve. See marked_ and pruned_ends_ below.
Pruning the graph at the same time is slower but not by too much (< 2x) and seems worth doing. Note that when the lower matrix is dense, most of the graph will likely be pruned. As a result, the symbolic phase will be negligible compared to the numerical phase so we don't really need a dense version of PermutedLowerSparseSolve().
The goal is to find which rows of the working column we will need to look at in PermutedLowerSparseSolve() when solving P^{-1}.L.P.x = rhs, 'P' being a row permutation, 'L' a lower triangular matrix and 'this' being 'P^{-1}.L'.
Let G denote the graph G = (V,E) of the column-to-row adjacency of A:
Let S denote the set of nodes i such that rhs_i != 0. Let R denote the set of all accessible nodes from S in G. x_k is possibly non-zero iff k is in R, i.e. if k is not in R then x_k = 0 for sure, and there is no need to look a the row k during the solve.
So, to solve P^{-1}.L.P.x = rhs, only rows corresponding to P.R have to be considered (ignoring the one that map to identity column of L). A topological sort of P.R is used to decide in which order one should iterate on them. This will be given by upper_column_rows_ and it will be populated in reverse order.
Topological sort based on Depth-First-Search. A few notes:
If the depth-first search from the current node is finished (i.e. there is a sentinel on the stack), we store the node (which is just before on the stack). This will store the nodes in reverse topological order.
Unmark and prune the nodes that are already unmarked. See the header comment on marked_ for the algorithm description.
Complexity note: The only difference with the "normal" DFS doing no pruning is this extra loop here and the marked_[entry_row] = true in the loop later in this function. On an already pruned graph, this is probably between 1 and 2 times slower than the "normal" DFS.
If the node is already stored, skip.
Expand only if we are not on a kNonPivotal row. Otherwise we can store the node right away.
Go one level forward in the depth-first search, and store the 'adjacent' node on nodes_to_explore_ for further processing.
The graph contains cycles? this is not supposed to happen.
Clear stored_.
void operations_research::glop::TriangularMatrix::PermutedLowerSolve | ( | const SparseColumn & | rhs, |
const RowPermutation & | row_perm, | ||
const RowMapping & | partial_inverse_row_perm, | ||
SparseColumn * | lower, | ||
SparseColumn * | upper ) const |
This is currently only used for testing. It achieves the same result as PermutedLowerSparseSolve() below, but the latter exploits the sparsity of rhs and is thus faster for our use case.
IMPORTANT: lower will contain all the "symbolic" non-zero entries. A "symbolic" zero entry is one that will be zero whatever the coefficients of the rhs entries. That is it only depends on the position of its entries, not on their values. Thus, some of its coefficients may be zero. This fact is exploited by the LU factorization code. The zero coefficients of upper will be cleaned, however.
void operations_research::glop::TriangularMatrix::PermutedLowerSparseSolve | ( | const ColumnView & | rhs, |
const RowPermutation & | row_perm, | ||
SparseColumn * | lower, | ||
SparseColumn * | upper ) |
This solves a lower triangular system with only ones on the diagonal where the matrix and the input rhs are permuted by the inverse of row_perm. Note that the output will also be permuted by the inverse of row_perm. The function also supports partial permutation. That is if row_perm[i] < 0 then column row_perm[i] is assumed to be an identity column.
The output is given as follow:
Compute the set of rows that will be non zero in the result (lower_column, upper_column).
Copy rhs into initially_all_zero_scratchpad_.
We clear lower_column first in case upper_column and lower_column point to the same underlying SparseColumn.
rows_to_consider_ contains the row to process in reverse order. Note in particular that each "permuted_row" will never be touched again and so its value is final. We copy the result in (lower_column, upper_column) and clear initially_all_zero_scratchpad_ at the same time.
void operations_research::glop::TriangularMatrix::PopulateFromTranspose | ( | const TriangularMatrix & | input | ) |
Only a subset of the functions from CompactSparseMatrix are exposed (note the private inheritance). They are extended to deal with diagonal coefficients properly.
This takes care of the triangular special case.
The elimination structure of the transpose is not the same.
Compute first_non_identity_column_. Note that this is not necessarily the same as input.first_non_identity_column_ for an upper triangular matrix.
void operations_research::glop::TriangularMatrix::PopulateFromTriangularSparseMatrix | ( | const SparseMatrix & | input | ) |
Constructs a triangular matrix from the given SparseMatrix. The input is assumed to be lower or upper triangular without any permutations. This is checked in debug mode.
void operations_research::glop::TriangularMatrix::Reset | ( | RowIndex | num_rows, |
ColIndex | col_capacity ) |
On top of the CompactSparseMatrix functionality, TriangularMatrix::Reset() also pre-allocates space of size col_size for a number of internal vectors. This helps reduce costly push_back operations for large problems.
Non-zero entries in the first column always have an offset of 0.
void operations_research::glop::TriangularMatrix::Swap | ( | TriangularMatrix * | other | ) |
void operations_research::glop::TriangularMatrix::TransposeHyperSparseSolve | ( | DenseColumn * | rhs, |
RowIndexVector * | non_zero_rows ) const |
void operations_research::glop::TriangularMatrix::TransposeHyperSparseSolveWithReversedNonZeros | ( | DenseColumn * | rhs, |
RowIndexVector * | non_zero_rows ) const |
void operations_research::glop::TriangularMatrix::TransposeLowerSolve | ( | DenseColumn * | rhs | ) | const |
void operations_research::glop::TriangularMatrix::TransposeUpperSolve | ( | DenseColumn * | rhs | ) | const |
void operations_research::glop::TriangularMatrix::UpperSolve | ( | DenseColumn * | rhs | ) | const |