Add test functions substitutable_columns and substitutable_bidirectional

docs/src/dev.md

@@ -50,6 +50,8 @@ SparseMatrixColorings.directly_recoverable_columns
 SparseMatrixColorings.symmetrically_orthogonal_columns
 SparseMatrixColorings.structurally_orthogonal_columns
 SparseMatrixColorings.structurally_biorthogonal
+SparseMatrixColorings.substitutable_columns
+SparseMatrixColorings.substitutable_bidirectional
 SparseMatrixColorings.valid_dynamic_order
 ```
 

src/check.jl

@@ -86,11 +86,15 @@ A partition of the columns of a symmetric matrix `A` is _symmetrically orthogona
 1. the group containing the column `A[:, j]` has no other column with a nonzero in row `i`
 2. the group containing the column `A[:, i]` has no other column with a nonzero in row `j`
 
+It is equivalent to a __star coloring__.
+
 !!! warning
     This function is not coded with efficiency in mind, it is designed for small-scale tests.
 
 # References
 
+> [_On the Estimation of Sparse Hessian Matrices_](https://doi.org/10.1137/0716078), Powell and Toint (1979)
+> [_Estimation of sparse hessian matrices and graph coloring problems_](https://doi.org/10.1007/BF02612334), Coleman and Moré (1984)
 > [_What Color Is Your Jacobian? Graph Coloring for Computing Derivatives_](https://epubs.siam.org/doi/10.1137/S0036144504444711), Gebremedhin et al. (2005)
 """
 function symmetrically_orthogonal_columns(
@@ -102,7 +106,7 @@ function symmetrically_orthogonal_columns(
     end
     issymmetric(A) || return false
     group = group_by_color(color)
-    for i in axes(A, 2), j in axes(A, 2)
+    for i in axes(A, 1), j in axes(A, 2)
         iszero(A[i, j]) && continue
         ci, cj = color[i], color[j]
         check = _bilateral_check(
@@ -126,6 +130,8 @@ A bipartition of the rows and columns of a matrix `A` is _structurally biorthogo
 1. the group containing the column `A[:, j]` has no other column with a nonzero in row `i`
 2. the group containing the row `A[i, :]` has no other row with a nonzero in column `j`
 
+It is equivalent to an __star bicoloring__.
+
 !!! warning
     This function is not coded with efficiency in mind, it is designed for small-scale tests.
 """
@@ -138,8 +144,8 @@ function structurally_biorthogonal(
     if !proper_length_bicoloring(A, row_color, column_color; verbose)
         return false
     end
-    column_group = group_by_color(column_color)
     row_group = group_by_color(row_color)
+    column_group = group_by_color(column_color)
     for i in axes(A, 1), j in axes(A, 2)
         iszero(A[i, j]) && continue
         ci, cj = row_color[i], column_color[j]
@@ -261,6 +267,174 @@ function directly_recoverable_columns(
     return true
 end
 
+"""
+    substitutable_columns(
+        A::AbstractMatrix, order_nonzeros::AbstractMatrix, color::AbstractVector{<:Integer};
+        verbose=false
+    )
+
+Return `true` if coloring the columns of the symmetric matrix `A` with the vector `color` results in a partition that is substitutable, and `false` otherwise.
+For all nonzeros `A[i, j]`, `order_nonzeros[i, j]` provides its order of recovery.
+
+A partition of the columns of a symmetric matrix `A` is _substitutable_ if, for every nonzero element `A[i, j]`, either of the following statements holds:
+
+1. the group containing the column `A[:, j]` has all nonzeros in row `i` ordered before `A[i, j]`
+2. the group containing the column `A[:, i]` has all nonzeros in row `j` ordered before `A[i, j]`
+
+It is equivalent to an __acyclic coloring__.
+
+!!! warning
+    This function is not coded with efficiency in mind, it is designed for small-scale tests.
+
+# References
+
+> [_On the Estimation of Sparse Hessian Matrices_](https://doi.org/10.1137/0716078), Powell and Toint (1979)
+> [_The Cyclic Coloring Problem and Estimation of Sparse Hessian Matrices_](https://doi.org/10.1137/0607026), Coleman and Cai (1986)
+> [_What Color Is Your Jacobian? Graph Coloring for Computing Derivatives_](https://epubs.siam.org/doi/10.1137/S0036144504444711), Gebremedhin et al. (2005)
+"""
+function substitutable_columns(
+    A::AbstractMatrix,
+    order_nonzeros::AbstractMatrix,
+    color::AbstractVector{<:Integer};
+    verbose::Bool=false,
+)
+    checksquare(A)
+    if !proper_length_coloring(A, color; verbose)
+        return false
+    end
+    issymmetric(A) || return false
+    group = group_by_color(color)
+    for i in axes(A, 1), j in axes(A, 2)
+        iszero(A[i, j]) && continue
+        ci, cj = color[i], color[j]
+        check = _substitutable_check(
+            A, order_nonzeros; i, j, ci, cj, row_group=group, column_group=group, verbose
+        )
+        !check && return false
+    end
+    return true
+end
+
+"""
+    substitutable_bidirectional(
+        A::AbstractMatrix, order_nonzeros::AbstractMatrix, row_color::AbstractVector{<:Integer}, column_color::AbstractVector{<:Integer};
+        verbose=false
+    )
+
+Return `true` if bicoloring of the matrix `A` with the vectors `row_color` and `column_color` results in a bipartition that is substitutable, and `false` otherwise.
+For all nonzeros `A[i, j]`, `order_nonzeros[i, j]` provides its order of recovery.
+
+A bipartition of the rows and columns of a matrix `A` is _substitutable_ if, for every nonzero element `A[i, j]`, either of the following statements holds:
+
+1. the group containing the column `A[:, j]` has all nonzeros in row `i` ordered before `A[i, j]`
+2. the group containing the row `A[i, :]` has all nonzeros in column `j` ordered before `A[i, j]`
+
+It is equivalent to an __acyclic bicoloring__.
+
+!!! warning
+    This function is not coded with efficiency in mind, it is designed for small-scale tests.
+"""
+function substitutable_bidirectional(
+    A::AbstractMatrix,
+    order_nonzeros::AbstractMatrix,
+    row_color::AbstractVector{<:Integer},
+    column_color::AbstractVector{<:Integer};
+    verbose::Bool=false,
+)
+    if !proper_length_bicoloring(A, row_color, column_color; verbose)
+        return false
+    end
+    row_group = group_by_color(row_color)
+    column_group = group_by_color(column_color)
+    for i in axes(A, 1), j in axes(A, 2)
+        iszero(A[i, j]) && continue
+        ci, cj = row_color[i], column_color[j]
+        check = _substitutable_check(
+            A, order_nonzeros; i, j, ci, cj, row_group, column_group, verbose
+        )
+        !check && return false
+    end
+    return true
+end
+
+function _substitutable_check(
+    A::AbstractMatrix,
+    order_nonzeros::AbstractMatrix;
+    i::Integer,
+    j::Integer,
+    ci::Integer,
+    cj::Integer,
+    row_group::AbstractVector,
+    column_group::AbstractVector,
+    verbose::Bool,
+)
+    order_ij = order_nonzeros[i, j]
+    k_row = 0
+    k_column = 0
+    if ci != 0
+        for k in row_group[ci]
+            (k == i) && continue
+            if !iszero(A[k, j])
+                order_kj = order_nonzeros[k, j]
+                @assert !iszero(order_kj)
+                if order_kj > order_ij
+                    k_row = k
+                end
+            end
+        end
+    end
+    if cj != 0
+        for k in column_group[cj]
+            (k == j) && continue
+            if !iszero(A[i, k])
+                order_ik = order_nonzeros[i, k]
+                @assert !iszero(order_ik)
+                if order_ik > order_ij
+                    k_column = k
+                end
+            end
+        end
+    end
+    if ci == 0 && cj == 0
+        if verbose
+            @warn """
+            For coefficient (i=$i, j=$j) with colors (ci=$ci, cj=$cj):
+            - Row color ci=$ci is neutral.
+            - Column color cj=$cj is neutral.
+            """
+        end
+        return false
+    elseif ci == 0 && !iszero(k_column)
+        if verbose
+            @warn """
+            For coefficient (i=$i, j=$j) with colors (ci=$ci, cj=$cj):
+            - Row color ci=$ci is neutral.
+            - For the column $k_column in column color cj=$cj, A[$i, $k_column] is ordered after A[$i, $j].
+            """
+        end
+        return false
+    elseif cj == 0 && !iszero(k_row)
+        if verbose
+            @warn """
+            For coefficient (i=$i, j=$j) with colors (ci=$ci, cj=$cj):
+            - For the row $k_row in row color ci=$ci, A[$k_row, $j] is ordered after A[$i, $j].
+            - Column color cj=$cj is neutral.
+            """
+        end
+        return false
+    elseif !iszero(k_row) && !iszero(k_column)
+        if verbose
+            @warn """
+            For coefficient (i=$i, j=$j) with colors (ci=$ci, cj=$cj):
+            - For the row $k_row in row color ci=$ci, A[$k_row, $j] is ordered after A[$i, $j].
+            - For the column $k_column in column color cj=$cj, A[$i, $k_column] is ordered after A[$i, $j].
+            """
+        end
+        return false
+    end
+    return true
+end
+
 """
     valid_dynamic_order(g::AdjacencyGraph, π::AbstractVector{<:Integer}, order::DynamicDegreeBasedOrder)
     valid_dynamic_order(bg::AdjacencyGraph, ::Val{side}, π::AbstractVector{<:Integer}, order::DynamicDegreeBasedOrder)

test/check.jl

@@ -4,6 +4,8 @@ using SparseMatrixColorings:
     symmetrically_orthogonal_columns,
     structurally_biorthogonal,
     directly_recoverable_columns,
+    substitutable_columns,
+    substitutable_bidirectional,
     what_fig_41,
     efficient_fig_1
 using Test
@@ -184,3 +186,125 @@ For coefficient (i=1, j=1) with colors (ci=0, cj=0):
 """,
     ) structurally_biorthogonal(A, [0, 2, 2, 3], [0, 2, 2, 2, 3], verbose=true)
 end
+
+@testset "Substitutable columns" begin
+    A1 = [
+        1 1 1 1 1
+        1 1 0 0 0
+        1 0 1 0 0
+        1 0 0 1 0
+        1 0 0 0 1
+    ]
+    B1 = [
+        1 6 7 8 9
+        6 2 0 0 0
+        7 0 3 0 0
+        8 0 0 4 0
+        9 0 0 0 5
+    ]
+    A2 = [
+        1 1 0 0 0
+        1 1 1 0 0
+        0 1 1 1 0
+        0 0 1 1 1
+        0 0 0 1 1
+    ]
+    B2 = [
+        5 1 0 0 0
+        1 6 2 0 0
+        0 2 7 3 0
+        0 0 3 8 4
+        0 0 0 4 9
+    ]
+    A3 = [
+        0 1 1 1 1
+        1 0 1 1 1
+        1 1 0 1 1
+        1 1 1 0 1
+        1 1 1 1 0
+    ]
+    B3 = [
+        0 1 2 3 4
+        1 0 5 6 7
+        2 5 0 8 9
+        3 6 8 0 10
+        4 7 9 10 0
+    ]
+
+    # success
+
+    substitutable_columns(A1, B1, [1, 2, 2, 2, 2])
+    substitutable_columns(A2, B2, [1, 2, 3, 1, 2])
+    substitutable_columns(A3, B3, [1, 2, 3, 4, 0])
+
+    # failure
+
+    @test !substitutable_columns(A1, B1, [1, 1, 1, 1])
+    log = (:warn, "4 colors provided for 5 columns.")
+    @test_logs log substitutable_columns(A1, B1, [1, 1, 1, 1]; verbose=true)
+
+    @test !substitutable_columns(A1, B1, [1, 1, 1, 1, 1])
+    @test_logs (
+        :warn,
+        """
+For coefficient (i=1, j=1) with colors (ci=1, cj=1):
+- For the row 5 in row color ci=1, A[5, 1] is ordered after A[1, 1].
+- For the column 5 in column color cj=1, A[1, 5] is ordered after A[1, 1].
+""",
+    ) substitutable_columns(A1, B1, [1, 1, 1, 1, 1]; verbose=true)
+
+    @test !substitutable_columns(A2, B2, [1, 2, 0, 1, 2])
+    @test_logs (
+        :warn,
+        """
+For coefficient (i=3, j=3) with colors (ci=0, cj=0):
+- Row color ci=0 is neutral.
+- Column color cj=0 is neutral.
+""",
+    ) substitutable_columns(A2, B2, [1, 2, 0, 1, 2]; verbose=true)
+
+    @test !substitutable_columns(A3, B3, [0, 1, 2, 3, 3])
+    @test_logs (
+        :warn,
+        """
+For coefficient (i=1, j=4) with colors (ci=0, cj=3):
+- Row color ci=0 is neutral.
+- For the column 5 in column color cj=3, A[1, 5] is ordered after A[1, 4].
+""",
+    ) substitutable_columns(A3, B3, [0, 1, 2, 3, 3]; verbose=true)
+
+    @test !substitutable_columns(A3, B3, [1, 2, 3, 3, 0])
+    @test_logs (
+        :warn,
+        """
+For coefficient (i=3, j=5) with colors (ci=3, cj=0):
+- For the row 4 in row color ci=3, A[4, 5] is ordered after A[3, 5].
+- Column color cj=0 is neutral.
+""",
+    ) substitutable_columns(A3, B3, [1, 2, 3, 3, 0]; verbose=true)
+end
+
+@testset "Substitutable bidirectional" begin
+    A = [
+        1 0 0
+        0 1 0
+        0 0 1
+    ]
+    B = [
+        1 0 0
+        0 2 0
+        0 0 3
+    ]
+
+    # success
+
+    substitutable_bidirectional(A, B, [1, 0, 0], [0, 1, 1])
+
+    # failure
+
+    log = (:warn, "2 colors provided for 3 columns.")
+    @test_logs log !substitutable_bidirectional(A, B, [1, 0, 0], [0, 1]; verbose=true)
+
+    log = (:warn, "4 colors provided for 3 rows.")
+    @test_logs log !substitutable_bidirectional(A, B, [1, 0, 0, 1], [0, 1, 1]; verbose=true)
+end