What does `--method` mean?

[vasdommes]

Hi, I want to compute a RREF where the leftmost (independent) columns are always chosen as pivots.
From what I saw in a couple of test cases, --method=1 does the job.
Is it true? If so, I'd love to see it mentioned somewhere in the docs.

sparserref --help is not very helpful:

-m, --method                     method of RREF  [default: 0]

And this comment in the code also does not go into detail:

SparseRREF/sparse_rref.h

Lines 159 to 162 in 3f1aae6

	// rref_option
	// method 0: right and left search
	// method 1: only right search
	// method 2: hibrid

[Li-Xiang-Ideal]

Hi @vasdommes, this piece of code was not written by me, but I will try to answer your question based on my understanding.

Roughly speaking, the left search in opt->method means searching for pivots column by column, while the right search means searching for pivots row by row. More specifically, opt->method determines which pivot-searching function is called (pivots_search or pivots_search_right) through the following code:

SparseRREF/sparse_mat.h

Lines 1053 to 1067 in 3f1aae6

	bool only_right_search = true;
	sparse_mat<bool, index_t> tranmat;
	if (opt->method != 1) {
	only_right_search = false;
	tranmat = sparse_mat<bool, index_t>(mat.ncol, mat.nrow);
	sparse_mat_subview<T, index_t> matview(mat, leftrows);
	sparse_mat_transpose_replace(tranmat, matview, &pool);
	// sort pivots by nnz, it will be faster
	std::ranges::stable_sort(leftcols, std::less{},
	[&tranmat](size_t r) {
	return tranmat[r].nnz();
	});
	}

and

SparseRREF/sparse_mat.h

Lines 1082 to 1087 in 3f1aae6

	if (only_right_search) {
	ps = pivots_search_right(mat, leftrows, leftcols, opt->col_weight);
	}
	else {
	ps = pivots_search(mat, tranmat, leftrows, leftcols, opt->col_weight);
	}

The pivot-searching functions are defined as follows:

SparseRREF/sparse_mat.h

Lines 238 to 421 in 3f1aae6

	// first choose the pivot with minimal col_weight
	// if the col_weight is negative, we do not choose it
	template <typename T, typename index_t>
	std::vector<pivot_t<index_t>> pivots_search(
	const sparse_mat<T, index_t>& mat, const sparse_mat<bool, index_t>& tranmat,
	const std::vector<size_t>& leftrows, const std::vector<index_t>& leftcols,
	const std::function<int64_t(int64_t)>& col_weight = [](int64_t i) { return i; }) {
	std::list<pivot_t<index_t>> pivots;
	std::unordered_set<index_t> dict;
	dict.reserve((size_t)4096);
	// leftlook first
	for (auto col : leftcols) {
	if (tranmat[col].nnz() == 0)
	continue;
	// negative weight means that we do not want to select this column
	if (col_weight(col) < 0)
	continue;
	index_t row;
	size_t mnnz = SIZE_MAX;
	bool flag = true;
	for (auto r : tranmat[col].index_span()) {
	flag = (dict.count(r) == 0);
	if (!flag)
	break;
	size_t newnnz = mat[r].nnz();
	if (newnnz < mnnz) {
	row = r;
	mnnz = newnnz;
	}
	// make the result stable
	else if (newnnz == mnnz) {
	if (r < row)
	row = r;
	}
	}
	if (!flag)
	continue;
	if (mnnz != SIZE_MAX) {
	pivots.emplace_front(row, col);
	dict.insert(row);
	}
	}
	// rightlook then
	dict.clear();
	for (auto [r, c] : pivots)
	dict.insert(c);
	for (auto row : leftrows) {
	index_t col = 0;
	size_t mnnz = SIZE_MAX;
	bool flag = true;
	for (auto c : mat[row].index_span()) {
	// negative weight means that we do not want to select this column
	if (col_weight(c) < 0)
	continue;
	flag = (dict.count(c) == 0);
	if (!flag)
	break;
	if (tranmat[c].nnz() < mnnz) {
	mnnz = tranmat[c].nnz();
	col = c;
	}
	// make the result stable
	else if (tranmat[c].nnz() == mnnz) {
	if (col_weight(c) < col_weight(col))
	col = c;
	else if (col_weight(c) == col_weight(col) && c < col)
	col = c;
	}
	}
	if (!flag)
	continue;
	if (mnnz != SIZE_MAX) {
	pivots.emplace_back(row, col);
	dict.insert(col);
	}
	}
	std::vector<pivot_t<index_t>> result(pivots.begin(), pivots.end());
	return result;
	}
	// first choose the pivot with minimal col_weight
	// if the col_weight is negative, we do not choose it
	// only right search version, do not need the full tranpose
	template <typename T, typename index_t>
	std::vector<pivot_t<index_t>> pivots_search_right(const sparse_mat<T, index_t>& mat,
	const std::vector<size_t>& leftrows, const std::vector<index_t>& leftcols,
	const std::function<int64_t(int64_t)>& col_weight = [](int64_t i) { return i; }) {
	std::vector<pivot_t<index_t>> pivots;
	std::unordered_set<index_t> c_dict;
	c_dict.reserve((size_t)4096);
	for (auto row : leftrows) {
	index_t col = 0;
	size_t mnnz = SIZE_MAX;
	bool flag = true;
	for (auto c : mat[row].index_span()) {
	// negative weight means that we do not want to select this column
	if (col_weight(c) < 0)
	continue;
	flag = (c_dict.count(c) == 0);
	if (!flag)
	break;
	if (mnnz > 0) {
	mnnz = 0;
	col = c;
	}
	// make the result stable
	else if (mnnz == 0) {
	if (col_weight(c) < col_weight(col))
	col = c;
	else if (col_weight(c) == col_weight(col) && c < col)
	col = c;
	}
	}
	if (!flag)
	continue;
	if (mnnz != SIZE_MAX) {
	pivots.emplace_back(row, col);
	c_dict.insert(col);
	}
	}
	return pivots;
	}
	template <typename T, typename index_t>
	std::vector<pivot_t<index_t>> pivots_search_left(const sparse_mat<T, index_t>& mat,
	const sparse_mat<bool, index_t>& tranmat, const std::vector<index_t>& rowpivs,
	const std::function<int64_t(int64_t)>& col_weight = [](int64_t i) { return i; }) {
	std::vector<pivot_t<index_t>> pivots;
	std::unordered_set<index_t> r_dict;
	r_dict.reserve((size_t)4096);
	for (auto col = 0; col < tranmat.nrow; col++) {
	if (tranmat[col].nnz() == 0 || col_weight(col) < 0)
	continue;
	index_t row = 0;
	size_t mnnz = SIZE_MAX;
	bool flag = true;
	for (auto r : tranmat[col].index_span()) {
	if (rowpivs[r] != index_sval<index_t>())
	continue;
	flag = !r_dict.contains(r);
	if (!flag)
	break;
	if (mnnz > mat[r].nnz()) {
	mnnz = mat[r].nnz();
	row = r;
	}
	// make the result stable
	else if (mnnz == mat[r].nnz()) {
	if (r < row)
	row = r;
	}
	}
	if (!flag)
	continue;
	if (mnnz != SIZE_MAX) {
	pivots.emplace_back(row, col);
	r_dict.insert(row);
	}
	}
	// reverse the ordering of pivots
	for (size_t i = 0; i < pivots.size() / 2; i++) {
	std::swap(pivots[i], pivots[pivots.size() - 1 - i]);
	}
	return pivots;
	}

Now back to your question. If you mean that you want the leftmost (independent) columns to be always chosen as pivots, then you should use opt->method = 1, i.e., only use right search. In principle, column-wise pivot search and row-wise pivot search should give the same result. However, for performance reasons, during left search we sort leftcols according to nnz (number of nonzeros) of each column (see lines 1062-1066 in the first code block). This means that during left search, instead of choosing the leftmost (independent) column as the pivot, we choose the (independent) column with the fewest nonzeros. In contrast, right search still chooses the leftmost (independent) column as the pivot.

If there is anything incorrect in my understanding, please feel free to correct me @munuxi. Hope this helps.

[vasdommes]

@Li-Xiang-Ideal thanks for the detailed explanation!
Yes, opt->method = 1 is what I want indeed.