From 016e0b087f9e7e1c8636291ba3d30148ddc60363 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Beno=C3=AEt=20Legat?= Date: Mon, 19 Oct 2026 11:37:06 +0200 Subject: [PATCH 1/6] Speed up CSC substitution decompression for one-edge trees --- src/decompression.jl | 104 ++++++++++++++++++++++++++++-------------- test/decompression.jl | 89 +++++++++++++++++++++++++++++++++++- 2 files changed, 158 insertions(+), 35 deletions(-) diff --git a/src/decompression.jl b/src/decompression.jl index 35e0e71e..27534d53 100644 --- a/src/decompression.jl +++ b/src/decompression.jl @@ -610,6 +610,46 @@ function decompress!( return A end +@inline function set_symmetric_csc!( + nzA, A_colptr, lower_triangle_offsets, upper_triangle_offsets, counter, i, j, val, uplo +) + #! format: off + # A[i,j] is in the lower triangular part of A + if in_triangle(i, j, :L) + # uplo = :L or uplo = :F + # A[i,j] is stored at index_ij = (A.colptr[j+1] - offset_L) in A.nzval + if uplo != :U + nzind = A_colptr[j + 1] - lower_triangle_offsets[counter] + nzA[nzind] = val + end + + # uplo = :U or uplo = :F + # A[j,i] is stored at index_ji = (A.colptr[i] + offset_U) in A.nzval + if uplo != :L + nzind = A_colptr[i] + upper_triangle_offsets[counter] + nzA[nzind] = val + end + + # A[i,j] is in the upper triangular part of A + else + # uplo = :U or uplo = :F + # A[i,j] is stored at index_ij = (A.colptr[j] + offset_U) in A.nzval + if uplo != :L + nzind = A_colptr[j] + upper_triangle_offsets[counter] + nzA[nzind] = val + end + + # uplo = :L or uplo = :F + # A[j,i] is stored at index_ji = (A.colptr[i+1] - offset_L) in A.nzval + if uplo != :U + nzind = A_colptr[i + 1] - lower_triangle_offsets[counter] + nzA[nzind] = val + end + end + #! format: on + return nothing +end + function decompress_csc!( nzA::AbstractVector{R}, A_colptr::AbstractVector{<:Integer}, @@ -667,6 +707,25 @@ function decompress_csc!( first = tree_edge_indices[k] last = tree_edge_indices[k + 1] - 1 + # A one-edge tree needs no substitution or scratch-buffer updates. + if first == last + (i, j) = reverse_bfs_orders[first] + counter += 1 + val = B[i, color[j]] + set_symmetric_csc!( + nzA, + A_colptr, + lower_triangle_offsets, + upper_triangle_offsets, + counter, + i, + j, + val, + uplo, + ) + continue + end + # Reset the buffer to zero for all vertices in the tree (except the root) for pos in first:last (vertex, _) = reverse_bfs_orders[pos] @@ -682,40 +741,17 @@ function decompress_csc!( val = B[i, color[j]] - buffer_right_type[i] buffer_right_type[j] = buffer_right_type[j] + val - #! format: off - # A[i,j] is in the lower triangular part of A - if in_triangle(i, j, :L) - # uplo = :L or uplo = :F - # A[i,j] is stored at index_ij = (A.colptr[j+1] - offset_L) in A.nzval - if uplo != :U - nzind = A_colptr[j + 1] - lower_triangle_offsets[counter] - nzA[nzind] = val - end - - # uplo = :U or uplo = :F - # A[j,i] is stored at index_ji = (A.colptr[i] + offset_U) in A.nzval - if uplo != :L - nzind = A_colptr[i] + upper_triangle_offsets[counter] - nzA[nzind] = val - end - - # A[i,j] is in the upper triangular part of A - else - # uplo = :U or uplo = :F - # A[i,j] is stored at index_ij = (A.colptr[j] + offset_U) in A.nzval - if uplo != :L - nzind = A_colptr[j] + upper_triangle_offsets[counter] - nzA[nzind] = val - end - - # uplo = :L or uplo = :F - # A[j,i] is stored at index_ji = (A.colptr[i+1] - offset_L) in A.nzval - if uplo != :U - nzind = A_colptr[i + 1] - lower_triangle_offsets[counter] - nzA[nzind] = val - end - end - #! format: on + set_symmetric_csc!( + nzA, + A_colptr, + lower_triangle_offsets, + upper_triangle_offsets, + counter, + i, + j, + val, + uplo, + ) end end return nothing diff --git a/test/decompression.jl b/test/decompression.jl index 86efa122..8d57a6f2 100644 --- a/test/decompression.jl +++ b/test/decompression.jl @@ -1,7 +1,7 @@ using LinearAlgebra using SparseArrays using SparseMatrixColorings -using SparseMatrixColorings: decompress_csc! +using SparseMatrixColorings: SparsityPatternCSC, decompress_csc! using Test function test_star_csc_decompression(T, Ti, order, postprocessing) @@ -76,3 +76,90 @@ end test_star_csc_decompression(T, Ti, order, postprocessing) end end + +function test_tree_csc_decompression(T, Ti, order, postprocessing) + problem = ColoringProblem(; structure=:symmetric) + algorithm = GreedyColoringAlgorithm(order; decompression=:substitution, postprocessing) + @testset "$pattern" for pattern in (:blocks, :mixed, :diagonal, :empty) + A = spzeros(T, Ti, 10, 10) + if pattern in (:blocks, :mixed) + for j in 1:4, i in 1:(j - 1) + A[i, j] = A[j, i] = T(10i + j) + end + edges = if pattern == :blocks + [(6, 7), (6, 8), (7, 8)] + else + [(4, 5), (5, 6), (6, 7), (8, 9)] + end + for (i, j) in edges + A[i, j] = A[j, i] = T(10i + j) + end + end + if pattern != :empty + # Missing diagonals, a structural zero, and an empty last column. + A[2, 2] = T(2) + A[8, 8] = T(8) + A[3, 3] = one(T) + A[3, 3] = zero(T) + end + result = coloring(SparsityPatternCSC(A), problem, algorithm; decompression_eltype=T) + sizes = diff(result.tree_edge_indices) + if pattern == :blocks + @test !isempty(sizes) + @test all(==(1), sizes) + elseif pattern == :mixed + @test any(==(1), sizes) + @test any(>(1), sizes) + end + B = compress(A, result) + @testset "$uplo" for uplo in (:F, :L, :U) + expected = if uplo == :F + A + elseif uplo == :L + tril(A) + else + triu(A) + end + values = fill(T(NaN), nnz(expected)) + # Reuse the result with changed compressed values and stale scratch data. + for scale in (one(T), T(-2)) + compressed = scale .* B + fill!(result.buffer, T(NaN)) + @test decompress_csc!(values, expected.colptr, compressed, result, uplo) === + nothing + @test values == scale .* nonzeros(expected) + @test (@allocated decompress_csc!( + values, expected.colptr, compressed, result, uplo + )) == 0 + + dest = copy(expected) + fill!(nonzeros(dest), T(NaN)) + @test decompress!(dest, compressed, result, uplo) === dest + @test dest == scale .* expected + @test (@allocated decompress!(dest, compressed, result, uplo)) == 0 + end + + # The CSC entry point also accepts views of its input/output buffers. + storage = fill(T(NaN), 2nnz(expected)) + values_view = @view storage[1:2:end] + compressed_view = @view B[:, :] + @test decompress_csc!( + values_view, expected.colptr, compressed_view, result, uplo + ) === nothing + @test values_view == nonzeros(expected) + @test (@allocated decompress_csc!( + values_view, expected.colptr, compressed_view, result, uplo + )) == 0 + end + end + return nothing +end + +@testset "Tree CSC decompression" begin + @testset "$T, $Ti, $order, postprocessing=$postprocessing" for T in (Float32, Float64), + Ti in (Int32, Int64), order in (NaturalOrder(), LargestFirst()), + postprocessing in (false, true) + + test_tree_csc_decompression(T, Ti, order, postprocessing) + end +end From 7fe7c242a9e27731b7350a6be73438f48ad6e123 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Beno=C3=AEt=20Legat?= Date: Thu, 22 Oct 2026 12:10:34 +0200 Subject: [PATCH 2/6] No substitution for stars --- src/decompression.jl | 35 +++++++++++++++++++---------------- src/result.jl | 4 +++- test/decompression.jl | 28 ++++++++++++++++++++++++++-- 3 files changed, 48 insertions(+), 19 deletions(-) diff --git a/src/decompression.jl b/src/decompression.jl index 27534d53..1fab6266 100644 --- a/src/decompression.jl +++ b/src/decompression.jl @@ -661,6 +661,7 @@ function decompress_csc!( ag, color, reverse_bfs_orders, + is_star, tree_edge_indices, nt, diagonal_indices, @@ -707,22 +708,24 @@ function decompress_csc!( first = tree_edge_indices[k] last = tree_edge_indices[k + 1] - 1 - # A one-edge tree needs no substitution or scratch-buffer updates. - if first == last - (i, j) = reverse_bfs_orders[first] - counter += 1 - val = B[i, color[j]] - set_symmetric_csc!( - nzA, - A_colptr, - lower_triangle_offsets, - upper_triangle_offsets, - counter, - i, - j, - val, - uplo, - ) + # A star needs no substitution or scratch-buffer updates. + if is_star[k] + for pos in first:last + (i, j) = reverse_bfs_orders[pos] + counter += 1 + val = B[i, color[j]] + set_symmetric_csc!( + nzA, + A_colptr, + lower_triangle_offsets, + upper_triangle_offsets, + counter, + i, + j, + val, + uplo, + ) + end continue end diff --git a/src/result.jl b/src/result.jl index f7120670..66060a9f 100644 --- a/src/result.jl +++ b/src/result.jl @@ -387,6 +387,7 @@ struct TreeSetColoringResult{ color::Vector{T} group::GT reverse_bfs_orders::Vector{Tuple{T,T}} + is_star::Vector{Bool} tree_edge_indices::Vector{T} nt::T diagonal_indices::Vector{T} @@ -403,7 +404,7 @@ function TreeSetColoringResult( tree_set::TreeSet{<:Integer}, decompression_eltype::Type{R}, ) where {T<:Integer,R} - (; reverse_bfs_orders, tree_edge_indices, nt) = tree_set + (; reverse_bfs_orders, is_star, tree_edge_indices, nt) = tree_set (; S, nb_self_loops) = ag nvertices = length(color) group = group_by_color(T, color) @@ -484,6 +485,7 @@ function TreeSetColoringResult( color, group, reverse_bfs_orders, + is_star, tree_edge_indices, nt, diagonal_indices, diff --git a/test/decompression.jl b/test/decompression.jl index 8d57a6f2..02ec5419 100644 --- a/test/decompression.jl +++ b/test/decompression.jl @@ -80,10 +80,13 @@ end function test_tree_csc_decompression(T, Ti, order, postprocessing) problem = ColoringProblem(; structure=:symmetric) algorithm = GreedyColoringAlgorithm(order; decompression=:substitution, postprocessing) - @testset "$pattern" for pattern in (:blocks, :mixed, :diagonal, :empty) + @testset "$pattern" for pattern in ( + :blocks, :star_low, :star_middle, :star_high, :two_edges, :mixed, :diagonal, :empty + ) A = spzeros(T, Ti, 10, 10) if pattern in (:blocks, :mixed) - for j in 1:4, i in 1:(j - 1) + clique_size = pattern == :blocks ? 4 : 3 + for j in 1:clique_size, i in 1:(j - 1) A[i, j] = A[j, i] = T(10i + j) end edges = if pattern == :blocks @@ -95,6 +98,19 @@ function test_tree_csc_decompression(T, Ti, order, postprocessing) A[i, j] = A[j, i] = T(10i + j) end end + if pattern in (:star_low, :star_middle, :star_high, :two_edges) + hub = if pattern == :star_low + 1 + elseif pattern == :star_high + 9 + else + 5 + end + leaves = pattern == :two_edges ? (2, 7) : (2, 4, 7) + for leaf in leaves + A[leaf, hub] = A[hub, leaf] = T(10leaf + hub) + end + end if pattern != :empty # Missing diagonals, a structural zero, and an empty last column. A[2, 2] = T(2) @@ -107,9 +123,14 @@ function test_tree_csc_decompression(T, Ti, order, postprocessing) if pattern == :blocks @test !isempty(sizes) @test all(==(1), sizes) + @test all(result.is_star) + elseif pattern in (:star_low, :star_middle, :star_high, :two_edges) + @test sizes == [pattern == :two_edges ? 2 : 3] + @test all(result.is_star) elseif pattern == :mixed @test any(==(1), sizes) @test any(>(1), sizes) + @test any(!, result.is_star) end B = compress(A, result) @testset "$uplo" for uplo in (:F, :L, :U) @@ -128,6 +149,9 @@ function test_tree_csc_decompression(T, Ti, order, postprocessing) @test decompress_csc!(values, expected.colptr, compressed, result, uplo) === nothing @test values == scale .* nonzeros(expected) + if all(result.is_star) + @test all(isnan, result.buffer) + end @test (@allocated decompress_csc!( values, expected.colptr, compressed, result, uplo )) == 0 From d2e0a0e110da7285fc8daeaecf42dd036afb50a4 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Beno=C3=AEt=20Legat?= Date: Thu, 22 Oct 2026 15:34:12 +0200 Subject: [PATCH 3/6] else instead of continue --- src/decompression.jl | 53 ++++++++++++++++++++++---------------------- 1 file changed, 26 insertions(+), 27 deletions(-) diff --git a/src/decompression.jl b/src/decompression.jl index 1fab6266..ab2217ec 100644 --- a/src/decompression.jl +++ b/src/decompression.jl @@ -726,35 +726,34 @@ function decompress_csc!( uplo, ) end - continue - end - - # Reset the buffer to zero for all vertices in the tree (except the root) - for pos in first:last - (vertex, _) = reverse_bfs_orders[pos] - buffer_right_type[vertex] = zero(R) - end - # Reset the buffer to zero for the root vertex - (_, root) = reverse_bfs_orders[last] - buffer_right_type[root] = zero(R) + else + # Reset the buffer to zero for all vertices in the tree (except the root) + for pos in first:last + (vertex, _) = reverse_bfs_orders[pos] + buffer_right_type[vertex] = zero(R) + end + # Reset the buffer to zero for the root vertex + (_, root) = reverse_bfs_orders[last] + buffer_right_type[root] = zero(R) - for pos in first:last - (i, j) = reverse_bfs_orders[pos] - counter += 1 - val = B[i, color[j]] - buffer_right_type[i] - buffer_right_type[j] = buffer_right_type[j] + val + for pos in first:last + (i, j) = reverse_bfs_orders[pos] + counter += 1 + val = B[i, color[j]] - buffer_right_type[i] + buffer_right_type[j] = buffer_right_type[j] + val - set_symmetric_csc!( - nzA, - A_colptr, - lower_triangle_offsets, - upper_triangle_offsets, - counter, - i, - j, - val, - uplo, - ) + set_symmetric_csc!( + nzA, + A_colptr, + lower_triangle_offsets, + upper_triangle_offsets, + counter, + i, + j, + val, + uplo, + ) + end end end return nothing From 09331278140b641896819a588d8e92a098c83070 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Beno=C3=AEt=20Legat?= Date: Thu, 22 Oct 2026 15:45:34 +0200 Subject: [PATCH 4/6] Speed up stars in bicoloring as well --- src/decompression.jl | 116 ++++++++++++++++++++++++++---------------- src/result.jl | 5 +- test/decompression.jl | 78 ++++++++++++++++++++++++++++ 3 files changed, 155 insertions(+), 44 deletions(-) diff --git a/src/decompression.jl b/src/decompression.jl index ab2217ec..8165d627 100644 --- a/src/decompression.jl +++ b/src/decompression.jl @@ -869,6 +869,7 @@ function decompress!( symmetric_to_row, symmetric_to_column, reverse_bfs_orders, + is_star, tree_edge_indices, nt, buffer, @@ -888,26 +889,39 @@ function decompress!( first = tree_edge_indices[k] last = tree_edge_indices[k + 1] - 1 - # Reset the buffer to zero for all vertices in the tree (except the root) - for pos in first:last - (vertex, _) = reverse_bfs_orders[pos] - buffer_right_type[vertex] = zero(R) - end - # Reset the buffer to zero for the root vertex - (_, root) = reverse_bfs_orders[last] - buffer_right_type[root] = zero(R) + # A star needs no substitution or scratch-buffer updates. + if is_star[k] + for pos in first:last + (i, j) = reverse_bfs_orders[pos] + cj = symmetric_color[j] + if in_triangle(i, j, :L) + A[i - n, j] = Bc[i - n, symmetric_to_column[cj]] + else + A[j - n, i] = Br[symmetric_to_row[cj], i] + end + end + else + # Reset the buffer to zero for all vertices in the tree (except the root) + for pos in first:last + (vertex, _) = reverse_bfs_orders[pos] + buffer_right_type[vertex] = zero(R) + end + # Reset the buffer to zero for the root vertex + (_, root) = reverse_bfs_orders[last] + buffer_right_type[root] = zero(R) - for pos in first:last - (i, j) = reverse_bfs_orders[pos] - cj = symmetric_color[j] - if in_triangle(i, j, :L) - val = Bc[i - n, symmetric_to_column[cj]] - buffer_right_type[i] - buffer_right_type[j] = buffer_right_type[j] + val - A[i - n, j] = val - else - val = Br[symmetric_to_row[cj], i] - buffer_right_type[i] - buffer_right_type[j] = buffer_right_type[j] + val - A[j - n, i] = val + for pos in first:last + (i, j) = reverse_bfs_orders[pos] + cj = symmetric_color[j] + if in_triangle(i, j, :L) + val = Bc[i - n, symmetric_to_column[cj]] - buffer_right_type[i] + buffer_right_type[j] = buffer_right_type[j] + val + A[i - n, j] = val + else + val = Br[symmetric_to_row[cj], i] - buffer_right_type[i] + buffer_right_type[j] = buffer_right_type[j] + val + A[j - n, i] = val + end end end end @@ -925,6 +939,7 @@ function decompress!( symmetric_to_column, symmetric_to_row, reverse_bfs_orders, + is_star, tree_edge_indices, nt, A_indices, @@ -946,36 +961,51 @@ function decompress!( first = tree_edge_indices[k] last = tree_edge_indices[k + 1] - 1 - # Reset the buffer to zero for all vertices in the tree (except the root) - for pos in first:last - (vertex, _) = reverse_bfs_orders[pos] - buffer_right_type[vertex] = zero(R) - end - # Reset the buffer to zero for the root vertex - (_, root) = reverse_bfs_orders[last] - buffer_right_type[root] = zero(R) + # A star needs no substitution or scratch-buffer updates. + if is_star[k] + for pos in first:last + (i, j) = reverse_bfs_orders[pos] + cj = symmetric_color[j] + counter += 1 + val = if in_triangle(i, j, :L) + Bc[i - n, symmetric_to_column[cj]] + else + Br[symmetric_to_row[cj], i] + end + nzA[A_indices[counter]] = val + end + else + # Reset the buffer to zero for all vertices in the tree (except the root) + for pos in first:last + (vertex, _) = reverse_bfs_orders[pos] + buffer_right_type[vertex] = zero(R) + end + # Reset the buffer to zero for the root vertex + (_, root) = reverse_bfs_orders[last] + buffer_right_type[root] = zero(R) - for pos in first:last - (i, j) = reverse_bfs_orders[pos] - cj = symmetric_color[j] - counter += 1 + for pos in first:last + (i, j) = reverse_bfs_orders[pos] + cj = symmetric_color[j] + counter += 1 - #! format: off - if in_triangle(i, j, :L) - val = Bc[i - n, symmetric_to_column[cj]] - buffer_right_type[i] - buffer_right_type[j] = buffer_right_type[j] + val + #! format: off + if in_triangle(i, j, :L) + val = Bc[i - n, symmetric_to_column[cj]] - buffer_right_type[i] + buffer_right_type[j] = buffer_right_type[j] + val - # A[i-n,j] is stored at A_indices[counter] in nonzeros(A) - nzA[A_indices[counter]] = val + # A[i-n,j] is stored at A_indices[counter] in nonzeros(A) + nzA[A_indices[counter]] = val - else - val = Br[symmetric_to_row[cj], i] - buffer_right_type[i] - buffer_right_type[j] = buffer_right_type[j] + val + else + val = Br[symmetric_to_row[cj], i] - buffer_right_type[i] + buffer_right_type[j] = buffer_right_type[j] + val - # A[j-n,i] is stored at A_indices[counter] in nonzeros(A) - nzA[A_indices[counter]] = val + # A[j-n,i] is stored at A_indices[counter] in nonzeros(A) + nzA[A_indices[counter]] = val + end + #! format: on end - #! format: on end end return A diff --git a/src/result.jl b/src/result.jl index 66060a9f..e12c630e 100644 --- a/src/result.jl +++ b/src/result.jl @@ -894,6 +894,8 @@ struct TreeSetBicoloringResult{ A_indices::Vector{T} "storage for the edges of each tree in reverse BFS order" reverse_bfs_orders::Vector{Tuple{T,T}} + "whether each tree is a star and can be decompressed without substitution" + is_star::Vector{Bool} "internal storage for the positions of the trees inside `reverse_bfs_orders`" tree_edge_indices::Vector{T} "number of trees" @@ -914,7 +916,7 @@ function TreeSetBicoloringResult( symmetric_to_column::Vector{T}, decompression_eltype::Type{R}, ) where {T<:Integer,R} - (; reverse_bfs_orders, tree_edge_indices, nt) = tree_set + (; reverse_bfs_orders, is_star, tree_edge_indices, nt) = tree_set m, n = size(A) column_group = group_by_color(T, column_color) row_group = group_by_color(T, row_color) @@ -969,6 +971,7 @@ function TreeSetBicoloringResult( symmetric_to_row, A_indices, reverse_bfs_orders, + is_star, tree_edge_indices, nt, buffer, diff --git a/test/decompression.jl b/test/decompression.jl index 02ec5419..cb140f8a 100644 --- a/test/decompression.jl +++ b/test/decompression.jl @@ -187,3 +187,81 @@ end test_tree_csc_decompression(T, Ti, order, postprocessing) end end + +function test_tree_bicoloring_decompression(T, Ti, order, target) + problem = ColoringProblem(; structure=:nonsymmetric, partition=:bidirectional) + algorithm = GreedyColoringAlgorithm( + order; + decompression=:substitution, + postprocessing=target != :none, + postprocessing_minimizes=target == :none ? :all_colors : target, + ) + @testset "$pattern" for pattern in ( + :single_edge, :row_star, :column_star, :two_edges, :mixed, :empty + ) + A = spzeros(T, Ti, 8, 9) + edges = if pattern == :single_edge + [(2, 3)] + elseif pattern == :row_star + [(2, 1), (2, 3), (2, 5)] + elseif pattern == :column_star + [(1, 2), (3, 2), (5, 2)] + elseif pattern == :two_edges + [(2, 1), (2, 3)] + elseif pattern == :mixed + [(1, 1), (1, 2), (2, 2), (4, 4), (5, 5), (5, 6), (5, 7)] + else + Tuple{Int,Int}[] + end + for (i, j) in edges + A[i, j] = T <: Complex ? T(10i + j, i - j) : T(10i + j) + end + # Keep a structural zero, empty rows, and empty columns. + if !isempty(edges) + nonzeros(A)[end] = zero(T) + end + result = coloring(A, problem, algorithm; decompression_eltype=T) + sizes = diff(result.tree_edge_indices) + if pattern == :mixed + @test any(!, result.is_star) + @test any(result.is_star .& (sizes .> 1)) + @test any(==(1), sizes) + else + @test all(result.is_star) + if pattern in (:row_star, :column_star) + @test sizes == [3] + elseif pattern == :two_edges + @test sizes == [2] + elseif pattern == :single_edge + @test sizes == [1] + else + @test isempty(sizes) + end + end + Br, Bc = compress(A, result) + @testset "$(typeof(dest))" for dest in (copy(A), Matrix(A)) + for scale in (one(T), T(-2)) + rows, columns = scale .* Br, scale .* Bc + rows_view, columns_view = @view(rows[:, :]), @view(columns[:, :]) + fill!(dest isa SparseMatrixCSC ? nonzeros(dest) : dest, T(NaN)) + fill!(result.buffer, T(NaN)) + @test decompress!(dest, rows_view, columns_view, result) === dest + @test dest == scale .* A + if all(result.is_star) + @test all(isnan, result.buffer) + end + @test (@allocated decompress!(dest, rows_view, columns_view, result)) == 0 + end + end + end + return nothing +end + +@testset "Tree bicoloring decompression" begin + @testset "$T, $Ti, $order, target=$target" for T in (Float32, Float64, ComplexF64), + Ti in (Int32, Int64), order in (NaturalOrder(), LargestFirst()), + target in (:none, :all_colors, :row_colors, :column_colors) + + test_tree_bicoloring_decompression(T, Ti, order, target) + end +end From 715b6f274a1e5bf3732ce67fa753b1853157db64 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Beno=C3=AEt=20Legat?= Date: Mon, 26 Oct 2026 08:09:00 +0100 Subject: [PATCH 5/6] Skip substitution for stars in dense decompress! of TreeSetColoringResult --- src/decompression.jl | 61 +++++++++++++++++++++++++++++--------------- 1 file changed, 41 insertions(+), 20 deletions(-) diff --git a/src/decompression.jl b/src/decompression.jl index 8165d627..027d4a2f 100644 --- a/src/decompression.jl +++ b/src/decompression.jl @@ -559,9 +559,16 @@ end function decompress!( A::AbstractMatrix, B::AbstractMatrix, result::TreeSetColoringResult, uplo::Symbol=:F ) - (; ag, color, reverse_bfs_orders, tree_edge_indices, nt, diagonal_indices, buffer) = - result - (; S) = ag + (; + ag, + color, + reverse_bfs_orders, + is_star, + tree_edge_indices, + nt, + diagonal_indices, + buffer, + ) = result check_compatible_pattern(A, ag, uplo) R = eltype(A) fill!(A, zero(R)) @@ -585,25 +592,39 @@ function decompress!( first = tree_edge_indices[k] last = tree_edge_indices[k + 1] - 1 - # Reset the buffer to zero for all vertices in the tree (except the root) - for pos in first:last - (vertex, _) = reverse_bfs_orders[pos] - buffer_right_type[vertex] = zero(R) - end - # Reset the buffer to zero for the root vertex - (_, root) = reverse_bfs_orders[last] - buffer_right_type[root] = zero(R) + # A star needs no substitution or scratch-buffer updates. + if is_star[k] + for pos in first:last + (i, j) = reverse_bfs_orders[pos] + val = B[i, color[j]] + if in_triangle(i, j, uplo) + A[i, j] = val + end + if in_triangle(j, i, uplo) + A[j, i] = val + end + end + else + # Reset the buffer to zero for all vertices in the tree (except the root) + for pos in first:last + (vertex, _) = reverse_bfs_orders[pos] + buffer_right_type[vertex] = zero(R) + end + # Reset the buffer to zero for the root vertex + (_, root) = reverse_bfs_orders[last] + buffer_right_type[root] = zero(R) - for pos in first:last - (i, j) = reverse_bfs_orders[pos] - val = B[i, color[j]] - buffer_right_type[i] - buffer_right_type[j] = buffer_right_type[j] + val + for pos in first:last + (i, j) = reverse_bfs_orders[pos] + val = B[i, color[j]] - buffer_right_type[i] + buffer_right_type[j] = buffer_right_type[j] + val - if in_triangle(i, j, uplo) - A[i, j] = val - end - if in_triangle(j, i, uplo) - A[j, i] = val + if in_triangle(i, j, uplo) + A[i, j] = val + end + if in_triangle(j, i, uplo) + A[j, i] = val + end end end end From 9bb0e21e43a4aab89252293f16d01aece4cecd5c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Beno=C3=AEt=20Legat?= Date: Mon, 26 Oct 2026 08:09:06 +0100 Subject: [PATCH 6/6] Fix typos and clarify vertex offset comment in TreeSet --- src/coloring.jl | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/src/coloring.jl b/src/coloring.jl index 7a537c3e..22b1ea74 100644 --- a/src/coloring.jl +++ b/src/coloring.jl @@ -491,7 +491,7 @@ struct TreeSet{T} (BFS) traversal order of the `k`th tree starting from a depth-minimizing root. """ reverse_bfs_orders::Vector{Tuple{T,T}} - "For a tree index `1 <= k <= nt`, `is_star[k]` indicates whether the `k`th three is a star." + "For a tree index `1 <= k <= nt`, `is_star[k]` indicates whether the `k`th tree is a star." is_star::Vector{Bool} """ `tree_edge_indices[1]` is one and `tree_edge_indices[k+1] - tree_edge_indices[k]` is the number of edges in the `k`th tree. @@ -656,8 +656,9 @@ function TreeSet( queue_end = 0 # Positions of the first and last vertices in the current tree - # Note: tree_edge_indices contains the positions of the first and last edges, - # so we add to add an offset k-1 between edge indices and vertex indices + # Note: tree_edge_indices contains the positions of the first and last edges. + # Each tree has one more vertex than edges, so the k-1 previous trees shift + # the vertex positions by k-1 with respect to the edge positions. first_vertex = tree_edge_indices[k] + (k - 1) last_vertex = tree_edge_indices[k + 1] + (k - 1)