Merge pull request #15 from ArrogantGao/jg/polish-interface

Polish interfaces in OptimalBranchingCore

docs/src/index.md

@@ -65,7 +65,7 @@ SolverConfig
 
 
 # the result shows that the size of the maximum independent set is 9
-julia> reduce_and_branch(problem, config)
+julia> branch_and_reduce(problem, config)
 9
 
 # we can also use the EliminateGraphs package to verify the result
@@ -80,7 +80,7 @@ Furthermore, one can check the count of branches in the following way:
 julia> config = SolverConfig(MISReducer(), branching_strategy, MISCount)
 SolverConfig{MISReducer, BranchingStrategy{TensorNetworkSolver, IPSolver, EnvFilter, MinBoundarySelector, D3Measure}, MISCount}(MISReducer(), BranchingStrategy{TensorNetworkSolver, IPSolver, EnvFilter, MinBoundarySelector, D3Measure}(TensorNetworkSolver(), IPSolver(10), EnvFilter(), MinBoundarySelector(2), D3Measure()), MISCount)
 
-julia> reduce_and_branch(problem, config)
+julia> branch_and_reduce(problem, config)
 MISCount(9, 1)
 ```
 which shows that it takes only one branch to find the maximum independent set of size 9.

lib/OptimalBranchingCore/src/OptimalBranchingCore.jl

@@ -3,23 +3,32 @@ module OptimalBranchingCore
 using JuMP, HiGHS, SCIP
 using BitBasis
 
-export complexity_bv
-export Clause, BranchingTable, CandidateClause, DNF, Branch
-export BranchingStrategy
-export AbstractProblem, AbstractMeasure, AbstractReducer, AbstractSelector, AbstractTableSolver, AbstractSetCoverSolver
-export LPSolver, IPSolver
-export NoReducer
+# logic expressions
+export Clause, BranchingTable, DNF, booleans, ∨, ∧, ¬, covered_by, literals, is_true_literal, is_false_literal
+# weighted minimum set cover solvers and optimal branching rule
+export weighted_minimum_set_cover, AbstractSetCoverSolver, LPSolver, IPSolver
+export minimize_γ, optimal_branching_rule, OptimalBranchingResult
 
-export MaxSize, MaxSizeBranchCount
+##### interfaces #####
+# high-level interface
+export AbstractProblem, branch_and_reduce, BranchingStrategy
 
-export apply_branch, measure, reduce_problem, select, branching_table, weighted_minimum_set_cover
-export reduce_and_branch, optimal_branching_rule
+# variable selector interface
+export select_variable, AbstractSelector
+# branching table solver interface
+export branching_table, AbstractTableSolver
+# measure interface
+export measure, AbstractMeasure
+# reducer interface
+export reduce_problem, AbstractReducer, NoReducer
+# return type
+export MaxSize, MaxSizeBranchCount
 
 include("algebra.jl")
 include("bitbasis.jl")
 include("interfaces.jl")
-include("setcovering.jl")
 include("branching_table.jl")
+include("setcovering.jl")
 include("branch.jl")
 
 end

lib/OptimalBranchingCore/src/bitbasis.jl

@@ -54,6 +54,24 @@ function BitBasis.bdistance(c::Clause{INT}, b::INT) where INT <: Integer
     c1 = c.val & c.mask
     return bdistance(b1, c1)
 end
+"""
+    literals(c::Clause)
+
+Return all literals in the clause.
+"""
+literals(c::Clause) = [Clause(readbit(c.mask, i), readbit(c.val, i)) for i=1:bsizeof(c.mask) if readbit(c.mask, i) == 1]
+"""
+    is_true_literal(c::Clause)
+
+Check if the clause is a true literal.
+"""
+is_true_literal(c::Clause) = count_ones(c.mask) == 1 && all(i->readbit(c.val, i) == readbit(c.mask, i), 1:bsizeof(c.mask))
+"""
+    is_false_literal(c::Clause)
+
+Check if the clause is a false literal.
+"""
+is_false_literal(c::Clause) = count_ones(c.mask) == 1 && iszero(c.val)
 
 # Flip all bits in `b`, `n` is the number of bits
 function flip_all(n::Int, b::INT) where INT <: Integer

lib/OptimalBranchingCore/src/branch.jl

@@ -1,5 +1,5 @@
 """
-    optimal_branching_rule(table::BranchingTable, variables::Vector, problem::AbstractProblem, measure::AbstractMeasure, solver::AbstractSetCoverSolver)::DNF
+    optimal_branching_rule(table::BranchingTable, variables::Vector, problem::AbstractProblem, measure::AbstractMeasure, solver::AbstractSetCoverSolver)::OptimalBranchingResult
 
 Generate an optimal branching rule from a given branching table.
 
@@ -11,64 +11,12 @@ Generate an optimal branching rule from a given branching table.
 - `solver`: The solver used for the weighted minimum set cover problem, which can be either [`LPSolver`](@ref) or [`IPSolver`](@ref).
 
 ### Returns
-A [`DNF`](@ref) object representing the optimal branching rule.
+A [`OptimalBranchingResult`](@ref) object representing the optimal branching rule.
 """
 function optimal_branching_rule(table::BranchingTable, variables::Vector, problem::AbstractProblem, m::AbstractMeasure, solver::AbstractSetCoverSolver)
-    candidates = candidate_clauses(table)
-    size_reductions = [measure(problem, m) - measure(first(apply_branch(problem, candidate.clause, variables)), m) for candidate in candidates]
-    selection, _ = minimize_γ(length(table.table), candidates, size_reductions, solver; γ0=2.0)
-    return DNF(map(i->candidates[i].clause, selection))
-end
-
-# TODO: we need to extend this function to trim the candidate clauses
-"""
-    candidate_clauses(tbl::BranchingTable{INT}) where {INT}
-
-Generates candidate clauses from a branching table.
-
-### Arguments
-- `tbl::BranchingTable{INT}`: The branching table containing bit strings.
-
-### Returns
-- `Vector{CandidateClause{INT}}`: A vector of `CandidateClause` objects generated from the branching table.
-"""
-function candidate_clauses(tbl::BranchingTable{INT}) where {INT}
-    n, bss = tbl.bit_length, tbl.table
-    bs = vcat(bss...)
-    all_clauses = Set{Clause{INT}}()
-    temp_clauses = [Clause(bmask(INT, 1:n), bs[i]) for i in 1:length(bs)]
-    while !isempty(temp_clauses)
-        c = pop!(temp_clauses)
-        if !(c in all_clauses)
-            push!(all_clauses, c)
-            idc = Set(covered_items(bss, c))
-            for i in 1:length(bss)
-                if i ∉ idc                
-                    for b in bss[i]
-                        c_new = gather2(n, c, Clause(bmask(INT, 1:n), b))
-                        if (c_new != c) && c_new.mask != 0
-                            push!(temp_clauses, c_new)
-                        end
-                    end
-                end
-            end
-        end
-    end
-
-    allcovers = [CandidateClause(covered_items(bss, c), c) for c in all_clauses]
-    return allcovers
-end
-# Returns the indices of the bit strings that are covered by the clause.
-function covered_items(bitstrings, clause::Clause)
-    return findall(bs -> any(x->covered_by(x, clause), bs), bitstrings)
-end
-# merge two clauses, i.e. generate a new clause covering both
-function gather2(n::Int, c1::Clause{INT}, c2::Clause{INT}) where INT
-    b1 = c1.val & c1.mask
-    b2 = c2.val & c2.mask
-    mask = (b1 ⊻ flip_all(n, b2)) & c1.mask & c2.mask
-    val = b1 & mask
-    return Clause(mask, val)
+    candidates = collect(candidate_clauses(table))
+    size_reductions = [measure(problem, m) - measure(first(apply_branch(problem, candidate, variables)), m) for candidate in candidates]
+    return minimize_γ(table, candidates, size_reductions, solver; γ0=2.0)
 end
 
 """
@@ -99,7 +47,7 @@ BranchingStrategy
 """)
 
 """
-    reduce_and_branch(problem::AbstractProblem, config::BranchingStrategy; reducer::AbstractReducer=NoReducer(), result_type=Int)
+    branch_and_reduce(problem::AbstractProblem, config::BranchingStrategy; reducer::AbstractReducer=NoReducer(), result_type=Int)
 
 Branch the given problem using the specified solver configuration.
 
@@ -114,18 +62,18 @@ Branch the given problem using the specified solver configuration.
 ### Returns
 The resulting value, which may have different type depending on the `result_type`.
 """
-function reduce_and_branch(problem::AbstractProblem, config::BranchingStrategy, reducer::AbstractReducer, result_type)
+function branch_and_reduce(problem::AbstractProblem, config::BranchingStrategy, reducer::AbstractReducer, result_type)
     isempty(problem) && return zero(result_type)
     # reduce the problem
     rp, reducedvalue = reduce_problem(result_type, problem, reducer)
-    rp !== problem && return reduce_and_branch(rp, config, reducer, result_type) * reducedvalue
+    rp !== problem && return branch_and_reduce(rp, config, reducer, result_type) * reducedvalue
 
     # branch the problem
     variables = select_variables(rp, config.measure, config.selector)  # select a subset of variables
     tbl = branching_table(rp, config.table_solver, variables)      # compute the BranchingTable
-    rule = optimal_branching_rule(tbl, variables, rp, config.measure, config.set_cover_solver)  # compute the optimal branching rule
-    return sum(rule.clauses) do branch  # branch and recurse
+    result = optimal_branching_rule(tbl, variables, rp, config.measure, config.set_cover_solver)  # compute the optimal branching rule
+    return sum(result.optimal_rule.clauses) do branch  # branch and recurse
         subproblem, localvalue = apply_branch(rp, branch, variables)
-        reduce_and_branch(subproblem, config, reducer, result_type) * result_type(localvalue) * reducedvalue
+        branch_and_reduce(subproblem, config, reducer, result_type) * result_type(localvalue) * reducedvalue
     end
 end

lib/OptimalBranchingCore/src/setcovering.jl

@@ -92,64 +92,136 @@ function bisect_solve(f, a, fa, b, fb)
 end
 
 """
-    minimize_γ(candidate_clauses::AbstractVector{CandidateClause{INT}}, Δρ::Vector{TF}, solver) where{INT, TF}
+    OptimalBranchingResult{INT <: Integer}
+
+The result type for the optimal branching rule.
+
+### Fields
+- `selected_ids::Vector{Int}`: The indices of the selected rows in the branching table.
+- `optimal_rule::DNF{INT}`: The optimal branching rule.
+- `branching_vector::Vector{T<:Real}`: The branching vector that records the size reduction in each subproblem.
+- `γ::Float64`: The optimal γ value (the complexity of the branching rule).
+"""
+struct OptimalBranchingResult{INT <: Integer, T <: Real}
+    selected_ids::Vector{Int}
+    optimal_rule::DNF{INT}
+    branching_vector::Vector{T}
+    γ::Float64
+end
+
+"""
+    minimize_γ(table::BranchingTable, candidates::Vector{Clause}, Δρ::Vector, solver)
 
 Finds the optimal cover based on the provided vector of problem size reduction.
 This function implements a cover selection algorithm using an iterative process.
 It utilizes an integer programming solver to optimize the selection of sub-covers based on their complexity.
 
 ### Arguments
-- `candidate_clauses::AbstractVector{CandidateClause{INT}}`: A vector of CandidateClause structures.
-- `Δρ::Vector{TF}`: A vector of problem size reduction for each CandidateClause.
+- `table::BranchingTable`: A branching table containing clauses that need to be covered, a table entry is covered by a clause if one of its bit strings satisfies the clause. Please refer to [`covered_by`](@ref) for more details.
+- `candidates::Vector{Clause}`: A vector of candidate clauses to form the branching rule (in the form of [`DNF`](@ref)).
+- `Δρ::Vector`: A vector of problem size reduction for each candidate clause.
 - `solver`: The solver to be used. It can be an instance of `LPSolver` or `IPSolver`.
 
 ### Keyword Arguments
 - `γ0::Float64`: The initial γ value.
 
 ### Returns
-A tuple of two elements: (indices of selected clauses, γ)
+A tuple of two elements: (indices of selected subsets, γ)
 """
-function minimize_γ(num_items::Int, candidate_clauses::AbstractVector{CandidateClause{INT}}, Δρ::Vector{TF}, solver::AbstractSetCoverSolver; γ0::Float64 = 2.0) where{INT, TF}
-    @debug "solver = $(solver), sets = $(candidate_clauses), γ0 = $γ0"
+function minimize_γ(table::BranchingTable, candidates::Vector{Clause{INT}}, Δρ::Vector, solver::AbstractSetCoverSolver; γ0::Float64 = 2.0) where {INT}
+    @debug "solver = $(solver), subsets = $(subsets), γ0 = $γ0"
+    subsets = [covered_items(table.table, c) for c in candidates]
+    num_items = length(table.table)
 
     # Note: the following instance is captured for time saving, and also for it may cause IP solver to fail
-    for (k, clause) in enumerate(candidate_clauses)
-        (length(clause.covered_items) == num_items) && return [k], 1.0
+    for (k, subset) in enumerate(subsets)
+        (length(subset) == num_items) && return OptimalBranchingResult([k], DNF([candidates[k]]), [Δρ[k]], 1.0)
     end
 
     cx_old = cx = γ0
     local picked_scs
     for i = 1:solver.max_itr
         weights = 1 ./ cx_old .^ Δρ
-        picked_scs = weighted_minimum_set_cover(solver, weights, candidate_clauses, num_items)
+        picked_scs = weighted_minimum_set_cover(solver, weights, subsets, num_items)
         cx = complexity_bv(Δρ[picked_scs])
-        @debug "Iteration $i, picked indices = $(picked_scs), clauses = $(candidate_clauses[picked_scs]), branching_vector = $(Δρ[picked_scs]), γ = $cx"
+        @debug "Iteration $i, picked indices = $(picked_scs), subsets = $(subsets[picked_scs]), branching_vector = $(Δρ[picked_scs]), γ = $cx"
         cx ≈ cx_old && break  # convergence
         cx_old = cx
     end
-    return picked_scs, cx
+    return OptimalBranchingResult(picked_scs, DNF([candidates[i] for i in picked_scs]), Δρ[picked_scs], cx)
 end
 
+# TODO: we need to extend this function to trim the candidate clauses
+"""
+    candidate_clauses(tbl::BranchingTable{INT}) where {INT}
+
+Generates candidate clauses from a branching table.
+
+### Arguments
+- `tbl::BranchingTable{INT}`: The branching table containing bit strings.
+
+### Returns
+- `Vector{Clause{INT}}`: A vector of `Clause` objects generated from the branching table.
+"""
+function candidate_clauses(tbl::BranchingTable{INT}) where {INT}
+    n, bss = tbl.bit_length, tbl.table
+    bs = vcat(bss...)
+    all_clauses = Set{Clause{INT}}()
+    temp_clauses = [Clause(bmask(INT, 1:n), bs[i]) for i in 1:length(bs)]
+    while !isempty(temp_clauses)
+        c = pop!(temp_clauses)
+        if !(c in all_clauses)
+            push!(all_clauses, c)
+            idc = Set(covered_items(bss, c))
+            for i in 1:length(bss)
+                if i ∉ idc                
+                    for b in bss[i]
+                        c_new = gather2(n, c, Clause(bmask(INT, 1:n), b))
+                        if (c_new != c) && c_new.mask != 0
+                            push!(temp_clauses, c_new)
+                        end
+                    end
+                end
+            end
+        end
+    end
+    return all_clauses
+end
+# Returns the indices of the bit strings that are covered by the clause.
+function covered_items(bitstrings, clause::Clause)
+    return findall(bs -> any(x->covered_by(x, clause), bs), bitstrings)
+end
+# merge two clauses, i.e. generate a new clause covering both
+function gather2(n::Int, c1::Clause{INT}, c2::Clause{INT}) where INT
+    b1 = c1.val & c1.mask
+    b2 = c2.val & c2.mask
+    mask = (b1 ⊻ flip_all(n, b2)) & c1.mask & c2.mask
+    val = b1 & mask
+    return Clause(mask, val)
+end
+
+
+
 """
-    weighted_minimum_set_cover(solver, weights::AbstractVector, candidate_clauses::AbstractVector{CandidateClause{INT}}, num_items::Int) where{INT, TF, T}
+    weighted_minimum_set_cover(solver, weights::AbstractVector, subsets::Vector{Vector{Int}}, num_items::Int)
 
 Solves the weighted minimum set cover problem.
 
 ### Arguments
 - `solver`: The solver to be used. It can be an instance of `LPSolver` or `IPSolver`.
-- `weights::AbstractVector`: The weights of the candidate clauses.
-- `candidate_clauses::AbstractVector{CandidateClause{INT}}`: A vector of CandidateClause structures.
+- `weights::AbstractVector`: The weights of the subsets.
+- `subsets::Vector{Vector{Int}}`: A vector of subsets.
 - `num_items::Int`: The number of elements to cover.
 
 ### Returns
-A vector of indices of selected clauses.
+A vector of indices of selected subsets.
 """
-function weighted_minimum_set_cover(solver::LPSolver, weights::AbstractVector, candidate_clauses::AbstractVector{CandidateClause{INT}}, num_items::Int) where{INT}
-    nsc = length(candidate_clauses)
+function weighted_minimum_set_cover(solver::LPSolver, weights::AbstractVector, subsets::Vector{Vector{Int}}, num_items::Int)
+    nsc = length(subsets)
 
     sets_id = [Vector{Int}() for _=1:num_items]
     for i in 1:nsc
-        for j in candidate_clauses[i].covered_items
+        for j in subsets[i]
             push!(sets_id[j], i)
         end
     end
@@ -165,15 +237,15 @@ function weighted_minimum_set_cover(solver::LPSolver, weights::AbstractVector, c
 
     optimize!(model)
     @assert is_solved_and_feasible(model)
-    return pick_sets(value.(x), candidate_clauses, num_items)
+    return pick_sets(value.(x), subsets, num_items)
 end
 
-function weighted_minimum_set_cover(solver::IPSolver, weights::AbstractVector, candidate_clauses::AbstractVector{CandidateClause{INT}}, num_items::Int) where{INT}
-    nsc = length(candidate_clauses)
+function weighted_minimum_set_cover(solver::IPSolver, weights::AbstractVector, subsets::Vector{Vector{Int}}, num_items::Int)
+    nsc = length(subsets)
 
     sets_id = [Vector{Int}() for _=1:num_items]
     for i in 1:nsc
-        for j in candidate_clauses[i].covered_items
+        for j in subsets[i]
             push!(sets_id[j], i)
         end
     end
@@ -191,20 +263,20 @@ function weighted_minimum_set_cover(solver::IPSolver, weights::AbstractVector, c
 
     optimize!(model)
     @assert is_solved_and_feasible(model)
-    return pick_sets(value.(x), candidate_clauses, num_items)
+    return pick_sets(value.(x), subsets, num_items)
 end
 
 # by viewing xs as the probability of being selected, we can use a random algorithm to pick the sets
-function pick_sets(xs::Vector{TF}, candidate_clauses::AbstractVector{CandidateClause{INT}}, num_items::Int) where{INT, TF}
+function pick_sets(xs::Vector, subsets::Vector{Vector{Int}}, num_items::Int)
     picked = Set{Int}()
     picked_ids = Set{Int}()
-    nsc = length(candidate_clauses)
+    nsc = length(subsets)
     flag = true
     while flag 
         for i in 1:nsc
             if (rand() < xs[i]) && i ∉ picked
                 push!(picked, i)
-                picked_ids = union(picked_ids, candidate_clauses[i].covered_items)
+                picked_ids = union(picked_ids, subsets[i])
             end
             if length(picked_ids) == num_items
                 flag = false