test/max_matching_test.cc
author Peter Kovacs <kpeter@inf.elte.hu>
Fri, 03 Apr 2009 18:59:15 +0200
changeset 655 6ac5d9ae1d3d
parent 339 91d63b8b1a4c
child 637 b61354458b59
permissions -rw-r--r--
Support real types + numerical stability fix in NS (#254)

- Real types are supported by appropriate inicialization.
- A feature of the XTI spanning tree structure is removed to ensure
numerical stability (could cause problems using integer types).
The node potentials are updated always on the lower subtree,
in order to prevent overflow problems.
The former method isn't notably faster during to our tests.
deba@338
     1
/* -*- mode: C++; indent-tabs-mode: nil; -*-
deba@338
     2
 *
deba@338
     3
 * This file is a part of LEMON, a generic C++ optimization library.
deba@338
     4
 *
alpar@463
     5
 * Copyright (C) 2003-2009
deba@338
     6
 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
deba@338
     7
 * (Egervary Research Group on Combinatorial Optimization, EGRES).
deba@338
     8
 *
deba@338
     9
 * Permission to use, modify and distribute this software is granted
deba@338
    10
 * provided that this copyright notice appears in all copies. For
deba@338
    11
 * precise terms see the accompanying LICENSE file.
deba@338
    12
 *
deba@338
    13
 * This software is provided "AS IS" with no warranty of any kind,
deba@338
    14
 * express or implied, and with no claim as to its suitability for any
deba@338
    15
 * purpose.
deba@338
    16
 *
deba@338
    17
 */
deba@338
    18
deba@338
    19
#include <iostream>
deba@339
    20
#include <sstream>
deba@338
    21
#include <vector>
deba@338
    22
#include <queue>
deba@338
    23
#include <lemon/math.h>
deba@338
    24
#include <cstdlib>
deba@338
    25
deba@339
    26
#include <lemon/max_matching.h>
deba@339
    27
#include <lemon/smart_graph.h>
deba@339
    28
#include <lemon/lgf_reader.h>
deba@339
    29
deba@338
    30
#include "test_tools.h"
deba@338
    31
deba@338
    32
using namespace std;
deba@338
    33
using namespace lemon;
deba@338
    34
deba@339
    35
GRAPH_TYPEDEFS(SmartGraph);
deba@339
    36
deba@339
    37
deba@339
    38
const int lgfn = 3;
deba@339
    39
const std::string lgf[lgfn] = {
deba@339
    40
  "@nodes\n"
deba@339
    41
  "label\n"
deba@339
    42
  "0\n"
deba@339
    43
  "1\n"
deba@339
    44
  "2\n"
deba@339
    45
  "3\n"
deba@339
    46
  "4\n"
deba@339
    47
  "5\n"
deba@339
    48
  "6\n"
deba@339
    49
  "7\n"
deba@339
    50
  "@edges\n"
deba@339
    51
  "     label  weight\n"
deba@339
    52
  "7 4  0      984\n"
deba@339
    53
  "0 7  1      73\n"
deba@339
    54
  "7 1  2      204\n"
deba@339
    55
  "2 3  3      583\n"
deba@339
    56
  "2 7  4      565\n"
deba@339
    57
  "2 1  5      582\n"
deba@339
    58
  "0 4  6      551\n"
deba@339
    59
  "2 5  7      385\n"
deba@339
    60
  "1 5  8      561\n"
deba@339
    61
  "5 3  9      484\n"
deba@339
    62
  "7 5  10     904\n"
deba@339
    63
  "3 6  11     47\n"
deba@339
    64
  "7 6  12     888\n"
deba@339
    65
  "3 0  13     747\n"
deba@339
    66
  "6 1  14     310\n",
deba@339
    67
deba@339
    68
  "@nodes\n"
deba@339
    69
  "label\n"
deba@339
    70
  "0\n"
deba@339
    71
  "1\n"
deba@339
    72
  "2\n"
deba@339
    73
  "3\n"
deba@339
    74
  "4\n"
deba@339
    75
  "5\n"
deba@339
    76
  "6\n"
deba@339
    77
  "7\n"
deba@339
    78
  "@edges\n"
deba@339
    79
  "     label  weight\n"
deba@339
    80
  "2 5  0      710\n"
deba@339
    81
  "0 5  1      241\n"
deba@339
    82
  "2 4  2      856\n"
deba@339
    83
  "2 6  3      762\n"
deba@339
    84
  "4 1  4      747\n"
deba@339
    85
  "6 1  5      962\n"
deba@339
    86
  "4 7  6      723\n"
deba@339
    87
  "1 7  7      661\n"
deba@339
    88
  "2 3  8      376\n"
deba@339
    89
  "1 0  9      416\n"
deba@339
    90
  "6 7  10     391\n",
deba@339
    91
deba@339
    92
  "@nodes\n"
deba@339
    93
  "label\n"
deba@339
    94
  "0\n"
deba@339
    95
  "1\n"
deba@339
    96
  "2\n"
deba@339
    97
  "3\n"
deba@339
    98
  "4\n"
deba@339
    99
  "5\n"
deba@339
   100
  "6\n"
deba@339
   101
  "7\n"
deba@339
   102
  "@edges\n"
deba@339
   103
  "     label  weight\n"
deba@339
   104
  "6 2  0      553\n"
deba@339
   105
  "0 7  1      653\n"
deba@339
   106
  "6 3  2      22\n"
deba@339
   107
  "4 7  3      846\n"
deba@339
   108
  "7 2  4      981\n"
deba@339
   109
  "7 6  5      250\n"
deba@339
   110
  "5 2  6      539\n",
deba@339
   111
};
deba@339
   112
deba@339
   113
void checkMatching(const SmartGraph& graph,
deba@339
   114
                   const MaxMatching<SmartGraph>& mm) {
deba@339
   115
  int num = 0;
deba@339
   116
deba@339
   117
  IntNodeMap comp_index(graph);
deba@339
   118
  UnionFind<IntNodeMap> comp(comp_index);
deba@339
   119
deba@339
   120
  int barrier_num = 0;
deba@339
   121
deba@339
   122
  for (NodeIt n(graph); n != INVALID; ++n) {
deba@339
   123
    check(mm.decomposition(n) == MaxMatching<SmartGraph>::EVEN ||
deba@339
   124
          mm.matching(n) != INVALID, "Wrong Gallai-Edmonds decomposition");
deba@339
   125
    if (mm.decomposition(n) == MaxMatching<SmartGraph>::ODD) {
deba@339
   126
      ++barrier_num;
deba@339
   127
    } else {
deba@339
   128
      comp.insert(n);
deba@339
   129
    }
deba@339
   130
  }
deba@339
   131
deba@339
   132
  for (EdgeIt e(graph); e != INVALID; ++e) {
deba@339
   133
    if (mm.matching(e)) {
deba@339
   134
      check(e == mm.matching(graph.u(e)), "Wrong matching");
deba@339
   135
      check(e == mm.matching(graph.v(e)), "Wrong matching");
deba@339
   136
      ++num;
deba@339
   137
    }
deba@339
   138
    check(mm.decomposition(graph.u(e)) != MaxMatching<SmartGraph>::EVEN ||
deba@339
   139
          mm.decomposition(graph.v(e)) != MaxMatching<SmartGraph>::MATCHED,
deba@339
   140
          "Wrong Gallai-Edmonds decomposition");
deba@339
   141
deba@339
   142
    check(mm.decomposition(graph.v(e)) != MaxMatching<SmartGraph>::EVEN ||
deba@339
   143
          mm.decomposition(graph.u(e)) != MaxMatching<SmartGraph>::MATCHED,
deba@339
   144
          "Wrong Gallai-Edmonds decomposition");
deba@339
   145
deba@339
   146
    if (mm.decomposition(graph.u(e)) != MaxMatching<SmartGraph>::ODD &&
deba@339
   147
        mm.decomposition(graph.v(e)) != MaxMatching<SmartGraph>::ODD) {
deba@339
   148
      comp.join(graph.u(e), graph.v(e));
deba@339
   149
    }
deba@339
   150
  }
deba@339
   151
deba@339
   152
  std::set<int> comp_root;
deba@339
   153
  int odd_comp_num = 0;
deba@339
   154
  for (NodeIt n(graph); n != INVALID; ++n) {
deba@339
   155
    if (mm.decomposition(n) != MaxMatching<SmartGraph>::ODD) {
deba@339
   156
      int root = comp.find(n);
deba@339
   157
      if (comp_root.find(root) == comp_root.end()) {
deba@339
   158
        comp_root.insert(root);
deba@339
   159
        if (comp.size(n) % 2 == 1) {
deba@339
   160
          ++odd_comp_num;
deba@339
   161
        }
deba@339
   162
      }
deba@339
   163
    }
deba@339
   164
  }
deba@339
   165
deba@339
   166
  check(mm.matchingSize() == num, "Wrong matching");
deba@339
   167
  check(2 * num == countNodes(graph) - (odd_comp_num - barrier_num),
deba@339
   168
         "Wrong matching");
deba@339
   169
  return;
deba@339
   170
}
deba@339
   171
deba@339
   172
void checkWeightedMatching(const SmartGraph& graph,
deba@339
   173
                   const SmartGraph::EdgeMap<int>& weight,
deba@339
   174
                   const MaxWeightedMatching<SmartGraph>& mwm) {
deba@339
   175
  for (SmartGraph::EdgeIt e(graph); e != INVALID; ++e) {
deba@339
   176
    if (graph.u(e) == graph.v(e)) continue;
deba@339
   177
    int rw = mwm.nodeValue(graph.u(e)) + mwm.nodeValue(graph.v(e));
deba@339
   178
deba@339
   179
    for (int i = 0; i < mwm.blossomNum(); ++i) {
deba@339
   180
      bool s = false, t = false;
deba@339
   181
      for (MaxWeightedMatching<SmartGraph>::BlossomIt n(mwm, i);
deba@339
   182
           n != INVALID; ++n) {
deba@339
   183
        if (graph.u(e) == n) s = true;
deba@339
   184
        if (graph.v(e) == n) t = true;
deba@339
   185
      }
deba@339
   186
      if (s == true && t == true) {
deba@339
   187
        rw += mwm.blossomValue(i);
deba@339
   188
      }
deba@339
   189
    }
deba@339
   190
    rw -= weight[e] * mwm.dualScale;
deba@339
   191
deba@339
   192
    check(rw >= 0, "Negative reduced weight");
deba@339
   193
    check(rw == 0 || !mwm.matching(e),
deba@339
   194
          "Non-zero reduced weight on matching edge");
deba@339
   195
  }
deba@339
   196
deba@339
   197
  int pv = 0;
deba@339
   198
  for (SmartGraph::NodeIt n(graph); n != INVALID; ++n) {
deba@339
   199
    if (mwm.matching(n) != INVALID) {
deba@339
   200
      check(mwm.nodeValue(n) >= 0, "Invalid node value");
deba@339
   201
      pv += weight[mwm.matching(n)];
deba@339
   202
      SmartGraph::Node o = graph.target(mwm.matching(n));
deba@339
   203
      check(mwm.mate(n) == o, "Invalid matching");
deba@339
   204
      check(mwm.matching(n) == graph.oppositeArc(mwm.matching(o)),
deba@339
   205
            "Invalid matching");
deba@339
   206
    } else {
deba@339
   207
      check(mwm.mate(n) == INVALID, "Invalid matching");
deba@339
   208
      check(mwm.nodeValue(n) == 0, "Invalid matching");
deba@339
   209
    }
deba@339
   210
  }
deba@339
   211
deba@339
   212
  int dv = 0;
deba@339
   213
  for (SmartGraph::NodeIt n(graph); n != INVALID; ++n) {
deba@339
   214
    dv += mwm.nodeValue(n);
deba@339
   215
  }
deba@339
   216
deba@339
   217
  for (int i = 0; i < mwm.blossomNum(); ++i) {
deba@339
   218
    check(mwm.blossomValue(i) >= 0, "Invalid blossom value");
deba@339
   219
    check(mwm.blossomSize(i) % 2 == 1, "Even blossom size");
deba@339
   220
    dv += mwm.blossomValue(i) * ((mwm.blossomSize(i) - 1) / 2);
deba@339
   221
  }
deba@339
   222
deba@339
   223
  check(pv * mwm.dualScale == dv * 2, "Wrong duality");
deba@339
   224
deba@339
   225
  return;
deba@339
   226
}
deba@339
   227
deba@339
   228
void checkWeightedPerfectMatching(const SmartGraph& graph,
deba@339
   229
                          const SmartGraph::EdgeMap<int>& weight,
deba@339
   230
                          const MaxWeightedPerfectMatching<SmartGraph>& mwpm) {
deba@339
   231
  for (SmartGraph::EdgeIt e(graph); e != INVALID; ++e) {
deba@339
   232
    if (graph.u(e) == graph.v(e)) continue;
deba@339
   233
    int rw = mwpm.nodeValue(graph.u(e)) + mwpm.nodeValue(graph.v(e));
deba@339
   234
deba@339
   235
    for (int i = 0; i < mwpm.blossomNum(); ++i) {
deba@339
   236
      bool s = false, t = false;
deba@339
   237
      for (MaxWeightedPerfectMatching<SmartGraph>::BlossomIt n(mwpm, i);
deba@339
   238
           n != INVALID; ++n) {
deba@339
   239
        if (graph.u(e) == n) s = true;
deba@339
   240
        if (graph.v(e) == n) t = true;
deba@339
   241
      }
deba@339
   242
      if (s == true && t == true) {
deba@339
   243
        rw += mwpm.blossomValue(i);
deba@339
   244
      }
deba@339
   245
    }
deba@339
   246
    rw -= weight[e] * mwpm.dualScale;
deba@339
   247
deba@339
   248
    check(rw >= 0, "Negative reduced weight");
deba@339
   249
    check(rw == 0 || !mwpm.matching(e),
deba@339
   250
          "Non-zero reduced weight on matching edge");
deba@339
   251
  }
deba@339
   252
deba@339
   253
  int pv = 0;
deba@339
   254
  for (SmartGraph::NodeIt n(graph); n != INVALID; ++n) {
deba@339
   255
    check(mwpm.matching(n) != INVALID, "Non perfect");
deba@339
   256
    pv += weight[mwpm.matching(n)];
deba@339
   257
    SmartGraph::Node o = graph.target(mwpm.matching(n));
deba@339
   258
    check(mwpm.mate(n) == o, "Invalid matching");
deba@339
   259
    check(mwpm.matching(n) == graph.oppositeArc(mwpm.matching(o)),
deba@339
   260
          "Invalid matching");
deba@339
   261
  }
deba@339
   262
deba@339
   263
  int dv = 0;
deba@339
   264
  for (SmartGraph::NodeIt n(graph); n != INVALID; ++n) {
deba@339
   265
    dv += mwpm.nodeValue(n);
deba@339
   266
  }
deba@339
   267
deba@339
   268
  for (int i = 0; i < mwpm.blossomNum(); ++i) {
deba@339
   269
    check(mwpm.blossomValue(i) >= 0, "Invalid blossom value");
deba@339
   270
    check(mwpm.blossomSize(i) % 2 == 1, "Even blossom size");
deba@339
   271
    dv += mwpm.blossomValue(i) * ((mwpm.blossomSize(i) - 1) / 2);
deba@339
   272
  }
deba@339
   273
deba@339
   274
  check(pv * mwpm.dualScale == dv * 2, "Wrong duality");
deba@339
   275
deba@339
   276
  return;
deba@339
   277
}
deba@339
   278
deba@339
   279
deba@338
   280
int main() {
deba@338
   281
deba@339
   282
  for (int i = 0; i < lgfn; ++i) {
deba@339
   283
    SmartGraph graph;
deba@339
   284
    SmartGraph::EdgeMap<int> weight(graph);
deba@338
   285
deba@339
   286
    istringstream lgfs(lgf[i]);
deba@339
   287
    graphReader(graph, lgfs).
deba@339
   288
      edgeMap("weight", weight).run();
deba@338
   289
deba@339
   290
    MaxMatching<SmartGraph> mm(graph);
deba@339
   291
    mm.run();
deba@339
   292
    checkMatching(graph, mm);
deba@338
   293
deba@339
   294
    MaxWeightedMatching<SmartGraph> mwm(graph, weight);
deba@339
   295
    mwm.run();
deba@339
   296
    checkWeightedMatching(graph, weight, mwm);
deba@338
   297
deba@339
   298
    MaxWeightedPerfectMatching<SmartGraph> mwpm(graph, weight);
deba@339
   299
    bool perfect = mwpm.run();
deba@338
   300
deba@339
   301
    check(perfect == (mm.matchingSize() * 2 == countNodes(graph)),
deba@339
   302
          "Perfect matching found");
deba@338
   303
deba@339
   304
    if (perfect) {
deba@339
   305
      checkWeightedPerfectMatching(graph, weight, mwpm);
deba@338
   306
    }
deba@338
   307
  }
deba@338
   308
deba@338
   309
  return 0;
deba@338
   310
}