Bicoloring postprocessing seems to favor removing column colors and not row colors

[gdalle]

... which is bad for autodiff since reverse passes are more expensive than forward passes.

I wonder if this could be due to the order in which we iterate through the stars / trees? It doesn't seem to be influenced by the order in which we greedily color vertices.

A = diagm(ones(Bool, 3))
problem = ColoringProblem{:nonsymmetric,:bidirectional}()
algo = GreedyColoringAlgorithm(RandomOrder(); postprocessing=true)
for _ in 1:10
    result = coloring(A, problem, algo)
    @show column_colors(result), row_colors(result)
end

Perhaps in this line

SparseMatrixColorings.jl/src/coloring.jl

Line 646 in 7efc6bc

for j in axes(S, 2)

we could use the computed vertex order instead, or its reverse order? Not sure which one makes more sense?

I noticed this because in JuliaDiff/DifferentiationInterface.jl#864 I sometimes got empty forward groups but never empty reverse groups.

@amontoison

[amontoison]

Sorry for not answering earlier @gdalle!
I have an answer for you, when we have a trivial star (or trivial tree), we use max(i,j) of a edge (i,j) as the hub h: https://github.com/gdalle/SparseMatrixColorings.jl/blob/main/src/coloring.jl#L187

Reminder: we use the first indices of the columns of the augmented matrix H for the columns J and the last indices of the columns of H for the rows of J.

So it means that we try to neutralize the color of the columns of J when we have the possibility and not the rows of J (because index_rows_J > index_col_J and thus index_rows_J is used as hub by default):
https://github.com/gdalle/SparseMatrixColorings.jl/blob/main/src/coloring.jl#L656-L662

We should add an option to specify what we want to prefer to neutralize (rows or columns) and by default neutralize the rows (most common case).

[amontoison]

Guillaume, I explained it in the paper but we never implemented it in SMC.jl.
In the case of symmetric coloring we don't care if we neutralize the color of a row or a column of a Hessian (same cost) but it is quite different for bicoloring of a Jacobian.

In the procedure post_processing (line 21 and 22):

• i > j and we prefer to neutralize columns: set i as hub
• i > j and we prefer to neutralize rows: set j as hub
• i > j => the column i of H is a row of J and the column j of H is a column of J

[gdalle]

I think there are two aspects we tweak in the postprocessing:

• Given a trivial star, which vertex do we pick as the hub? This is what you discuss above, and for that we can prioritize either row or column vertices. We should make this configurable.
• In which order do we enumerate the trivial stars? This is what I mentioned in my first comment, and at the moment we don't do it in a principled way. The intuition for vertex orderings like LargestFirst is "make the hardest, most constraining choices first". Similarly, we could go from the most constraining vertex/star (as defined by our ordering) to the least constraining one when deciding which colors are used or not. Indeed, while this loop is order-agnostic, that loop is not.

[amontoison]

For the second point, it doesn't impact the number of colors in bicoloring if we prioritize row OR column vertices.

But if we want to minimize the total of colors instead (bicoloring or symmetric coloring) then the ordering of the trivial edges is relevant.

[gdalle]

For the second point, it doesn't impact the number of colors in bicoloring if we prioritize row OR column vertices.

Everything can impact the number of colors since we make decisions sequentially

[amontoison]

@gdalle I explained you why on Messenger. Please find a counter-example if you think I am wrong.

[amontoison]

For the case where we want to minimize ncolors_row + ncolors_column and not only ncolors_row or ncolors_col I think we can do something optimal.

What we can do is to do a first pass on all trivial edges to mark the ones where both vertices can be used as hub.
For the other ones, we just update the hub if needed but we don't need them for the second pass because the hub can be imposed (one vertex of the trivial edge is already a hub in a normal tree or star).

In the second pass we only have the trivial edges where the hub must be choosen by us and what we want to do is to use as hub the vertices that appear in the maximum number of trivial edges still listed in the second pass.

If we compute a "degree" for each vertex in terms of apperance in trivial edges (where both vertices can be hub) then using the vertex with the biggest degree in each trivial edge as hub should give us the best number of colors possible after a post-processing.

Degree is probably not the best word, I can reformulate as "in how many trivial edges of the second pass each vertex is part of".

Once again, I may be wrong so don't hesitate to share your point of view 🙂

[gdalle]

Well, if you claim we can do something "optimal" in an NP-hard problem, the burden of proof is on you ^^ This remains a greedy coloring and postprocessing algorithm, so my opinion is that we can only hope for a "good enough" outcome, and not an optimal one.

About the importance of the order, let us assume that we always select the row $i$ as the hub for a trivial star (edge) $s = (i, j)$.
The point I'm trying to make is that when we enumerate $s_1, ..., s_k$ and we make a decision for star $s_\ell = (i_\ell, j_\ell)$, we may set used[color[i_l]] = true. This in turns influences what we do for $s_{\ell + 1}$ and all the stars after it, because if a color is already used and another isn't, the hub must follow that first. Thus, the order in which we enumerate trivial stars when removing colors matters, just as much as the order in which we enumerate vertices when adding colors.
My suggestion is to figure out if we can deduce a reasonable order over stars from the order over vertices that we already have.

Your suggestion in the message above is that we can set the hubs over all of these stars "in parallel", without going through them in a sequence and inducing dependencies.
It may be better than what we do now but I don't think it is optimal either. When you pick the hub for a star, you make a local decision that seems good, it is a greedy procedure. But the real problem you need to solve is global: what is the best set of colors to neutralize, or here the best set of vertices to pick as hubs. Just because these selected vertices have a larger degree individually doesn't preclude the existence of another set of candidate hubs which looks worse individually but works better as a whole?

[amontoison]

Yes exactly! We don't have any dependency if we only minimize ncolors_row or ncolors_col because in bicoloring, all edges (i,j) has exactly one vertex that is a row of J and one vertex that is a column of J (no overlap).

We can always impose as hub the vertices in one partition whatever is the order of treatment of the trivial edges.

We must use vertices in one partition as hub to minimize ncolors_row or ncolors_col and whatever is the order, we can always do that.

Another way to see it is:
If we want to minimize ncolors_row, we must always use vertices that are in the partition of columns of the Jacobian in J as hub.
And because we exactly have one vertex per trivial edge in this partition, we know which one to take and we are sure that the "row vertex" was never used as a hub previously before whatever is the order of treatment of the trivial edges.
It was never marked as hub because it is in the other partition.
Conclusion, we can always use the "column vertex" as hub.

[gdalle]

I think I have a clearer picture now. In the general symmetric case, after the first pass where we handle imposed hubs, the problem can be rephrased as: "select the set $H$ of hubs with the smallest number of distinct new colors $C(H)$ such that each trivial edge without an imposed hub is incident to at least one of them". It is a variant of the vertex cover problem on a subgraph of $G$, which is NP-hard, so a degree-based heuristic like the one you suggest may be clever but not optimal.
In the bicoloring case where we want to minimize either ncolors_row or ncolors_col, you noticed that the graph is bipartite which makes the problem easy. However, do we agree that in the symmetric case it is not obvious what the optimal $H$ looks like?