search

lanl-ansi / MathProgIncidence.jl

Tools for constructing and analyzing the incidence graph or matrix of variables and constraints in a JuMP model
Other
9 stars 1 forks source link

Performance issue in Dulmage-Mendelsohn #8

Closed Robbybp closed 1 year ago

Robbybp commented 1 year ago

dulmage_mendelsohn.jl uses the following:

filter = cat(vec1, vec2, ...; dims = 1)
other = [n for n in nodes if !(n in filter)]

As filter is a vector, this is quadratic time and will become slow if nodes, the vector of nodes in the bipartite graph, becomes large. As we do expect this set of nodes to become large, we should implement filter as a set.

Robbybp commented 1 year ago

This is not the current performance bottleneck. Testing on some medium-sized instances from pglib, a 2k-bus ACOPF mode takes 70 s to compute the DM partition. In a cursory test, the dominant computations cost appears to be in computing the maximum matching. However, BipartiteMatching.jl reports excellent performance on much larger graphs than this. The bottlneck may be translating our Graphs.jl graph into the format required by BipartiteMatching.jl?

Robbybp commented 1 year ago

BipartiteMatching.jl accepts as input a 2-dimensional bit array. The quadratic loop over both edge sets used to construct this array is likely the bottleneck. To get around this, we need an implementation that takes a sparse graph or matrix as an input.

Robbybp commented 1 year ago

The BlossomV url appears to be back up, so I can test with the GraphsMatching Hungarian algorithm implementation:

import Graphs
import GraphsMatching as GM

function maximum_matching(graph::Graphs.Graph, set1::Set)
    if !_is_valid_bipartition(graph, set1)
        throw(Exception)
    end 
    nvert = Graphs.nv(graph)
    weights = SparseArrays.spzeros(nvert, nvert)
    for e in Graphs.edges(graph)
        weights[Graphs.src(e), Graphs.dst(e)] = 1.0 
    end 
    println("Beginning maximum weight maximal matching")
    result = GM.maximum_weight_maximal_matching(
        graph,
        weights;
        algorithm = GM.HungarianAlgorithm(),
    )   
    println("Done with maximum weight maximal matching")
    matching = Dict(
        # The GraphsMatching convention is that mate[n] is -1 if n is unmatched.
        # Calling functions need a map from set1 nodes to set2 (other) nodes.
        n1 => result.mate[n1] for n1 in set1 if result.mate[n1] != -1
    )   
    return matching
end

Despite accepting sparse data structures, this is actually significantly slower than the BipartiteMatching.jl implementation. The time complexity of the Hungarian algorithm is O(n^3), although I'm not sure if this is O(nm) or a true O(n^3). We may need a Hopcroft-Karm implementation?

Robbybp commented 1 year ago

The GraphsMatching LPAlgorithm option with HiGHS as the solver is faster than HungarianAlgorithm (computes DM for IEEE-118 in 12 s as opposed to 17 s), although it seems to not be reliable. My current test suite fails with an infeasible model sent to HiGHS at some point.

Both of these are significantly slower than BipartiteMatching (DM for IEEE-118 in 0.7 s).

Robbybp commented 1 year ago

Timing results on small PGLib instances, with an initial implementation of Hopcroft-Karp in Graphs.jl.

Case Graph build time Dulmage-Mendelsohn time N. variables
pglib_opf_case3_lmbd.m 2.96e-04 3.07e-04 19
pglib_opf_case5_pjm.m 4.28e-04 4.40e-04 35
pglib_opf_case14_ieee.m 6.60e-04 8.38e-04 109
pglib_opf_case24_ieee_rts.m 1.06e-03 2.17e-03 201
pglib_opf_case30_as.m 1.07e-03 2.18e-03 225
pglib_opf_case30_ieee.m 1.06e-03 2.14e-03 225
pglib_opf_case39_epri.m 1.05e-03 2.72e-03 263
pglib_opf_case57_ieee.m 1.73e-03 6.40e-03 435
pglib_opf_case60_c.m 1.84e-03 7.44e-03 473
pglib_opf_case73_ieee_rts.m 2.51e-03 1.48e-02 627
pglib_opf_case118_ieee.m 4.67e-01 8.59e-01 981
pglib_opf_case89_pegase.m 3.81e-03 4.64e-02 1019
pglib_opf_case200_activ.m 4.95e-03 7.05e-02 1381
pglib_opf_case179_goc.m 5.05e-03 6.47e-02 1411
pglib_opf_case162_ieee_dtc.m 5.22e-03 6.77e-02 1461
pglib_opf_case300_ieee.m 7.93e-03 1.63e-01 2245
pglib_opf_case240_pserc.m 8.22e-03 1.66e-01 2273
pglib_opf_case500_goc.m 1.35e-02 4.68e-01 3913
pglib_opf_case588_sdet.m 1.34e-02 4.59e-01 3921
pglib_opf_case793_goc.m 1.80e-02 7.87e-01 5239
pglib_opf_case1354_pegase.m 5.22e-02 3.16e+00 10673
pglib_opf_case1888_rte.m 6.75e-02 5.43e+00 13901
pglib_opf_case1951_rte.m 6.80e-02 5.84e+00 14287
pglib_opf_case2000_goc.m 6.73e-02 9.79e+00 18533
pglib_opf_case2736sp_k.m 8.83e-02 1.01e+01 18549
pglib_opf_case2737sop_k.m 6.65e-02 9.86e+00 18551
pglib_opf_case2746wp_k.m 6.76e-02 1.04e+01 18609
pglib_opf_case2746wop_k.m 6.70e-02 1.04e+01 18721
pglib_opf_case2848_rte.m 7.75e-02 1.28e+01 20801
pglib_opf_case2868_rte.m 7.62e-02 1.29e+01 20969
pglib_opf_case2853_sdet.m 8.11e-02 1.39e+01 21391
pglib_opf_case2869_pegase.m 9.01e-02 1.69e+01 24067
pglib_opf_case2742_goc.m 8.60e-02 1.67e+01 24177

The IEEE-118 case is comparable to BipartiteMatching. On GOC-2742, BipartiteMatching takes 37.4 s.

Robbybp commented 1 year ago

The time spent computing the matching seems small. HK on GOC-2742 takes 0.15 s. The rest of the time is probably the quadratic loop mentioned above.

When I fix this quadratic loop, GOC-2742 takes about 1 s.

Robbybp commented 1 year ago

Results on all PGLib instances with HK matching and fixed quadratic loop:

Case Graph build time Dulmage-Mendelsohn time N. variables
pglib_opf_case3_lmbd.m 2.55e-04 3.10e-04 19
pglib_opf_case5_pjm.m 3.60e-04 3.91e-04 35
pglib_opf_case14_ieee.m 6.64e-04 6.90e-04 109
pglib_opf_case24_ieee_rts.m 1.08e-03 1.02e-03 201
pglib_opf_case30_as.m 1.01e-03 1.08e-03 225
pglib_opf_case30_ieee.m 1.05e-03 1.02e-03 225
pglib_opf_case39_epri.m 1.14e-03 1.15e-03 263
pglib_opf_case57_ieee.m 1.76e-03 1.98e-03 435
pglib_opf_case60_c.m 1.99e-03 1.93e-03 473
pglib_opf_case73_ieee_rts.m 2.60e-03 2.22e-03 627
pglib_opf_case118_ieee.m 4.08e-03 1.99e-02 981
pglib_opf_case89_pegase.m 4.50e-03 4.87e-03 1019
pglib_opf_case200_activ.m 5.10e-03 6.18e-03 1381
pglib_opf_case179_goc.m 4.71e-03 5.91e-03 1411
pglib_opf_case162_ieee_dtc.m 5.22e-03 7.15e-03 1461
pglib_opf_case300_ieee.m 4.20e-02 1.04e-02 2245
pglib_opf_case240_pserc.m 3.68e-02 1.05e-02 2273
pglib_opf_case500_goc.m 1.44e-02 1.64e-02 3913
pglib_opf_case588_sdet.m 1.46e-02 1.72e-02 3921
pglib_opf_case793_goc.m 1.92e-02 2.38e-02 5239
pglib_opf_case1354_pegase.m 4.31e-02 7.09e-02 10673
pglib_opf_case1888_rte.m 7.86e-02 6.63e-02 13901
pglib_opf_case1951_rte.m 5.41e-02 7.50e-02 14287
pglib_opf_case2383wp_k.m 5.95e-02 8.32e-02 16351
pglib_opf_case2312_goc.m 9.13e-02 9.26e-02 16677
pglib_opf_case2000_goc.m 6.96e-02 1.71e-01 18533
pglib_opf_case2736sp_k.m 6.88e-02 1.08e-01 18549
pglib_opf_case2737sop_k.m 1.02e-01 1.16e-01 18551
pglib_opf_case2746wp_k.m 6.79e-02 8.92e-02 18609
pglib_opf_case2746wop_k.m 6.90e-02 1.42e-01 18721
pglib_opf_case3012wp_k.m 7.29e-02 1.12e-01 20313
pglib_opf_case2848_rte.m 1.09e-01 1.13e-01 20801
pglib_opf_case2868_rte.m 1.11e-01 1.17e-01 20969
pglib_opf_case3120sp_k.m 1.11e-01 1.28e-01 21013
pglib_opf_case2853_sdet.m 1.13e-01 9.85e-02 21391
pglib_opf_case3022_goc.m 8.30e-02 1.13e-01 22585
pglib_opf_case3375wp_k.m 1.20e-01 1.20e-01 23393
pglib_opf_case2869_pegase.m 1.26e-01 1.11e-01 24067
pglib_opf_case2742_goc.m 1.23e-01 1.97e-01 24177
pglib_opf_case4661_sdet.m 1.83e-01 1.89e-01 33311
pglib_opf_case3970_goc.m 1.30e-01 3.50e-01 34505
pglib_opf_case4020_goc.m 1.35e-01 3.99e-01 35993
pglib_opf_case4917_goc.m 1.95e-01 2.19e-01 36739
pglib_opf_case4601_goc.m 1.45e-01 3.80e-01 37999
pglib_opf_case4837_goc.m 2.10e-01 3.72e-01 40735
pglib_opf_case4619_goc.m 2.07e-01 4.06e-01 41839
pglib_opf_case6468_rte.m 2.40e-01 3.11e-01 48937
pglib_opf_case6470_rte.m 2.37e-01 3.73e-01 48961
pglib_opf_case6495_rte.m 2.34e-01 3.65e-01 49067
pglib_opf_case6515_rte.m 1.91e-01 3.82e-01 49179
pglib_opf_case10000_goc.m 7.81e-01 1.58e+00 72773
pglib_opf_case8387_pegase.m 3.49e-01 6.23e-01 75019
pglib_opf_case9241_pegase.m 4.45e-01 5.95e-01 82679
pglib_opf_case9591_goc.m 3.74e-01 1.09e+00 82843
pglib_opf_case10480_goc.m 4.09e-01 1.27e+00 95197
pglib_opf_case13659_pegase.m 5.19e-01 8.03e-01 109187
pglib_opf_case19402_goc.m 8.62e-01 2.83e+00 177621
pglib_opf_case24464_goc.m 9.98e-01 2.52e+00 200193
pglib_opf_case30000_goc.m 9.87e-01 2.83e+00 201573

This is more like what I would expect.

Robbybp commented 1 year ago

Now I need to:

Robbybp commented 1 year ago

Fixed in #15