// Copyright 2004 The Trustees of Indiana University.

// Use, modification and distribution is subject to the Boost Software

// License, Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at

// http://www.boost.org/LICENSE_1_0.txt)

//  Authors: Douglas Gregor

//           Andrew Lumsdaine

/*

 * © Portions copyright (c) 2006-2007 Nokia Corporation.  All rights reserved.

*/

#include <boost/graph/betweenness_centrality.hpp>

#include <boost/graph/adjacency_list.hpp>

#include <vector>

#include <stack>

#include <queue>

#include <boost/property_map.hpp>

#include <boost/test/minimal.hpp>

#include <boost/random/uniform_01.hpp>

#include <boost/random/linear_congruential.hpp>

#ifdef __SYMBIAN32__

#include "std_log_result.h"

#define LOG_FILENAME_LINE __FILE__, __LINE__

#endif

using namespace boost;

const double error_tolerance = 0.001;

typedef property<edge_weight_t, double,

                 property<edge_index_t, std::size_t> > EdgeProperties;

struct weighted_edge 

{

  int source, target;

  double weight;

};

template<typename Graph>

void 

run_weighted_test(Graph*, int V, weighted_edge edge_init[], int E,

                  double correct_centrality[])

{

  Graph g(V);

  typedef typename graph_traits<Graph>::vertex_descriptor Vertex;

  typedef typename graph_traits<Graph>::vertex_iterator vertex_iterator;

  typedef typename graph_traits<Graph>::edge_descriptor Edge;

  std::vector<Vertex> vertices(V);

  {

    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);

      vertices[index] = *v;

    }

  }

  std::vector<Edge> edges(E);

  for (int e = 0; e < E; ++e) {

    edges[e] = add_edge(vertices[edge_init[e].source],

                        vertices[edge_init[e].target], 

                        g).first;

    put(edge_weight, g, edges[e], 1.0);

  }

  std::vector<double> centrality(V);

  brandes_betweenness_centrality(

    g,

    centrality_map(

      make_iterator_property_map(centrality.begin(), get(vertex_index, g),

                                 double()))

    .vertex_index_map(get(vertex_index, g)).weight_map(get(edge_weight, g)));

  for (int v = 0; v < V; ++v) {

    BOOST_CHECK(centrality[v] == correct_centrality[v]);

  }

}

struct unweighted_edge 

{

  int source, target;

};

template<typename Graph>

void 

run_unweighted_test(Graph*, int V, unweighted_edge edge_init[], int E,

                    double correct_centrality[],

                    double* correct_edge_centrality = 0)

{

  Graph g(V);

  typedef typename graph_traits<Graph>::vertex_descriptor Vertex;

  typedef typename graph_traits<Graph>::vertex_iterator vertex_iterator;

  typedef typename graph_traits<Graph>::edge_descriptor Edge;

  std::vector<Vertex> vertices(V);

  {

    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);

      vertices[index] = *v;

    }

  }

  std::vector<Edge> edges(E);

  for (int e = 0; e < E; ++e) {

    edges[e] = add_edge(vertices[edge_init[e].source],

                        vertices[edge_init[e].target], 

                        g).first;

    put(edge_weight, g, edges[e], 1.0);

    put(edge_index, g, edges[e], e);

  }

  std::vector<double> centrality(V);

  std::vector<double> edge_centrality1(E);

  brandes_betweenness_centrality(

    g,

    centrality_map(

      make_iterator_property_map(centrality.begin(), get(vertex_index, g),

                                 double()))

    .edge_centrality_map(

       make_iterator_property_map(edge_centrality1.begin(), 

                                  get(edge_index, g), double()))

    .vertex_index_map(get(vertex_index, g)));

  std::vector<double> centrality2(V);

  std::vector<double> edge_centrality2(E);

  brandes_betweenness_centrality(

    g,

    vertex_index_map(get(vertex_index, g)).weight_map(get(edge_weight, g))

    .centrality_map(

       make_iterator_property_map(centrality2.begin(), get(vertex_index, g),

                                  double()))

    .edge_centrality_map(

       make_iterator_property_map(edge_centrality2.begin(), 

                                  get(edge_index, g), double())));

  std::vector<double> edge_centrality3(E);

  brandes_betweenness_centrality(

    g,

    edge_centrality_map(

      make_iterator_property_map(edge_centrality3.begin(), 

                                 get(edge_index, g), double())));

  for (int v = 0; v < V; ++v) {

    BOOST_CHECK(centrality[v] == centrality2[v]);

    double relative_error = 

      correct_centrality[v] == 0.0? centrality[v]

      : (centrality[v] - correct_centrality[v]) / correct_centrality[v];

    if (relative_error < 0) relative_error = -relative_error;

    BOOST_CHECK(relative_error < error_tolerance);

  }  

  for (int e = 0; e < E; ++e) {

    BOOST_CHECK(edge_centrality1[e] == edge_centrality2[e]);

    BOOST_CHECK(edge_centrality1[e] == edge_centrality3[e]);

    if (correct_edge_centrality) {

      double relative_error = 

        correct_edge_centrality[e] == 0.0? edge_centrality1[e]

        : (edge_centrality1[e] - correct_edge_centrality[e]) 

        / correct_edge_centrality[e];

      if (relative_error < 0) relative_error = -relative_error;

      BOOST_CHECK(relative_error < error_tolerance);

      if (relative_error >= error_tolerance) {

        std::cerr << "Edge " << e << " has edge centrality " 

                  << edge_centrality1[e] << ", should be " 

                  << correct_edge_centrality[e] << std::endl;

      }

    }

  }

}

template<typename Graph>

void

run_wheel_test(Graph*, int V)

{

  typedef typename graph_traits<Graph>::vertex_descriptor Vertex;

  typedef typename graph_traits<Graph>::vertex_iterator vertex_iterator;

  typedef typename graph_traits<Graph>::edge_descriptor Edge;

  Graph g(V);

  Vertex center = *boost::vertices(g).first;

  std::vector<Vertex> vertices(V);

  {

    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);

      vertices[index] = *v;

      if (*v != center) {

        Edge e = add_edge(*v, center, g).first;

        put(edge_weight, g, e, 1.0);

      }

    }

  }

  std::vector<double> centrality(V);

  brandes_betweenness_centrality(

    g,

    make_iterator_property_map(centrality.begin(), get(vertex_index, g),

                               double()));

  std::vector<double> centrality2(V);

  brandes_betweenness_centrality(

    g,

    centrality_map(

      make_iterator_property_map(centrality2.begin(), get(vertex_index, g),

                                 double()))

    .vertex_index_map(get(vertex_index, g)).weight_map(get(edge_weight, g)));

  relative_betweenness_centrality(

    g,

    make_iterator_property_map(centrality.begin(), get(vertex_index, g),

                               double()));

  relative_betweenness_centrality(

    g,

    make_iterator_property_map(centrality2.begin(), get(vertex_index, g),

                               double()));

  for (int v = 0; v < V; ++v) {

    BOOST_CHECK(centrality[v] == centrality2[v]);

    BOOST_CHECK((v == 0 && centrality[v] == 1)

               || (v != 0 && centrality[v] == 0));

  } 

  double dominance = 

    central_point_dominance(

      g, 

      make_iterator_property_map(centrality2.begin(), get(vertex_index, g),

                                 double()));

  BOOST_CHECK(dominance == 1.0);

}

template<typename MutableGraph>

void randomly_add_edges(MutableGraph& g, double edge_probability)

{

  typedef typename graph_traits<MutableGraph>::directed_category

    directed_category;

  const bool is_undirected = 

    is_same<directed_category, undirected_tag>::value;

  minstd_rand gen;

  uniform_01<minstd_rand, double> rand_gen(gen);

  typedef typename graph_traits<MutableGraph>::vertex_descriptor vertex;

  typename graph_traits<MutableGraph>::vertex_iterator vi, vi_end;

  for (tie(vi, vi_end) = vertices(g); vi != vi_end; ++vi) {

    vertex v = *vi;

    typename graph_traits<MutableGraph>::vertex_iterator wi 

      = is_undirected? vi : vertices(g).first;

    while (wi != vi_end) {

      vertex w = *wi++;

      if (v != w) {

        if (rand_gen() < edge_probability) add_edge(v, w, g);

      }

    }

  }

}

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);

#ifdef __SYMBIAN32__

	std_log(LOG_FILENAME_LINE,"[End Test Case ]");

	testResultXml("betweenness_centrality_test");

	close_log_file();