CMR  1.3.0
Namespaces | Classes | Typedefs | Enumerations | Functions
tu Namespace Reference

Namespaces

 detail
 
 util
 

Classes

struct  abs_greater
 
struct  articulation_edge_filter
 
class  binary_linear_space
 
struct  bipartite_graph_bfs_node
 
struct  bipartite_graph_dimensions
 
class  bipartite_r10_graphs
 
class  combination
 
class  decomposed_matroid
 
class  decomposed_matroid_leaf
 
class  decomposed_matroid_separator
 
class  DenseMatrix
 Dense matrix with entries of type V. More...
 
struct  elaborate_extension_matrix_modifier
 
class  ghouila_houri_enumerator
 
struct  is_all_ones
 
struct  is_non_zero
 
class  last_combination_exception
 
class  logger
 
class  Matrix
 Base class for matrices whose entries have type V. More...
 
class  matrix_modified
 
class  matrix_permuted
 
struct  matrix_reorder_row_less
 
class  matrix_transposed
 
class  matroid
 
class  matroid_permuted
 
class  matroid_transposed
 
class  nested_minor_sequence
 
class  nested_minor_sequence_transposed
 
class  permutation
 
class  permutation_enumerator
 
class  permutation_shrink_exception
 
class  separation
 
class  SparseMatrix
 Sparse matrix with entries of type V. More...
 
class  Submatrix
 
struct  submatrix_indices
 
class  SubmatrixIndices
 Indexes a submatrix. More...
 
class  TransposedMatrix
 Matrix proxy for the transpose of a matrix. More...
 
struct  vector_less
 
class  vector_three_connectivity
 
struct  zero_block_matrix_modifier
 

Typedefs

typedef boost::numeric::ublas::matrix< long long > integer_matrix
 
typedef boost::numeric::ublas::matrix_indirect< const integer_matrix, submatrix_indices::indirect_array_typeinteger_submatrix
 
typedef matroid< int > integer_matroid
 A matroid with integer names. More...
 
typedef std::set< int > matroid_element_set
 
typedef boost::property< tu::edge_matroid_element_t, int > matroid_element_property
 
typedef boost::adjacency_list< boost::vecS, boost::vecS, boost::undirectedS, boost::no_property, matroid_element_propertymatroid_graph
 
typedef boost::property_map< matroid_graph, edge_matroid_element_t >::const_type const_matroid_element_map
 
typedef boost::property_map< matroid_graph, edge_matroid_element_t >::type matroid_element_map
 
typedef std::pair< size_t, size_t > size_pair_t
 

Enumerations

enum  log_level { LOG_QUIET, LOG_PROGRESSIVE, LOG_VERBOSE }
 
enum  edge_matroid_element_t { edge_matroid_element }
 
enum  rank_distribution { RANK_TOO_HIGH = 0, RANK_BL_1_TR_1 = 1, RANK_BL_2_TR_0 = 2, RANK_BL_0_TR_2 = 3 }
 

Functions

template<typename MatroidType , typename MatrixType >
std::pair< bool, decomposed_matroid * > decompose_binary_matroid (MatroidType &matroid, MatrixType &matrix, matroid_element_set extra_elements, bool construct_decomposition, logger &log)
 Forward declaration. More...
 
template<typename MatroidType , typename MatrixType >
std::pair< bool, decomposed_matroid * > decompose_with_minor_sequence (matroid_permuted< MatroidType > &permuted_matroid, matrix_permuted< MatrixType > &permuted_matrix, nested_minor_sequence &nested_minors, matroid_element_set extra_elements, bool construct_decomposition, logger &log)
 
template<typename MatroidType , typename MatrixType >
matroid_graphconstruct_small_matroid_graph (MatroidType &matroid, MatrixType &matrix)
 
std::ostream & operator<< (std::ostream &stream, const binary_linear_space &space)
 
template<typename MatrixType , typename StartNodesContainerType , typename EndNodesContainerType >
bool bipartite_graph_bfs (const MatrixType &matrix, const bipartite_graph_dimensions &dimensions, const StartNodesContainerType &start_nodes, const EndNodesContainerType &end_nodes, bool reach_all, std::vector< bipartite_graph_bfs_node > &result)
 
std::ostream & operator<< (std::ostream &stream, const combination &c)
 
int submatrix_determinant (const integer_matrix &matrix, const submatrix_indices &submatrix)
 
bool submatrix_camion (const integer_matrix &matrix, const submatrix_indices &submatrix)
 
bool determinant_is_totally_unimodular (const integer_matrix &matrix)
 
bool determinant_is_totally_unimodular (const integer_matrix &matrix, submatrix_indices &violator)
 
template<typename MatroidType , typename MatrixType , typename NestedMinorSequence >
separation enumerate_separations (MatroidType &matroid, MatrixType &matrix, const NestedMinorSequence &nested_minors, matroid_element_set &extra_elements, logger &log)
 
template<typename MatroidType , typename MatrixType , typename NestedMinorSequenceType , typename RowThreeConnectivity , typename ColumnThreeConnectivity >
bool find_simple_row_extension (MatroidType &matroid, MatrixType &matrix, NestedMinorSequenceType &nested_minors, RowThreeConnectivity &row_three_connectivity, ColumnThreeConnectivity &column_three_connectivity)
 
template<typename MatroidType , typename MatrixType , typename NestedMinorSequenceType , typename RowThreeConnectivity , typename ColumnThreeConnectivity >
char find_parallel_or_unit_vector (MatroidType &matroid, MatrixType &matrix, NestedMinorSequenceType &nested_minors, RowThreeConnectivity &row_three_connectivity, ColumnThreeConnectivity &column_three_connectivity, size_t &index)
 
template<typename MatroidType , typename MatrixType , typename NestedMinorSequenceType , typename RowThreeConnectivity , typename ColumnThreeConnectivity >
separation find_elaborate_extension (MatroidType &matroid, MatrixType &matrix, NestedMinorSequenceType &nested_minors, RowThreeConnectivity &row_three_connectivity, ColumnThreeConnectivity &column_three_connectivity, size_t index, matroid_element_set &extra_elements)
 
template<typename MatroidType , typename MatrixType >
separation find_minor_sequence (MatroidType &matroid, MatrixType &matrix, nested_minor_sequence &nested_minors, matroid_element_set &extra_elements, logger &log)
 
template<typename MatroidType , typename MatrixType >
separation find_wheel_minor (matroid_permuted< MatroidType > &permuted_matroid, matrix_permuted< MatrixType > &permuted_matrix, matroid_element_set &extra_elements)
 
template<typename T >
gcd_impl (T a, T b, T &s, T &t)
 
template<typename T >
gcd (T a, T b, T &s, T &t)
 
template<typename T >
gcd (T a, T b)
 
bool ghouila_houri_is_totally_unimodular_enum_rows (const integer_matrix &matrix)
 
bool ghouila_houri_is_totally_unimodular_enum_columns (const integer_matrix &matrix)
 
bool ghouila_houri_is_totally_unimodular (const integer_matrix &matrix)
 
template<typename MatroidType , typename MatrixType >
int find_parallel_to_row (const MatroidType &matroid, const MatrixType &matrix, const size_t minor_height, const size_t minor_width, const size_t row)
 
template<typename Graph , typename Vertex , typename EdgeSet >
struct articulation_edge_filter< Graph > make_articulation_edge_filter (const Graph *graph, const Vertex *articulation_vertex, const EdgeSet *evil_edges)
 
template<typename MatroidType , typename MatrixType >
bool extend_graph (matroid_graph &graph, const MatroidType &matroid, const MatrixType &matrix, const size_t minor_height, const size_t minor_width, const nested_minor_sequence::extension_type extension_type)
 
template<typename MatroidType , typename MatrixType , typename NestedMinorSequenceType >
matroid_graphconstruct_matroid_graph (const MatroidType &matroid, const MatrixType &matrix, const NestedMinorSequenceType &nested_minors, size_t &largest_graphic_minor)
 
template<typename Matrix >
void matrix_row_combine (Matrix &matrix, size_t row1, size_t row2, long long ul, long long ur, long long ll, long long lr)
 
template<typename Matrix >
bool matrix_row_gcd (Matrix &matrix, size_t row1, size_t row2, size_t column)
 
template<typename Matrix >
bool matrix_column_gcd (Matrix &matrix, size_t column1, size_t column2, size_t row)
 
template<typename InputMatrix , typename OutputMatrix >
size_t matrix_find_column_basis_and_transform_integral (const InputMatrix &input_matrix, OutputMatrix &output_matrix, std::vector< size_t > &column_basis)
 
std::ostream & operator<< (std::ostream &stream, logger &log)
 
template<typename M >
TransposedMatrix< M > transposedMatrix (M &matrix)
 Returns a TransposedMatrix for matrix. More...
 
template<typename Matrix >
bool find_smallest_nonzero_matrix_entry (const Matrix &matrix, size_t row_first, size_t row_beyond, size_t column_first, size_t column_beyond, size_t &row, size_t &column)
 
template<typename Matrix >
bool matrix_row_zero (const Matrix &matrix, size_t row, size_t column_first, size_t column_beyond)
 
template<typename Matrix >
bool matrix_column_zero (const Matrix &matrix, size_t column, size_t row_first, size_t row_beyond)
 
template<typename Matrix1 , typename Matrix2 >
bool equals (const Matrix1 &first, const Matrix2 &second)
 
template<typename MatrixType >
void matrix_permute1 (MatrixType &matrix, size_t index1, size_t index2)
 
template<typename MatrixType >
void matrix_permute2 (MatrixType &matrix, size_t index1, size_t index2)
 
template<typename MatrixType >
void matrix_set_value (matrix_permuted< MatrixType > &matrix, size_t row, size_t column, typename MatrixType::value_type value)
 
template<typename MatrixType >
void matrix_set_value (matrix_permuted< const MatrixType > &matrix, size_t row, size_t column, typename MatrixType::value_type value)
 
template<typename MatrixType >
void matrix_permute1 (matrix_permuted< MatrixType > &matrix, size_t index1, size_t index2)
 
template<typename MatrixType >
void matrix_permute2 (matrix_permuted< MatrixType > &matrix, size_t index1, size_t index2)
 
template<typename MatrixType >
void matrix_binary_pivot (matrix_permuted< MatrixType > &matrix, size_t i, size_t j)
 
template<typename MatrixType >
permutation matrix_get_perm1 (matrix_permuted< MatrixType > &matrix)
 
template<typename MatrixType >
permutation matrix_get_perm2 (matrix_permuted< MatrixType > &matrix)
 
template<typename MatrixType >
permutation matrix_get_perm1 (matrix_transposed< MatrixType > &matrix)
 
template<typename MatrixType >
permutation matrix_get_perm2 (matrix_transposed< MatrixType > &matrix)
 
template<typename MatrixType >
void matrix_set_perm1 (matrix_permuted< MatrixType > &matrix, const permutation &permutation)
 
template<typename MatrixType >
void matrix_set_perm2 (matrix_permuted< MatrixType > &matrix, const permutation &permutation)
 
template<typename MatrixType >
void matrix_set_perm1 (matrix_transposed< MatrixType > &matrix, const permutation &permutation)
 
template<typename MatrixType >
void matrix_set_perm2 (matrix_transposed< MatrixType > &matrix, const permutation &permutation)
 
template<typename MatrixType >
void matrix_apply_row_permutation (MatrixType &matrix, const permutation &perm)
 
template<typename MatrixType >
void matrix_apply_column_permutation (MatrixType &matrix, const permutation &perm)
 
template<typename MatrixType , typename ElementLess >
void matrix_reorder_rows (const MatrixType &matrix, size_t row_first, size_t row_beyond, size_t column_first, size_t column_beyond, ElementLess element_less, permutation &result_permutation)
 
template<typename MatrixType , typename ElementLess >
void matrix_reorder_columns (const MatrixType &matrix, size_t row_first, size_t row_beyond, size_t column_first, size_t column_beyond, ElementLess element_less, permutation &result_permutation)
 
template<typename MatrixType , typename ElementLess >
void matrix_reorder_rows (MatrixType &matrix, size_t row_first, size_t row_beyond, size_t column_first, size_t column_beyond, ElementLess element_less)
 
template<typename MatrixType , typename ElementLess >
void matrix_reorder_columns (MatrixType &matrix, size_t row_first, size_t row_beyond, size_t column_first, size_t column_beyond, ElementLess element_less)
 
template<typename MatrixType >
void matrix_binary_pivot (MatrixType &matrix, size_t i, size_t j)
 
template<typename MatrixType >
matrix_transposed< MatrixType > make_transposed_matrix (MatrixType &matrix)
 
template<typename MatrixType >
void matrix_set_value (matrix_transposed< MatrixType > &matrix, size_t row, size_t column, typename MatrixType::value_type value)
 
template<typename MatrixType >
void matrix_set_value (matrix_transposed< const MatrixType > &matrix, size_t row, size_t column, typename MatrixType::value_type value)
 
template<typename MatrixType >
void matrix_binary_pivot (matrix_transposed< MatrixType > &matrix, size_t i, size_t j)
 
template<typename MatrixType >
void matrix_permute1 (matrix_transposed< MatrixType > &matrix, size_t index1, size_t index2)
 
template<typename MatrixType >
void matrix_permute2 (matrix_transposed< MatrixType > &matrix, size_t index1, size_t index2)
 
template<typename NameType >
void matroid_permute1 (matroid< NameType > &matroid, size_t index1, size_t index2)
 
template<typename MatroidType , typename MatrixType >
void matroid_permute1 (MatroidType &matroid, MatrixType &matrix, size_t index1, size_t index2)
 
template<typename NameType >
void matroid_permute2 (matroid< NameType > &matroid, size_t index1, size_t index2)
 
template<typename MatroidType , typename MatrixType >
void matroid_permute2 (MatroidType &matroid, MatrixType &matrix, size_t index1, size_t index2)
 
template<typename NameType >
void matroid_binary_pivot (matroid< NameType > &matroid, size_t i, size_t j)
 
template<typename MatroidType , typename MatrixType >
void matroid_binary_pivot (MatroidType &matroid, MatrixType &matrix, size_t i, size_t j)
 
template<typename MatroidType , typename MatrixType >
void matroid_print (const MatroidType &matroid, const MatrixType &matrix)
 
template<typename MatroidType , typename MatrixType >
void matroid_print_minor (const MatroidType &matroid, const MatrixType &matrix, size_t height, size_t width)
 
template<typename MatroidType >
std::set< int > matroid_elements (const MatroidType &matroid)
 
std::ostream & operator<< (std::ostream &stream, const tu::matroid_graph &graph)
 
template<typename MatroidType >
void matroid_permute1 (matroid_permuted< MatroidType > &matroid, size_t index1, size_t index2)
 
template<typename MatroidType >
void matroid_permute2 (matroid_permuted< MatroidType > &matroid, size_t index1, size_t index2)
 
template<typename MatroidType >
void matroid_binary_pivot (matroid_permuted< MatroidType > &matroid, size_t i, size_t j)
 
template<typename MatroidType >
void matroid_set_perm1 (matroid_permuted< MatroidType > &matroid, const permutation &permutation)
 
template<typename MatroidType >
void matroid_set_perm2 (matroid_permuted< MatroidType > &matroid, const permutation &permutation)
 
template<typename MatroidType >
void matroid_set_perm1 (matroid_transposed< MatroidType > &matroid, const permutation &permutation)
 
template<typename MatroidType >
void matroid_set_perm2 (matroid_transposed< MatroidType > &matroid, const permutation &permutation)
 
template<typename MatroidType , typename MatrixType >
void matroid_apply_row_permutation (MatroidType &matroid, MatrixType &matrix, const permutation &perm)
 
template<typename MatroidType , typename MatrixType >
void matroid_apply_column_permutation (MatroidType &matroid, MatrixType &matrix, const permutation &perm)
 
template<typename MatroidType , typename MatrixType , typename ElementLess >
void matroid_reorder_rows (MatroidType &matroid, MatrixType &matrix, size_t row_first, size_t row_beyond, size_t column_first, size_t column_beyond, ElementLess element_less)
 
template<typename MatroidType , typename MatrixType , typename ElementLess >
void matroid_reorder_columns (MatroidType &matroid, MatrixType &matrix, size_t row_first, size_t row_beyond, size_t column_first, size_t column_beyond, ElementLess element_less)
 
template<typename MatroidType , typename MatrixType >
void matroid_separate (MatroidType &matroid, MatrixType &matrix, const int type, const std::pair< size_t, size_t > &split, integer_matroid &upper_left_matroid, integer_matrix &upper_left_matrix, integer_matroid &lower_right_matroid, integer_matrix &lower_right_matrix)
 
template<typename MatroidType >
matroid_transposed< MatroidType > make_transposed_matroid (MatroidType &matroid)
 
template<typename MatroidType >
void matroid_permute1 (matroid_transposed< MatroidType > &matroid, size_t index1, size_t index2)
 
template<typename MatroidType >
void matroid_permute2 (matroid_transposed< MatroidType > &matroid, size_t index1, size_t index2)
 
template<typename MatroidType >
void matroid_binary_pivot (matroid_transposed< MatroidType > &matroid, size_t i, size_t j)
 
template<typename NestedMinorSequenceType >
nested_minor_sequence_transposed< NestedMinorSequenceType > make_transposed_nested_minor_sequence (NestedMinorSequenceType &sequence)
 
template<typename MatrixType >
rank_distribution partition (matrix_permuted< const MatrixType > &matrix, size_t &top_left_height, size_t &top_left_width, size_t &bottom_right_height, size_t &bottom_right_width)
 
bool operator== (const permutation &p, const permutation &q)
 
std::ostream & operator<< (std::ostream &stream, const permutation &p)
 
template<class Less >
void sort (permutation &permutation, size_t first, size_t beyond, Less &less)
 
template<class Less >
void sort (permutation &permutation, Less &less)
 
template<typename MatrixType >
bool is_r10 (MatrixType matrix)
 
template<typename MatrixType >
bool find_nonzero_column (MatrixType &matrix, size_t column_first, size_t column_beyond, size_t row_first, size_t row_beyond, size_t target_column)
 
template<typename MatrixType >
void check_sign (const MatrixType &matrix, const std::vector< bipartite_graph_bfs_node > &spanning_tree, const bipartite_graph_dimensions &dim, const std::set< size_t > &nodes, size_t current_index, size_t column, std::map< size_t, bool > &changes)
 
template<typename M >
bool sign_matrix (M &matrix, submatrix_indices *violator)
 
template<typename Matrix >
void smith_normal_form_diagonal (Matrix &input_matrix, std::vector< int > &diagonal)
 
decomposed_matroiddecompose_binary_matroid (const integer_matrix &matrix, log_level level)
 
bool is_totally_unimodular (const integer_matrix &matrix, log_level level)
 
bool is_totally_unimodular (const integer_matrix &matrix, decomposed_matroid *&decomposition, log_level level)
 
bool is_totally_unimodular (const integer_matrix &matrix, decomposed_matroid *&decomposition, submatrix_indices &violator, log_level level)
 
bool is_totally_unimodular (const integer_matrix &matrix, submatrix_indices &violator, log_level level)
 
bool is_signed_matrix (const integer_matrix &matrix)
 
bool is_signed_matrix (const integer_matrix &matrix, submatrix_indices &violator)
 
bool sign_matrix (integer_matrix &matrix)
 
void support_matrix (integer_matrix &matrix)
 
CMR_EXPORT bool is_zero_plus_minus_one_matrix (const integer_matrix &matrix)
 
CMR_EXPORT bool is_zero_plus_minus_one_matrix (const integer_matrix &matrix, std::pair< integer_matrix::size_type, integer_matrix::size_type > &position)
 
CMR_EXPORT bool is_zero_one_matrix (const integer_matrix &matrix)
 
CMR_EXPORT bool is_zero_one_matrix (const integer_matrix &matrix, std::pair< integer_matrix::size_type, integer_matrix::size_type > &position)
 
template<typename Matrix >
bool test_k_modularity (const Matrix &matrix, size_t &rank, unsigned int *pk, bool enforce_unimodularity, log_level level)
 
bool is_unimodular (const integer_matrix &matrix, size_t &rank, log_level level)
 
bool is_strongly_unimodular (const integer_matrix &matrix, size_t &rank, log_level level)
 
bool is_k_modular (const integer_matrix &matrix, size_t &rank, log_level level)
 
bool is_k_modular (const integer_matrix &matrix, size_t &rank, unsigned int &k, log_level level)
 
bool is_strongly_k_modular (const integer_matrix &matrix, size_t &rank, log_level level)
 
bool is_strongly_k_modular (const integer_matrix &matrix, size_t &rank, unsigned int &k, log_level level)
 
unsigned int get_k_modular_integrality (const integer_matrix &matrix, const integer_matrix &rhs)
 
bool is_k_modular_integral (const integer_matrix &matrix, const integer_matrix &rhs)
 
bool is_complement_total_unimodular (const integer_matrix &matrix, std::size_t &complementedRow, std::size_t &complementedColumn, log_level level)
 

Detailed Description

Copyright Matthias Walter 2010. Distributed under the Boost Software License, Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)

Typedef Documentation

typedef boost::property_map<matroid_graph, edge_matroid_element_t>::const_type tu::const_matroid_element_map
typedef boost::numeric::ublas::matrix<long long> tu::integer_matrix

Integer matrix

A matroid with integer names.

typedef boost::numeric::ublas::matrix_indirect<const integer_matrix, submatrix_indices::indirect_array_type> tu::integer_submatrix

Indirect integer matrix

typedef boost::property_map<matroid_graph, edge_matroid_element_t>::type tu::matroid_element_map
typedef std::set<int> tu::matroid_element_set

A set of matroid elements. Here, the original row elements are numbered -1, ..., -h and the original column elements +1, ..., +w

typedef boost::adjacency_list<boost::vecS, boost::vecS, boost::undirectedS, boost::no_property, matroid_element_property> tu::matroid_graph
typedef std::pair<size_t, size_t> tu::size_pair_t

Enumeration Type Documentation

Matroid element property

Enumerator
edge_matroid_element 
Enumerator
LOG_QUIET 
LOG_PROGRESSIVE 
LOG_VERBOSE 

Types of distributions of the rank 2

Enumerator
RANK_TOO_HIGH 
RANK_BL_1_TR_1 
RANK_BL_2_TR_0 
RANK_BL_0_TR_2 

Function Documentation

template<typename MatrixType , typename StartNodesContainerType , typename EndNodesContainerType >
bool tu::bipartite_graph_bfs ( const MatrixType &  matrix,
const bipartite_graph_dimensions dimensions,
const StartNodesContainerType &  start_nodes,
const EndNodesContainerType &  end_nodes,
bool  reach_all,
std::vector< bipartite_graph_bfs_node > &  result 
)
inline

Does a breadth-first-search on the bipartite graph defined by a submatrix of a given matrix. Running time: O(height * width)

Parameters
matrixGiven matrix
dimensionsObject for row/column-to-index calculations
start_nodesSet of starting nodes for search
end_nodesSet of nodes to be found
reach_allWhether we'd like to find all end nodes, or just one
resultWhether all resp. one end node was found.
Returns
true iff enough of the target nodes were found.

Initialize datastructure

Start bfs

row -> columns

column -> rows

template<typename MatrixType >
void tu::check_sign ( const MatrixType &  matrix,
const std::vector< bipartite_graph_bfs_node > &  spanning_tree,
const bipartite_graph_dimensions dim,
const std::set< size_t > &  nodes,
size_t  current_index,
size_t  column,
std::map< size_t, bool > &  changes 
)

Takes a spanning tree in the bipartite graph of a matrix and a set of nodes to be signed.

Parameters
matrixThe given matrix
spanning_treeA spanning tree
dimIndex-mapping for bipartite graph
nodesSet of nodes
current_indexCurrent index in the spanning tree
columnThe column to be signed
changesNecessary changes to the entries in the column

Root does never change.

Search for ancestors until reaching one of the given nodes.

If the ancestor is not yet processed, we recurse.

If sum (modulo 4) is not 0, we'd like to change the current one

template<typename MatroidType , typename MatrixType , typename NestedMinorSequenceType >
matroid_graph* tu::construct_matroid_graph ( const MatroidType &  matroid,
const MatrixType &  matrix,
const NestedMinorSequenceType &  nested_minors,
size_t &  largest_graphic_minor 
)

Either constructs a graph whose forest matroid is the given matroid or detects that the given matroid is non-graphic.

Parameters
matroidThe given matroid
matrixRepresentation matrix of the given matroid
nested_minorsSequence of nested minors
Returns
The constructed graph or NULL if the matroid is not graphic.

Initialize W3

template<typename MatroidType , typename MatrixType >
matroid_graph* tu::construct_small_matroid_graph ( MatroidType &  matroid,
MatrixType &  matrix 
)

Create the graph for 2 x w or h x 2 matrices.

Parameters
matroid
matrix
Returns
Created graph

0 1

1 0

1 1

template<typename MatroidType , typename MatrixType >
std::pair< bool, decomposed_matroid * > tu::decompose_binary_matroid ( MatroidType &  matroid,
MatrixType &  matrix,
matroid_element_set  extra_elements,
bool  construct_decomposition,
logger log 
)

Forward declaration.

Decomposes a given binary matroid.

Parameters
matroid
matrix
construct_decomposition
Returns

Identifies a W_3 minor in the upper left corner or finds a separation

Separation case

Filter or copy extra_elements depending on type of separation

Identifies a W_3 minor in the upper left corner or finds a separation

Separation case

Filter or copy extra_elements depending on type of separation

CMR_EXPORT decomposed_matroid * tu::decompose_binary_matroid ( const integer_matrix matrix,
log_level  level 
)

Returns a decomposition of a given binary matroid into a k-sum-decomposition (k=1,2,3) in graphic, cographic, R10 and maybe irregular components.

Parameters
matrixRepresentation matrix of a binary matroid
levelLog level
Returns
Root of decomposition tree
template<typename MatroidType , typename MatrixType >
std::pair<bool, decomposed_matroid*> tu::decompose_with_minor_sequence ( matroid_permuted< MatroidType > &  permuted_matroid,
matrix_permuted< MatrixType > &  permuted_matrix,
nested_minor_sequence nested_minors,
matroid_element_set  extra_elements,
bool  construct_decomposition,
logger log 
)

Filter or copy extra_elements depending on type of separation

CMR_EXPORT bool tu::determinant_is_totally_unimodular ( const integer_matrix matrix)

Checks all subdeterminants to test a given matrix for total unimodularity.

Parameters
matrixThe given matrix
Returns
true if and only if this matrix is totally unimdular
CMR_EXPORT bool tu::determinant_is_totally_unimodular ( const integer_matrix matrix,
submatrix_indices violator 
)

Checks all subdeterminants to test a given matrix for total unimodularity. If this is not the case, violator describes a violating submatrix.

Parameters
matrixThe given matrix
violatorThe violating submatrix, if the result is false
Returns
true if and only if the this matrix is totally unimodular

Examine the determinant

template<typename MatroidType , typename MatrixType , typename NestedMinorSequence >
separation tu::enumerate_separations ( MatroidType &  matroid,
MatrixType &  matrix,
const NestedMinorSequence &  nested_minors,
matroid_element_set extra_elements,
logger log 
)
inline

Enumerates the partitions along the sequence of nested minors.

Parameters
matroidGiven matroid
matrixRepresentation matrix of the given matroid
nested_minorsSequence of nested minors
extra_elementsSet of matroid-elements to be filled with pivot-elements
logLogger
Returns
true if and only if the partition algorithm was successful

Every regular 3-connected matroid which is non-graphic, non-cographic and not isomorphic to R10 must contain R12

The first minors must be enumerated fully until we have one with at least 7 elements.

preparations

Calculate number of enumerations

Full enumeration

Transform the mappings into row/column permutations.

Clever enumeration along the sequence

template<typename Matrix1 , typename Matrix2 >
bool tu::equals ( const Matrix1 &  first,
const Matrix2 &  second 
)
template<typename MatroidType , typename MatrixType >
bool tu::extend_graph ( matroid_graph graph,
const MatroidType &  matroid,
const MatrixType &  matrix,
const size_t  minor_height,
const size_t  minor_width,
const nested_minor_sequence::extension_type  extension_type 
)

Tries to extend a 3-connected graph of a graphic minor or detects that the extended minor is not graphic.

Parameters
graphThe graph of the graphic minor
matroidThe matroid
matrixRepresentation matrix of the matroid
minor_heightNumber of row elements of the minor
minor_widthNumber of column elements of the minor
extension_typeType of the extension
Returns
true if and only if the extension was possible.

Distinct extension cases

Find vertices of corresponding edges

Remove old edges

Create new vertices

Create new edges

Find vertices of corresponding edges

Create new vertex

Create new edges

Find vertices of corresponding edges

Remove old edges

Create new vertex

Create new edges

Check if edges form a path

Collect all 1-edges

Create the new vertex and connect it with the_vertex.

Reconnect the 1-edges

Count for each vertex the number of paths that use it

Find articulation points for the graph without 1-edges

Filter the chosen articulation point and the one-edges and check the remaining graph

We should really have articulation point + 2 further components

There cannot be a 1-edge inside one component.

template<typename MatroidType , typename MatrixType , typename NestedMinorSequenceType , typename RowThreeConnectivity , typename ColumnThreeConnectivity >
separation tu::find_elaborate_extension ( MatroidType &  matroid,
MatrixType &  matrix,
NestedMinorSequenceType &  nested_minors,
RowThreeConnectivity &  row_three_connectivity,
ColumnThreeConnectivity &  column_three_connectivity,
size_t  index,
matroid_element_set extra_elements 
)

Searches for elaborate extensions or a 2-separation by BFS in the bipartite graph, starting at a parallel row or unit column.

Parameters
matroidA given matroid
matrixRepresentation matrix for the given matroid
nested_minorsSequence of nested minors
row_three_connectivityConnectivity of the rows
column_three_connectivityConnectivity of the columns
indexIndex of the row in the minor
extra_elementsSet of matroid-elements to be filled with pivot elements.
Returns
Either a 2-separation or no separation

Initialize BFS

Sets row types.

Sets column types.

Call BFS on modified matrix.

If a path it found, we can achieve one of the 3 elaborate extensions via path shortening technique.

Do pivots

1st case

2nd case

This implies a 2-separation

template<typename MatroidType , typename MatrixType >
separation tu::find_minor_sequence ( MatroidType &  matroid,
MatrixType &  matrix,
nested_minor_sequence nested_minors,
matroid_element_set extra_elements,
logger log 
)

Constructs a either a sequence of nested minors, either up to the whole matroid (which is then 3-connected) or until we find a 1- or 2-separation. At first, we search for simple row / column extensions. If there are no, we try elaborate extensions.

Parameters
matroidA given matroid
matrixRepresentation matrix for the given matroid
nested_minorsSequence of nested minors
extra_elementsSet of matroid-elements to be filled with pivot elements
logLogger
Returns
Either a 2-separation or no separation, in case the matroid is 3-connected.

To search for parallel-columns we search for parallel rows in the transpose.

Simple row extension

Simple column extension

Elaborate extension with starting row

Elaborate extension with starting column

template<typename MatrixType >
bool tu::find_nonzero_column ( MatrixType &  matrix,
size_t  column_first,
size_t  column_beyond,
size_t  row_first,
size_t  row_beyond,
size_t  target_column 
)

Generic function to find a non-zero column and swap it to a given position.

Parameters
matrixThe given matrix
column_firstFirst index of a column range
column_beyondBeyond index of a column range
row_firstFirst index of a row range
row_beyondBeyond index of a row range
target_columnColumn to be swapped to
Returns
Whether a non-zero column was found.
template<typename MatroidType , typename MatrixType , typename NestedMinorSequenceType , typename RowThreeConnectivity , typename ColumnThreeConnectivity >
char tu::find_parallel_or_unit_vector ( MatroidType &  matroid,
MatrixType &  matrix,
NestedMinorSequenceType &  nested_minors,
RowThreeConnectivity &  row_three_connectivity,
ColumnThreeConnectivity &  column_three_connectivity,
size_t &  index 
)

Takes a represented matroid and a sequence of nested minors and searches for a row (resp. column) below (right to) the last minor which is either parallel to one of the minor's rows (columns) or a unit vector.

Parameters
matroidA given matroid
matrixRepresentation matrix for the given matroid
nested_minorsSequence of nested minors
row_three_connectivityConnectivities of all possible rows
column_three_connectivityConnectivities of all possible columns
indexReturns the index of the parallel row/column or the index of the unique one.
Returns
true if and only if something is found
template<typename MatroidType , typename MatrixType >
int tu::find_parallel_to_row ( const MatroidType &  matroid,
const MatrixType &  matrix,
const size_t  minor_height,
const size_t  minor_width,
const size_t  row 
)

Finds a matroid element which is parallel to a row element in minor of a given matroid or a unit vector in a column element in the minor.

Parameters
matroidThe given matroid
matrixRepresentation matrix of the given matroid
minor_heightNumber of rows of the minor
minor_widthNumber of columns of the minor
rowThe element must be parallel to this row (along the minor)
Returns
The matroid element that was found to be parallel / unit vector

Look for unit vector

Check for parallel

template<typename MatroidType , typename MatrixType , typename NestedMinorSequenceType , typename RowThreeConnectivity , typename ColumnThreeConnectivity >
bool tu::find_simple_row_extension ( MatroidType &  matroid,
MatrixType &  matrix,
NestedMinorSequenceType &  nested_minors,
RowThreeConnectivity &  row_three_connectivity,
ColumnThreeConnectivity &  column_three_connectivity 
)

Searches for an extension of a sequence of nested minors of a given matroid, consisting of a non-parallel, non-unit-vector and non-zero row.

Parameters
matroidA given matroid
matrixRepresenation matrix for the given matroid
nested_minorsSequence of nested minors
row_three_connectivityConnectivities of all possible rows
column_three_connectivityConnectivities of all possible columns
Returns
true if and only if the extension was found
template<typename Matrix >
bool tu::find_smallest_nonzero_matrix_entry ( const Matrix matrix,
size_t  row_first,
size_t  row_beyond,
size_t  column_first,
size_t  column_beyond,
size_t &  row,
size_t &  column 
)
template<typename MatroidType , typename MatrixType >
separation tu::find_wheel_minor ( matroid_permuted< MatroidType > &  permuted_matroid,
matrix_permuted< MatrixType > &  permuted_matrix,
matroid_element_set extra_elements 
)

Searches for a W3 minor of a given matroid. It might find a 1- or 2-separation instead. The matroid bust be at least 3 x 3.

Parameters
permuted_matroidGiven matroid
permuted_matrixRepresentation matrix for the given matroid
extra_elementsSet of matroid-elements to be filled with pivot elements
Returns
Separation, if one was found

Permute 1's in first row to the left.

Trivial 1-separation

1- or 2-separation

Ensure we have a 2x2 block of ones

Grow the block to a set-maximal one.

Search for a path in BFS graph

There must be a 2-separation.

Swap unreachable rows to top

Swap unreachable columns to left

Go back along the path to find the nearest endpoint.

Find indices for basic W3-parts.

Apply path shortening technique.

Collect more W3-indices.

Some normalization

template<typename T >
T tu::gcd ( a,
b,
T &  s,
T &  t 
)
template<typename T >
T tu::gcd ( a,
b 
)
template<typename T >
T tu::gcd_impl ( a,
b,
T &  s,
T &  t 
)
CMR_EXPORT unsigned int tu::get_k_modular_integrality ( const integer_matrix matrix,
const integer_matrix rhs 
)

In case a matrix A is k-modular, it may lead to q-integrality of the polyhedron A*x = b, x >= 0 if B*x = b*q for some basis B of A. This method finds a minimal integer q, assuming that A is k-modular. If q = 1, the polyhedron is integral, if q = 2 it is half-integral, etc.

Parameters
matrixThe matrix A, being k-modular.
rhsThe rhs, given as a m x 1 matrix.
Returns
Minimal q as defined above
CMR_EXPORT bool tu::ghouila_houri_is_totally_unimodular ( const integer_matrix matrix)

Tests a given matrix to be totally unimodular using ghouila-houri's criterion by enumeration of either row or column subsets by choosing the one which induces fewer enumerations.

Parameters
matrixThe given matrix
Returns
true if and only if this matrix is totally unimodular
CMR_EXPORT bool tu::ghouila_houri_is_totally_unimodular_enum_columns ( const integer_matrix matrix)

Tests a given matrix to be totally unimodular using ghouila-houri's criterion by enumeration of column subsets.

Parameters
matrixThe given matrix
Returns
true if and only if this matrix is totally unimodular
CMR_EXPORT bool tu::ghouila_houri_is_totally_unimodular_enum_rows ( const integer_matrix matrix)

Tests a given matrix to be totally unimodular using ghouila-houri's criterion by enumeration of row subsets.

Parameters
matrixThe given matrix
Returns
true if and only if this matrix is totally unimodular
CMR_EXPORT bool tu::is_complement_total_unimodular ( const integer_matrix matrix,
std::size_t &  complementedRow,
std::size_t &  complementedColumn,
log_level  level = LOG_QUIET 
)

Tests if a matrix A is complement totally unimodular (ctu).

Parameters
matrixThe matrix A.
complementedRowIf A is not ctu, indicates the complemented row; #rows(A) if no row was complemented.
complementedColumnIf A is not ctu, indicates the complemented column; #columns(A) if no column was complemented.
levelLog level
Returns
true if and only if the matrix is ctu.
CMR_EXPORT bool tu::is_k_modular ( const integer_matrix matrix,
size_t &  rank,
log_level  level 
)

Tests for k-modularity without certificates. A matrix of rank r is k-modular if and only if for every column basis B the gcd (greatest common divisor) of the determinants of regular r x r submatrices of B has the same absolute value.

Parameters
matrixThe matrix to be tested
rankReturns the rank r
levelLog level
Returns
true if and only if the matrix is k-modular
CMR_EXPORT bool tu::is_k_modular ( const integer_matrix matrix,
size_t &  rank,
unsigned int &  k,
log_level  level 
)

Tests for k-modularity without certificates. A matrix of rank r is k-modular if and only if for every column basis B the gcd (greatest common divisor) of the determinants of regular r x r submatrices of B has the same absolute value k. It also computes k.

Parameters
matrixThe matrix to be tested
rankReturns the rank r
kThe common absolute value
levelLog level
Returns
true if and only if the matrix is k-modular
CMR_EXPORT bool tu::is_k_modular_integral ( const integer_matrix matrix,
const integer_matrix rhs 
)

In case a matrix A is k-modular, it may lead to integrality of the polyhedron A*x = b, x >= 0 if B*x = b for some basis B of A. This method tests of that property.

Parameters
matrixThe matrix A, having the Dantzig property.
rhsThe rhs, given as a m x 1 matrix.
Returns
true if and only if the polyhedron is integral
template<typename MatrixType >
bool tu::is_r10 ( MatrixType  matrix)
inline

Tests a given matrix to be matrix-isomorphic to one of the R10-representing matrices by examining the corresponding bipartite graphs.

Parameters
matrixA given matrix
Returns
true if and only if this matrix is a represenation matrix for R10
CMR_EXPORT bool tu::is_signed_matrix ( const integer_matrix matrix)

Tests if a given matrix is a signed version of its support matrix already. Running time: O(height * width * min(height, width))

Parameters
matrixThe matrix to be tested
Returns
true if and only if the support matrix can be signed to the orignal
CMR_EXPORT bool tu::is_signed_matrix ( const integer_matrix matrix,
submatrix_indices violator 
)

Tests if a given matrix is a signed version of its support matrix already, returning the indices of a violating submatrix if this is not the case. Running time: O(height * width * min(height, width))

Parameters
matrixThe matrix to be tested
violatorReturns violator indices
Returns
true if and only if the support matrix can be signed to the original
CMR_EXPORT bool tu::is_strongly_k_modular ( const integer_matrix matrix,
size_t &  rank,
log_level  level 
)

Tests for strong k-modularity without certificates. A matrix of rank r is k-modular if and only if it and its transpose are k-modular.

Parameters
matrixThe matrix to be tested
rankReturns the rank r
levelLog level
Returns
true if and only if the matrix is k-modular
CMR_EXPORT bool tu::is_strongly_k_modular ( const integer_matrix matrix,
size_t &  rank,
unsigned int &  k,
log_level  level 
)

Tests for strong k-modularity without certificates. A matrix of rank r is k-modular if and only if it and its transpose are k-modular. It also computes k.

Parameters
matrixThe matrix to be tested
rankReturns the rank r
kThe common absolute value
levelLog level
Returns
true if and only if the matrix is k-modular
CMR_EXPORT bool tu::is_strongly_unimodular ( const integer_matrix matrix,
size_t &  rank,
log_level  level 
)

Tests for strong unimodularity without certificates. A matrix is strongly unimodular if and only if it and its transpose are unimodular.

Parameters
matrixThe matrix to be tested
rankReturns the rank r of the matrix
levelLog level
Returns
true if and only if the matrix is strongly unimodular
CMR_EXPORT bool tu::is_totally_unimodular ( const integer_matrix matrix,
log_level  level 
)

Tests for total unimodularity without certificates. A matrix is totally unimodular if and only if every square submatrix has a determinant of -1, 0 or +1.

Parameters
matrixThe matrix to be tested
levelLog level
Returns
true if and only if the matrix is totally unimodular

Test for being a -1,0,+1 matrix

Signing test

Decomposition of matroid represented by support matrix

CMR_EXPORT bool tu::is_totally_unimodular ( const integer_matrix matrix,
decomposed_matroid *&  decomposition,
log_level  level 
)

Tests for total unimodularity with a positive certificate. A matrix is totally unimodular if and only if every square submatrix has a determinant of -1, 0 or +1. If the matrix has this property, the routine returns a decomposition of the underlying binary matroid into a k-sum-decomposition (k=1,2,3) in graphic, cographic, R10 and maybe irregular components.

Parameters
matrixThe matrix to be tested
decompositionReturns the root of the decomposition tree
levelLog level
Returns
true if and only if the matrix is totally unimodular

Test for being a -1,0,+1 matrix

Signing test

Decomposition of matroid represented by support matrix

CMR_EXPORT bool tu::is_totally_unimodular ( const integer_matrix matrix,
decomposed_matroid *&  decomposition,
submatrix_indices violator,
log_level  level 
)

Tests for total unimodularity with certificates. A matrix is totally unimodular if and only if every square submatrix has a determinant of -1, 0 or +1. If the matrix has this property, the routine returns a decomposition of the underlying binary matroid into a k-sum-decomposition (k=1,2,3) in graphic, cographic, R10 and maybe irregular components. If the matrix does not have this property,the routine returns the indices of a violating submatrix.

Parameters
matrixThe matrix to be tested
decompositionReturns the root of the decomposition tree.
violatorReturns violator indices
levelLog level
Returns
true if and only if the matrix is totally unimodular

Test for being a -1,0,+1 matrix

Signing test

Decomposition of matroid represented by support matrix

CMR_EXPORT bool tu::is_totally_unimodular ( const integer_matrix matrix,
submatrix_indices violator,
log_level  level 
)

Tests for total unimodularity with a negative certificate. A matrix is totally unimodular if and only if every square submatrix has a determinant of -1, 0 or +1. If the matrix does not have this property,the routine returns the indices of a violating submatrix.

Parameters
matrixThe matrix to be tested
violatorReturns violator indices
levelLog level
Returns
true if and only if the matrix is totally unimodular
CMR_EXPORT bool tu::is_unimodular ( const integer_matrix matrix,
size_t &  rank,
log_level  level 
)

Tests for unimodularity without certificates. A matrix of rank r is unimodular if and only if for every column basis B the gcd (greatest common divisor) of the determinants of regular r x r submatrices of B is 1.

Parameters
matrixThe matrix to be tested
rankReturns the rank k
levelLog level
Returns
true if and only if the matrix is unimodular
bool tu::is_zero_one_matrix ( const integer_matrix matrix)

Tests if the given matrix contains only 0 or 1 entries.

Parameters
matrixThe given matrix
Returns
true if and only if it is a 0-1 matrix
bool tu::is_zero_one_matrix ( const integer_matrix matrix,
std::pair< integer_matrix::size_type, integer_matrix::size_type > &  position 
)

Tests if the given matrix contains only 0 or 1 entries, returning the violating position if this is not the case.

Parameters
matrixThe given matrix
positionReturns a violating entry
Returns
true if and only if it is a 0-1 matrix
bool tu::is_zero_plus_minus_one_matrix ( const integer_matrix matrix)

Tests if the given matrix contains only -1,0,+1 entries.

Parameters
matrixThe given matrix
Returns
true if and only if it is a -1,0,+1 matrix
bool tu::is_zero_plus_minus_one_matrix ( const integer_matrix matrix,
std::pair< integer_matrix::size_type, integer_matrix::size_type > &  position 
)

Tests if the given matrix contains only -1,0,+1 entries, returning the violating position if this is not the case.

Parameters
matrixThe given matrix
positionReturns a violating entry
Returns
true if and only if it is a -1,0,+1 matrix
template<typename Graph , typename Vertex , typename EdgeSet >
struct articulation_edge_filter< Graph > tu::make_articulation_edge_filter ( const Graph *  graph,
const Vertex *  articulation_vertex,
const EdgeSet *  evil_edges 
)

Constructs an articulation edge filter for a given graph.

Parameters
graphThe given graph
articulation_vertexVertex to be excluded
evil_edgesEdges to be excluded
Returns
The constructed filter
template<typename MatrixType >
matrix_transposed<MatrixType> tu::make_transposed_matrix ( MatrixType &  matrix)
inline

Creates a transpose proxy of a given matrix.

Parameters
matrixThe original matrix
Returns
A matrix_transposed proxy
template<typename MatroidType >
matroid_transposed<MatroidType> tu::make_transposed_matroid ( MatroidType &  matroid)
inline

Creates a transpose proxy of a given matroid, i.e. the dual matroid.

Parameters
matroidThe original matroid
Returns
A matroid_transposed proxy
template<typename NestedMinorSequenceType >
nested_minor_sequence_transposed<NestedMinorSequenceType> tu::make_transposed_nested_minor_sequence ( NestedMinorSequenceType &  sequence)
inline

Creates a transpose proxy of a given nested minor sequence.

Parameters
sequenceGiven sequence
Returns
The transpose proxy
template<typename MatrixType >
void tu::matrix_apply_column_permutation ( MatrixType &  matrix,
const permutation perm 
)
inline

Applies a columnn permutation to a given matrix, i.e. performs the necessary column swaps at the matrix directly.

Parameters
matrixThe matrix to be changed.
permThe given permutation
template<typename MatrixType >
void tu::matrix_apply_row_permutation ( MatrixType &  matrix,
const permutation perm 
)
inline

Applies a row permutation to a given matrix, i.e. performs the necessary row swaps at the matrix directly.

Parameters
matrixThe matrix to be changed.
permThe given permutation
template<typename MatrixType >
void tu::matrix_binary_pivot ( MatrixType &  matrix,
size_t  i,
size_t  j 
)
template<typename MatrixType >
void tu::matrix_binary_pivot ( matrix_transposed< MatrixType > &  matrix,
size_t  i,
size_t  j 
)
inline

Free function to perform a pivot on a transposed matrix.

Parameters
matrixThe transposed matrix
iRow index
jColumn index
template<typename MatrixType >
void tu::matrix_binary_pivot ( matrix_permuted< MatrixType > &  matrix,
size_t  i,
size_t  j 
)
inline

Free function to perform a pivot on a permuted matrix.

Parameters
matrixThe permuted matrix
iRow index
jColumn index
template<typename Matrix >
bool tu::matrix_column_gcd ( Matrix matrix,
size_t  column1,
size_t  column2,
size_t  row 
)
template<typename Matrix >
bool tu::matrix_column_zero ( const Matrix matrix,
size_t  column,
size_t  row_first,
size_t  row_beyond 
)
template<typename InputMatrix , typename OutputMatrix >
size_t tu::matrix_find_column_basis_and_transform_integral ( const InputMatrix &  input_matrix,
OutputMatrix &  output_matrix,
std::vector< size_t > &  column_basis 
)

Swap a non-zero column to the pivot column.

Swap the row to the pivot row (in the output matrix. This is necessary later anyway.

Reduce the entries below the pivot which is at (rank, rank).

Divide row by gcd of all entries.

Reduce column above entry

template<typename MatrixType >
permutation tu::matrix_get_perm1 ( matrix_permuted< MatrixType > &  matrix)
inline
template<typename MatrixType >
permutation tu::matrix_get_perm1 ( matrix_transposed< MatrixType > &  matrix)
inline
template<typename MatrixType >
permutation tu::matrix_get_perm2 ( matrix_permuted< MatrixType > &  matrix)
inline
template<typename MatrixType >
permutation tu::matrix_get_perm2 ( matrix_transposed< MatrixType > &  matrix)
inline
template<typename MatrixType >
void tu::matrix_permute1 ( matrix_transposed< MatrixType > &  matrix,
size_t  index1,
size_t  index2 
)
inline

Free function to permute two rows of a transposed matrix.

Parameters
matrixThe transposed matrix
index1First index
index2Second index
template<typename MatrixType >
void tu::matrix_permute1 ( matrix_permuted< MatrixType > &  matrix,
size_t  index1,
size_t  index2 
)
inline

Free function to permute two rows of a permuted matrix.

Parameters
matrixThe permuted matrix
index1First index
index2Second index
template<typename MatrixType >
void tu::matrix_permute1 ( MatrixType &  matrix,
size_t  index1,
size_t  index2 
)
inline

Free function to permute two rows of a permuted matrix.

Parameters
matrixThe permuted matrix
index1First index
index2Second index
template<typename MatrixType >
void tu::matrix_permute2 ( matrix_transposed< MatrixType > &  matrix,
size_t  index1,
size_t  index2 
)
inline

Free function to permute two columns of a transposed matrix.

Parameters
matrixThe transposed matrix
index1First index
index2Second index
template<typename MatrixType >
void tu::matrix_permute2 ( matrix_permuted< MatrixType > &  matrix,
size_t  index1,
size_t  index2 
)
inline

Free function to permute two columns of a permuted matrix.

Parameters
matrixThe permuted matrix
index1First index
index2Second index
template<typename MatrixType >
void tu::matrix_permute2 ( MatrixType &  matrix,
size_t  index1,
size_t  index2 
)
inline

Free function to permute two columns of a permuted matrix.

Parameters
matrixThe permuted matrix
index1First index
index2Second index
template<typename MatrixType , typename ElementLess >
void tu::matrix_reorder_columns ( const MatrixType &  matrix,
size_t  row_first,
size_t  row_beyond,
size_t  column_first,
size_t  column_beyond,
ElementLess  element_less,
permutation result_permutation 
)
inline

Reorders the columns of a matrix by setting up a permutation to view the ordering.

Parameters
matrixThe given matrix, which is unchanged.
row_firstFirst row for comparison
row_beyondBeyond row for comparison
column_firstFirst column to reorder
column_beyondBeyond column to reorder
element_lessCompare functor
result_permutationResulting permutation
template<typename MatrixType , typename ElementLess >
void tu::matrix_reorder_columns ( MatrixType &  matrix,
size_t  row_first,
size_t  row_beyond,
size_t  column_first,
size_t  column_beyond,
ElementLess  element_less 
)
inline

Reorders the columns of a matrix.

Parameters
matrixThe given matrix, which is unchanged.
row_firstFirst row for comparison
row_beyondBeyond row for comparison
column_firstFirst column to reorder
column_beyondBeyond column to reorder
element_lessCompare functor
template<typename MatrixType , typename ElementLess >
void tu::matrix_reorder_rows ( const MatrixType &  matrix,
size_t  row_first,
size_t  row_beyond,
size_t  column_first,
size_t  column_beyond,
ElementLess  element_less,
permutation result_permutation 
)
inline

Reorders the rows of a matrix by setting up a permutation to view the ordering.

Parameters
matrixThe given matrix, which is unchanged.
row_firstFirst row to reorder
row_beyondBeyond row to reorder
column_firstFirst column for comparison
column_beyondBeyond column for comparison
element_lessCompare functor
result_permutationResulting permutation
template<typename MatrixType , typename ElementLess >
void tu::matrix_reorder_rows ( MatrixType &  matrix,
size_t  row_first,
size_t  row_beyond,
size_t  column_first,
size_t  column_beyond,
ElementLess  element_less 
)
inline

Reorders the rows of a matrix.

Parameters
matrixThe given matrix, which is unchanged.
row_firstFirst row for comparison
row_beyondBeyond row for comparison
column_firstFirst column to reorder
column_beyondBeyond column to reorder
element_lessCompare functor
template<typename Matrix >
void tu::matrix_row_combine ( Matrix matrix,
size_t  row1,
size_t  row2,
long long  ul,
long long  ur,
long long  ll,
long long  lr 
)
template<typename Matrix >
bool tu::matrix_row_gcd ( Matrix matrix,
size_t  row1,
size_t  row2,
size_t  column 
)

Special case: If both are |.| = g then s should not be zero.

template<typename Matrix >
bool tu::matrix_row_zero ( const Matrix matrix,
size_t  row,
size_t  column_first,
size_t  column_beyond 
)
template<typename MatrixType >
void tu::matrix_set_perm1 ( matrix_permuted< MatrixType > &  matrix,
const permutation permutation 
)
inline
template<typename MatrixType >
void tu::matrix_set_perm1 ( matrix_transposed< MatrixType > &  matrix,
const permutation permutation 
)
inline
template<typename MatrixType >
void tu::matrix_set_perm2 ( matrix_permuted< MatrixType > &  matrix,
const permutation permutation 
)
inline
template<typename MatrixType >
void tu::matrix_set_perm2 ( matrix_transposed< MatrixType > &  matrix,
const permutation permutation 
)
inline
template<typename MatrixType >
void tu::matrix_set_value ( matrix_transposed< MatrixType > &  matrix,
size_t  row,
size_t  column,
typename MatrixType::value_type  value 
)
inline

Free function to set a matrix value of a permuted matrix.

Parameters
matrixThe transposed matrix
rowRow index
columnColumn index
valueNew value
template<typename MatrixType >
void tu::matrix_set_value ( matrix_permuted< MatrixType > &  matrix,
size_t  row,
size_t  column,
typename MatrixType::value_type  value 
)
inline

Free function to set a matrix value of a permuted matrix.

Parameters
matrixThe permuted matrix
rowRow index
columnColumn index
valueNew value
template<typename MatrixType >
void tu::matrix_set_value ( matrix_transposed< const MatrixType > &  matrix,
size_t  row,
size_t  column,
typename MatrixType::value_type  value 
)
inline

Dummy function for the version with writable orignal matrix.

template<typename MatrixType >
void tu::matrix_set_value ( matrix_permuted< const MatrixType > &  matrix,
size_t  row,
size_t  column,
typename MatrixType::value_type  value 
)
inline

Dummy function for the version with writable orignal matrix.

template<typename MatroidType , typename MatrixType >
void tu::matroid_apply_column_permutation ( MatroidType &  matroid,
MatrixType &  matrix,
const permutation perm 
)
inline

Applies a column permutation to a given matroid and its representation matrix, i.e. performs the necessary columns swaps at the matroid / matrix directly.

Parameters
matroidThe given matroid to be changed
matrixRepresentation matrix of the given matroid
permThe given permutation
template<typename MatroidType , typename MatrixType >
void tu::matroid_apply_row_permutation ( MatroidType &  matroid,
MatrixType &  matrix,
const permutation perm 
)
inline

Applies a row permutation to a given matroid and its representation matrix, i.e. performs the necessary row swaps at the matroid / matrix directly.

Parameters
matroidThe given matroid to be changed
matrixRepresentation matrix of the given matroid
permThe given permutation
template<typename MatroidType >
void tu::matroid_binary_pivot ( matroid_transposed< MatroidType > &  matroid,
size_t  i,
size_t  j 
)

Free function to perform a pivot on a matroid.

Parameters
matroidThe given matroid
iRow index
jColumn index
template<typename MatroidType >
void tu::matroid_binary_pivot ( matroid_permuted< MatroidType > &  matroid,
size_t  i,
size_t  j 
)

Free function to perform a pivot on a matroid.

Parameters
matroidThe given matroid
iRow index
jColumn index
template<typename NameType >
void tu::matroid_binary_pivot ( matroid< NameType > &  matroid,
size_t  i,
size_t  j 
)

Free function to perform a pivot on a matroid.

Parameters
matroidThe given matroid
iRow index
jColumn index
template<typename MatroidType , typename MatrixType >
void tu::matroid_binary_pivot ( MatroidType &  matroid,
MatrixType &  matrix,
size_t  i,
size_t  j 
)

Free function to perform a pivot on a matroid and its representation matrix.

Parameters
matroidThe given matroid
matrixRepresentation matrix of the given matroid
iRow index
jColumn index
template<typename MatroidType >
std::set<int> tu::matroid_elements ( const MatroidType &  matroid)

Returns all matroid elements.

Parameters
matroidThe given matroid
Returns
A set of matroid elements
template<typename MatroidType >
void tu::matroid_permute1 ( matroid_transposed< MatroidType > &  matroid,
size_t  index1,
size_t  index2 
)
inline

Free function to permute two rows.

Parameters
matroidThe given matroid
index1First row index
index2Second row index
template<typename MatroidType >
void tu::matroid_permute1 ( matroid_permuted< MatroidType > &  matroid,
size_t  index1,
size_t  index2 
)
inline

Free function to permute two rows.

Parameters
matroidThe given matroid
index1First row index
index2Second row index
template<typename NameType >
void tu::matroid_permute1 ( matroid< NameType > &  matroid,
size_t  index1,
size_t  index2 
)
inline

Free function to permute two rows.

Parameters
matroidThe given matroid
index1First row index
index2Second row index
template<typename MatroidType , typename MatrixType >
void tu::matroid_permute1 ( MatroidType &  matroid,
MatrixType &  matrix,
size_t  index1,
size_t  index2 
)
inline

Free function to permute two rows of a matroid and its representation matrix.

Parameters
matroidThe given matroid
matrixRepresentation matrix of the given matroid
index1First row index
index2Second row index
template<typename MatroidType >
void tu::matroid_permute2 ( matroid_transposed< MatroidType > &  matroid,
size_t  index1,
size_t  index2 
)
inline

Free function to permute two columns.

Parameters
matroidThe given matroid
index1First column index
index2Second column index
template<typename MatroidType >
void tu::matroid_permute2 ( matroid_permuted< MatroidType > &  matroid,
size_t  index1,
size_t  index2 
)
inline

Free function to permute two columns.

Parameters
matroidThe given matroid
index1First column index
index2Second column index
template<typename NameType >
void tu::matroid_permute2 ( matroid< NameType > &  matroid,
size_t  index1,
size_t  index2 
)
inline

Free function to permute two columns.

Parameters
matroidThe given matroid
index1First column index
index2Second column index
template<typename MatroidType , typename MatrixType >
void tu::matroid_permute2 ( MatroidType &  matroid,
MatrixType &  matrix,
size_t  index1,
size_t  index2 
)
inline

Free function to permute two columns of a matroid and its representation matrix.

Parameters
matroidThe given matroid
matrixRepresentation matrix of the given matroid
index1First column index
index2Second column index
template<typename MatroidType , typename MatrixType >
void tu::matroid_print ( const MatroidType &  matroid,
const MatrixType &  matrix 
)

Prints a matroid with its representation matrix.

Parameters
matroidThe given matrix
matrixRepresentation matrix of the given matroid
template<typename MatroidType , typename MatrixType >
void tu::matroid_print_minor ( const MatroidType &  matroid,
const MatrixType &  matrix,
size_t  height,
size_t  width 
)

Prints a matroid minor with its representation matrix.

Parameters
matroidThe given matrix
matrixRepresentation matrix of the given matroid
heightNumber of rows of the minor
widthNumber of columns of the minor
template<typename MatroidType , typename MatrixType , typename ElementLess >
void tu::matroid_reorder_columns ( MatroidType &  matroid,
MatrixType &  matrix,
size_t  row_first,
size_t  row_beyond,
size_t  column_first,
size_t  column_beyond,
ElementLess  element_less 
)
inline

Reorders the columns of a matroid and its representation matrix by setting up a permutation to view the ordering.

Parameters
matroidThe given matrix, which is unchanged.
matrixRepresentation matrix of the given matroid
row_firstFirst row for comparison
row_beyondBeyond row for comparison
column_firstFirst column to reorder
column_beyondBeyond column to reorder
element_lessCompare functor
template<typename MatroidType , typename MatrixType , typename ElementLess >
void tu::matroid_reorder_rows ( MatroidType &  matroid,
MatrixType &  matrix,
size_t  row_first,
size_t  row_beyond,
size_t  column_first,
size_t  column_beyond,
ElementLess  element_less 
)
inline

Reorders the rows of a matroid and its representation matrix by setting up a permutation to view the ordering.

Parameters
matroidThe given matrix, which is unchanged.
matrixRepresentation matrix of the given matroid
row_firstFirst row to reorder
row_beyondBeyond row to reorder
column_firstFirst column for comparison
column_beyondBeyond column for comparison
element_lessCompare functor
template<typename MatroidType , typename MatrixType >
void tu::matroid_separate ( MatroidType &  matroid,
MatrixType &  matrix,
const int  type,
const std::pair< size_t, size_t > &  split,
integer_matroid upper_left_matroid,
integer_matrix upper_left_matrix,
integer_matroid lower_right_matroid,
integer_matrix lower_right_matrix 
)

Upper left

Lower right

template<typename MatroidType >
void tu::matroid_set_perm1 ( matroid_permuted< MatroidType > &  matroid,
const permutation permutation 
)
inline
template<typename MatroidType >
void tu::matroid_set_perm1 ( matroid_transposed< MatroidType > &  matroid,
const permutation permutation 
)
inline
template<typename MatroidType >
void tu::matroid_set_perm2 ( matroid_permuted< MatroidType > &  matroid,
const permutation permutation 
)
inline
template<typename MatroidType >
void tu::matroid_set_perm2 ( matroid_transposed< MatroidType > &  matroid,
const permutation permutation 
)
inline
std::ostream & tu::operator<< ( std::ostream &  stream,
const tu::matroid_graph graph 
)

Output operator for a matroid graph

Parameters
streamOutput stream
graphThe given matroid graph
Returns
The output stream after writing to it
std::ostream & tu::operator<< ( std::ostream &  stream,
logger log 
)

Streams a line of a logger object and flushes the output stream.

Parameters
Outputstream
Loggerobject
Returns
Output stream
std::ostream& tu::operator<< ( std::ostream &  stream,
const combination c 
)
inline

Output operator for combinations.

Parameters
streamA given output stream
pA given permutation
Returns
The stream after writing the combination
std::ostream& tu::operator<< ( std::ostream &  stream,
const binary_linear_space space 
)
inline

Streams a binary linear vector space.

Parameters
streamOutput stream
spaceThe given vector space
Returns
The output stream after processing
std::ostream& tu::operator<< ( std::ostream &  stream,
const permutation p 
)
inline

Output operator for permutations.

Parameters
streamA given output stream
pA given permutation
Returns
The stream after writing the permutation
bool tu::operator== ( const permutation p,
const permutation q 
)
inline

Compares two permutations for equality.

Parameters
pFirst permutation
qSecond permutation
Returns
true if and only if the permutations are the same
template<typename MatrixType >
rank_distribution tu::partition ( matrix_permuted< const MatrixType > &  matrix,
size_t &  top_left_height,
size_t &  top_left_width,
size_t &  bottom_right_height,
size_t &  bottom_right_width 
)
inline

Partitioning routine. Takes a 3-separation of a minor and tries to enlarge it to a 3-separation of the whole matroid.

Parameters
matrixRepresentation matrix of the matroid
top_left_heightNumber of rows in the top-left separation part
top_left_widthNumber of columns in the top-left separation part
bottom_right_heightNumber of rows in the bottom-right separation part
bottom_right_widthNumber of columns in the bottom-right separation part
Returns
A rank distribution, if enlarging was possible or RANK_TOO_HIGH if not

Repeat until no rows/columns can be shifted.

Setup bottom left rows

Setup top right rows

Check we have rank sum of 2

Check all other rows

Setup top right columns

Setup bottom left columns

Assert we have rank sum of 2

Check all other columns

template<typename M >
bool tu::sign_matrix ( M &  matrix,
submatrix_indices violator 
)

Generic signing function of a matrix, which might also be const. Running time: O (height * width^2)

Parameters
matrixThe given matrix
violatorPointer to violator indices to be filled.
Returns
true if and only if the matrix is signed already.

Go trough column by column.

There is a non-zero column right of the already-handled submatrix.

Start a BFS on bipartite graph of the submatrix and look for shortest paths from first 1 to all others

Evaluate matrix-entries on the shortest paths

Checking changes

Find the violator, going along the path

Fill violator data

We are not just testing, so swap the sign on a one.

Augment submatrix by rows with 1 in the new column.

Handled upper-left submatrix and lower-right submatrix are disconnected

A zero column can be skipped, as it is handled by definition.

Found a nonzero column and swap ones to the top.

CMR_EXPORT bool tu::sign_matrix ( integer_matrix matrix)

Signes a given matrix to be a signed version of its support matrix. Running time: O(height * width * min(height, width))

Parameters
matrixThe matrix to be signed
Returns
true if and only if any change was necessary
template<typename Matrix >
void tu::smith_normal_form_diagonal ( Matrix input_matrix,
std::vector< int > &  diagonal 
)

Ensure it is positive

Reduce entries below pivot

Reduce entries right of pivot

Test remaining for gcd and move to column

template<class Less >
void tu::sort ( permutation permutation,
size_t  first,
size_t  beyond,
Less &  less 
)
inline

Arranges part of a given permutation such that it can be used to view a vector in a sorted way.

Parameters
permutationThe given permutation to be changed
firstFirst index of the part
beyondBeyond index of the part
lessFunctor to compare two elements
template<class Less >
void tu::sort ( permutation permutation,
Less &  less 
)
inline

Arranges a given permutation such that it can be used to view a vector in a sorted way.

Parameters
permutationThe given permutation to be changed
lessFunctor to compare two elements
bool tu::submatrix_camion ( const integer_matrix matrix,
const submatrix_indices submatrix 
)

Checks Camion's criterion for total unimodularity.

Parameters
matrixA given integer matrix
submatrixMatrix-indices describing a submatrix B of A
Returns
true iff B is not Eulerian or 1^T * B * 1 = 0 (mod 4)
CMR_EXPORT int tu::submatrix_determinant ( const integer_matrix matrix,
const submatrix_indices submatrix 
)

Calculates a subdeterminant of the given matrix.

Parameters
matrixA given integer matrix
submatrixMatrix-indices describing a submatrix
Returns
The submatrix' determinant

LU-factorize

Calculate the determinant as a product of the diagonal entries, taking the permutation matrix into account.

CMR_EXPORT void tu::support_matrix ( integer_matrix matrix)

Makes the matrix its own support matrix.

Parameters
matrixThe given matrix
template<typename Matrix >
bool tu::test_k_modularity ( const Matrix matrix,
size_t &  rank,
unsigned int *  pk,
bool  enforce_unimodularity,
log_level  level 
)

In case we need k, let's compute Smith Normal Form of submatrix under column basis.

Assert that Basis * Transformed = Original

Compute Smith Normal Form and check that the gcd of all full-rank submatrices is 1.

Create submatrix of non-basis.

template<typename M >
TransposedMatrix<M> tu::transposedMatrix ( M &  matrix)

Returns a TransposedMatrix for matrix.