betweenness_centrality_test.cpp

C++的一个好库。。。现在很流行

template<typename Graph, typename VertexIndexMap, typename CentralityMap>
void 
simple_unweighted_betweenness_centrality(const Graph& g, VertexIndexMap index,
                                         CentralityMap centrality)
{
  typedef typename boost::graph_traits<Graph>::vertex_descriptor vertex;
  typedef typename boost::graph_traits<Graph>::vertex_iterator vertex_iterator;
  typedef typename boost::graph_traits<Graph>::adjacency_iterator adjacency_iterator;
  typedef typename boost::graph_traits<Graph>::vertices_size_type vertices_size_type;
  typedef typename boost::property_traits<CentralityMap>::value_type centrality_type;

  vertex_iterator vi, vi_end;
  for (tie(vi, vi_end) = vertices(g); vi != vi_end; ++vi)
    put(centrality, *vi, 0);

  vertex_iterator si, si_end;
  for (tie(si, si_end) = vertices(g); si != si_end; ++si) {
    vertex s = *si;

    // S <-- empty stack
    std::stack<vertex> S;

    // P[w] <-- empty list, w \in V
    typedef std::vector<vertex> Predecessors;
    std::vector<Predecessors> predecessors(num_vertices(g));

    // sigma[t] <-- 0, t \in V
    std::vector<vertices_size_type> sigma(num_vertices(g), 0);

    // sigma[s] <-- 1
    sigma[get(index, s)] = 1;

    // d[t] <-- -1, t \in V
    std::vector<int> d(num_vertices(g), -1);

    // d[s] <-- 0
    d[get(index, s)] = 0;

    // Q <-- empty queue
    std::queue<vertex> Q;

    // enqueue s --> Q
    Q.push(s);

    while (!Q.empty()) {
      // dequeue v <-- Q
      vertex v = Q.front(); Q.pop();

      // push v --> S
      S.push(v);

      adjacency_iterator wi, wi_end;
      for (tie(wi, wi_end) = adjacent_vertices(v, g); wi != wi_end; ++wi) {
        vertex w = *wi;

        // w found for the first time?
        if (d[get(index, w)] < 0) {
          // enqueue w --> Q
          Q.push(w);
          
          // d[w] <-- d[v] + 1
          d[get(index, w)] = d[get(index, v)] + 1;
        }

        // shortest path to w via v?
        if (d[get(index, w)] == d[get(index, v)] + 1) {
          // sigma[w] = sigma[w] + sigma[v]
          sigma[get(index, w)] += sigma[get(index, v)];

          // append v --> P[w]
          predecessors[get(index, w)].push_back(v);
        }
      }
    }

    // delta[v] <-- 0, v \in V
    std::vector<centrality_type> delta(num_vertices(g), 0);

    // S returns vertices in order of non-increasing distance from s
    while (!S.empty()) {
      // pop w <-- S
      vertex w = S.top(); S.pop();

      const Predecessors& w_preds = predecessors[get(index, w)];
      for (typename Predecessors::const_iterator vi = w_preds.begin();
           vi != w_preds.end(); ++vi) {
        vertex v = *vi;
        // delta[v] <-- delta[v] + (sigma[v]/sigma[w])*(1 + delta[w])
        delta[get(index, v)] += 
          ((centrality_type)sigma[get(index, v)]/sigma[get(index, w)])
          * (1 + delta[get(index, w)]);
      }

      if (w != s) {
        // C_B[w] <-- C_B[w] + delta[w]
        centrality[w] += delta[get(index, w)];
      }
    }
  }

  typedef typename graph_traits<Graph>::directed_category directed_category;
  const bool is_undirected = 
    is_same<directed_category, undirected_tag>::value;
  if (is_undirected) {
    vertex_iterator v, v_end;
    for(tie(v, v_end) = vertices(g); v != v_end; ++v) {
      put(centrality, *v, get(centrality, *v) / centrality_type(2));
    }
  }
}

template<typename Graph>
void random_unweighted_test(Graph*, int n)
{
  Graph g(n);

  {
    typename graph_traits<Graph>::vertex_iterator v, v_end;
    int index = 0;
    for (tie(v, v_end) = boost::vertices(g); v != v_end; ++v, ++index) {
      put(vertex_index, g, *v, index);
    }
  }

  randomly_add_edges(g, 0.20);

  std::cout << "Random graph with " << n << " vertices and "
            << num_edges(g) << " edges.\n";

  std::cout << "  Direct translation of Brandes' algorithm...";
  std::vector<double> centrality(n);
  simple_unweighted_betweenness_centrality(g, get(vertex_index, g),
    make_iterator_property_map(centrality.begin(), get(vertex_index, g),
                               double()));
  std::cout << "DONE.\n";

  std::cout << "  Real version, unweighted...";
  std::vector<double> centrality2(n);
  brandes_betweenness_centrality(g, 
     make_iterator_property_map(centrality2.begin(), get(vertex_index, g),
                                double()));
  std::cout << "DONE.\n";

  if (!std::equal(centrality.begin(), centrality.end(),
                  centrality2.begin())) {
    for (std::size_t v = 0; v < centrality.size(); ++v) {
      double relative_error = 
        centrality[v] == 0.0? centrality2[v]
        : (centrality2[v] - centrality[v]) / centrality[v];
      if (relative_error < 0) relative_error = -relative_error;
      BOOST_CHECK(relative_error < error_tolerance);
    }
  }

  std::cout << "  Real version, weighted...";
  std::vector<double> centrality3(n);

  for (typename graph_traits<Graph>::edge_iterator ei = edges(g).first;
       ei != edges(g).second; ++ei)
    put(edge_weight, g, *ei, 1);

  brandes_betweenness_centrality(g, 
    weight_map(get(edge_weight, g))
    .centrality_map(
       make_iterator_property_map(centrality3.begin(), get(vertex_index, g),
                                  double())));
  std::cout << "DONE.\n";

  if (!std::equal(centrality.begin(), centrality.end(),
                  centrality3.begin())) {
    for (std::size_t v = 0; v < centrality.size(); ++v) {
      double relative_error = 
        centrality[v] == 0.0? centrality3[v]
        : (centrality3[v] - centrality[v]) / centrality[v];
      if (relative_error < 0) relative_error = -relative_error;
      BOOST_CHECK(relative_error < error_tolerance);
    }
  }
}

int test_main(int, char*[])
{
  typedef adjacency_list<listS, listS, undirectedS, 
                         property<vertex_index_t, int>, EdgeProperties> 
    Graph;
  typedef adjacency_list<listS, listS, directedS, 
                         property<vertex_index_t, int>, EdgeProperties> 
    Digraph;

  struct unweighted_edge ud_edge_init1[5] = { 
    { 0, 1 },
    { 0, 3 },
    { 1, 2 },
    { 3, 2 },
    { 2, 4 }
  };
  double ud_centrality1[5] = { 0.5, 1.0, 3.5, 1.0, 0.0 };
  run_unweighted_test((Graph*)0, 5, ud_edge_init1, 5, ud_centrality1);

  // Example borrowed from the JUNG test suite
  struct unweighted_edge ud_edge_init2[10] = { 
    { 0, 1 },
    { 0, 6 },
    { 1, 2 },
    { 1, 3 },
    { 2, 4 },
    { 3, 4 },
    { 4, 5 },
    { 5, 8 },
    { 7, 8 },
    { 6, 7 },
  };
  double ud_centrality2[9] = {
    0.2142 * 28, 
    0.2797 * 28,
    0.0892 * 28,
    0.0892 * 28,
    0.2797 * 28,
    0.2142 * 28,
    0.1666 * 28,
    0.1428 * 28,
    0.1666 * 28
  };
  double ud_edge_centrality2[10] = {
    10.66666,
    9.33333,
    6.5,
    6.5,
    6.5,
    6.5,
    10.66666,
    9.33333,
    8.0,
    8.0
  };

  run_unweighted_test((Graph*)0, 9, ud_edge_init2, 10, ud_centrality2,
                      ud_edge_centrality2);

  weighted_edge dw_edge_init1[6] = {
    { 0, 1, 1.0 },
    { 0, 3, 1.0 },
    { 1, 2, 0.5 },
    { 3, 1, 1.0 },
    { 3, 4, 1.0 },
    { 4, 2, 0.5 }
  };
  double dw_centrality1[5] = { 0.0, 1.5, 0.0, 1.0, 0.5 };
  run_weighted_test((Digraph*)0, 5, dw_edge_init1, 6, dw_centrality1);

  run_wheel_test((Graph*)0, 15);

  random_unweighted_test((Graph*)0, 300);

  return 0;
}