diff --git a/include/boost/graph/personalized_page_rank.hpp b/include/boost/graph/personalized_page_rank.hpp new file mode 100644 index 000000000..f312e64ba --- /dev/null +++ b/include/boost/graph/personalized_page_rank.hpp @@ -0,0 +1,248 @@ +// Copyright 2026 Emmanouil Krasanakis + +// Distributed under 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: Emmanouil Krasanakis + +#ifndef BOOST_GRAPH_PERSONALIZED_PAGE_RANK_HPP +#define BOOST_GRAPH_PERSONALIZED_PAGE_RANK_HPP + +#include +#include +#include +#include +#include +#include +#include +#include + +namespace boost +{ +namespace graph +{ + struct rank_convergence + { + explicit rank_convergence(std::size_t iters, double tol=0) : iters(iters), tol(tol) {} // allowing tolerance for early stopping + template < typename RankMap, typename RankMap2, typename Graph > + bool operator()(const RankMap& current, const RankMap2& previous, const Graph& g) + { + if (--iters == 0) + return true; + if (!tol) + return false; + using rank_type = typename property_traits< RankMap >::value_type; + rank_type sum_abs(0); + for (auto v : boost::make_iterator_range(vertices(g))) + sum_abs += std::abs(get(current, v) - get(previous, v)); + return sum_abs*num_vertices(g) + void personalized_page_rank_step( + const Graph& g, + WeightMap weight_map, + PersonalizationMap personalization_map, + RankMap from_rank, + RankMap2 to_rank, + typename property_traits< RankMap >::value_type damping, + incidence_graph_tag) + { + using rank_type = typename property_traits< RankMap >::value_type; + rank_type l1_norm(0); // Computing the norm simultaneously avoids an extra summing iteration. + + // Initialize the constant part of maps. + for (auto v : boost::make_iterator_range(vertices(g))) + { + auto v_constant = rank_type(1 - damping) * get(personalization_map, v); + put(to_rank, v, v_constant); + l1_norm += v_constant; + } + + // Accumulate from neighbors. + for (auto u : boost::make_iterator_range(vertices(g))) + { + rank_type u_rank_factor = damping * get(from_rank, u); + rank_type l1_accumulated_norm(0); // TBD: Consider making l1_norm volatile to reduce accumulation errors. + for (auto e : boost::make_iterator_range(out_edges(u, g))) + { + auto v = target(e, g); + rank_type u_rank_out = get(weight_map, e)*u_rank_factor; + put(to_rank, v, get(to_rank, v) + u_rank_out); + l1_accumulated_norm += u_rank_out; + } + l1_norm += l1_accumulated_norm; + } + // If there are negative edge weights, or if negative damping is used, l1_norm could be zero or near-zero. + // Division in those cases is conceptually correct for floating point weights, and actually expected behavior. + // That said, such edge cases are impossible to arise for all typical algorithm uses. + for (auto v : boost::make_iterator_range(vertices(g))) + put(to_rank, v, get(to_rank, v)/l1_norm); + } + + template < + typename Graph, + typename WeightMap, + typename PersonalizationMap, + typename RankMap, + typename RankMap2 > + void personalized_page_rank_step( + const Graph& g, + WeightMap weight_map, + PersonalizationMap personalization_map, + RankMap from_rank, + RankMap2 to_rank, + typename property_traits< RankMap >::value_type damping, + bidirectional_graph_tag) + { + using damping_type = typename property_traits< RankMap >::value_type; + damping_type l1_norm(0); // Computing the norm simultaneously avoids an extra summing iteration. + for (auto v : boost::make_iterator_range(vertices(g))) + { + damping_type rank(0); + for (auto e : boost::make_iterator_range(in_edges(v, g))) + rank += get(from_rank, source(e, g))*get(weight_map, e); + auto v_score = (damping_type(1) - damping) * get(personalization_map, v) + damping * rank; + put(to_rank, v, v_score); + l1_norm += v_score; + } + // See above function for potential division by zero comments. + for (auto v : boost::make_iterator_range(vertices(g))) + put(to_rank, v, get(to_rank, v)/l1_norm); + } + } // end namespace personalized_page_rank_detail + + template < + typename Graph, + typename WeightMap, + typename PersonalizationMap, + typename RankMap, + typename Done, + typename RankMap2 > + Done personalized_page_rank( + const Graph& g, + WeightMap weight_map, + PersonalizationMap personalization_map, + RankMap rank_map, + Done done, + typename property_traits< RankMap >::value_type damping, + RankMap2 rank_map2 + BOOST_GRAPH_ENABLE_IF_MODELS_PARM(Graph, vertex_list_graph_tag)) + { + using Vertex = typename graph_traits::vertex_descriptor; + using Edge = typename graph_traits::edge_descriptor; + BOOST_CONCEPT_ASSERT(( boost::IncidenceGraphConcept )); + BOOST_CONCEPT_ASSERT(( boost::VertexListGraphConcept )); + BOOST_CONCEPT_ASSERT(( boost::ReadablePropertyMapConcept )); + BOOST_CONCEPT_ASSERT(( boost::ReadablePropertyMapConcept)); + BOOST_CONCEPT_ASSERT(( boost::ReadWritePropertyMapConcept )); + + + assert (damping>=-1.0 && damping<1.0 && "Damping outside the closed-open range [-1.0,1.0) could induce numerical instability."); // non-inclussive upper limit is deliberate + + using rank_type = typename property_traits< PersonalizationMap >::value_type; + rank_type personalization_norm(0); + for (auto v : boost::make_iterator_range(vertices(g))) + personalization_norm += get(personalization_map, v); + + // TBD: This implementation couples iterators when possible under reduced L1 cache invalidation assumptions, + // but this is not necessarily the case because we may be grabbing 2x memory lanes each time to write there. + // Could investigate which pattern is faster in the future. + for (auto v : boost::make_iterator_range(vertices(g))) + { + rank_type value = get(personalization_map, v)/personalization_norm; + put(personalization_map, v, value); + put(rank_map, v, value); + } + + bool to_map_2 = true; + do + { + typedef typename graph_traits< Graph >::traversal_category category; + if (to_map_2) + personalized_page_rank_detail::personalized_page_rank_step(g, weight_map, personalization_map, rank_map, rank_map2, damping, category()); + else + personalized_page_rank_detail::personalized_page_rank_step(g, weight_map, personalization_map, rank_map2, rank_map, damping, category()); + to_map_2 = !to_map_2; + } + while ((to_map_2 && !done(rank_map, rank_map2, g)) || (!to_map_2 && !done(rank_map2, rank_map, g))); // Done may not be symmetric. + + // Now multiply the result with personalization_norm to restore the order of magnitude and store it in rank_map. + // Also restore the original personalization_map's magnitude for reuse (this is lossy up to numerical tolerance + // but leaner than making a copy). + if (!to_map_2) + { + for (auto v : boost::make_iterator_range(vertices(g))) + { + put(rank_map, v, get(rank_map2, v)*personalization_norm); + put(personalization_map, v, get(personalization_map, v)*personalization_norm); + } + } + else + { + for (auto v : boost::make_iterator_range(vertices(g))) + { + put(rank_map, v, get(rank_map, v)*personalization_norm); + put(personalization_map, v, get(personalization_map, v)*personalization_norm); + } + } + return done; + } + + template < + typename Graph, + typename WeightMap, + typename PersonalizationMap, + typename RankMap, + typename Done > + Done personalized_page_rank( + const Graph& g, + WeightMap weight_map, + PersonalizationMap personalization_map, + RankMap rank_map, + Done done, + typename property_traits< RankMap >::value_type damping) + { + using rank_type = typename property_traits< RankMap >::value_type; + std::vector< rank_type > ranks2(num_vertices(g)); + return personalized_page_rank(g, weight_map, personalization_map, rank_map, done, damping, + make_iterator_property_map(ranks2.begin(), get(vertex_index, g))); + } + + template < typename Graph, typename PersonalizationMap, typename RankMap > + rank_convergence personalized_page_rank( + const Graph& g, + PersonalizationMap personalization_map, + RankMap rank_map, + typename property_traits< RankMap >::value_type damping=0.85) + { + // This is the most traditional personalized PageRank implementation, with minimized signature. + using Edge = typename graph_traits::edge_descriptor; + using rank_type = typename property_traits< RankMap >::value_type; + std::vector< rank_type > ranks2(num_vertices(g)); + auto markovian_weights = make_function_property_map([&g](Edge e){ return 1.0 / out_degree(source(e, g), g); }); + return personalized_page_rank(g, + markovian_weights, + personalization_map, + rank_map, + rank_convergence(100, 1.E-9), + damping, + make_iterator_property_map(ranks2.begin(), get(vertex_index, g))); + } + +} +} // end namespace boost::graph + +#endif // BOOST_GRAPH_PERSONALIZED_PAGE_RANK_HPP diff --git a/test/Jamfile.v2 b/test/Jamfile.v2 index 0ef7513fe..faaa3e31b 100644 --- a/test/Jamfile.v2 +++ b/test/Jamfile.v2 @@ -175,6 +175,7 @@ alias graph_test_regular : [ run delete_edge.cpp ] [ run johnson-test.cpp ] [ run lvalue_pmap.cpp ] + [ run personalized_pagerank_test.cpp ] ; alias graph_test_with_filesystem : : diff --git a/test/personalized_pagerank_test.cpp b/test/personalized_pagerank_test.cpp new file mode 100644 index 000000000..1cade31f1 --- /dev/null +++ b/test/personalized_pagerank_test.cpp @@ -0,0 +1,77 @@ + +#include +#include +#include +#include +#include +#include + +using DirectedGraph = typename boost::adjacency_list; +using DirectedVertex = typename boost::graph_traits::vertex_descriptor; +using DirectedEdge = typename boost::graph_traits::edge_descriptor; + +struct custom_rank_convergence: public boost::graph::rank_convergence +{ + explicit custom_rank_convergence(std::size_t iters, double tol=0) : boost::graph::rank_convergence(iters,tol) {} + std::size_t get_remaining_iters() const { return iters; } +}; + +void directed_graph_tests(double damping, double renormalize) +{ + std::vector>, int>> graph_defs = { + // deliberately hard symmetric graph + {{ + {0,1},{1,0},{1,2},{2,1},{2,3},{3,2}, + {4,5},{5,4},{5,6},{6,5},{6,7},{7,6},{7,8},{8,7},{8,9},{9,8},{9,10},{10,9}, + {0,3},{3,0},{1,3},{3,1},{1,4},{4,1}, + {4,6},{6,4},{6,9},{9,6},{6,8},{8,6},{7,9},{9,7},{8,10},{10,8}, + {11,10},{10,11},{10,12},{12,10} + }, 13}, + // undirected circle with 0 and 2 being symmetric, and a non-symmetric directed one-way blocks 2 hops away from 0 ({7,6} blocked) and 2 ({4,5} blocked) + {{ {0,1},{1,0},{1,3},{3,1},{3,2},{2,3},{2,4},{4,2},{5,4},{5,6},{6,5},{6,7},{7,0},{0,7}}, 8}, + // fully undirected graph + {{ {0,1}, {2,1}, {0,3}, {2,3}, {3,4}, {4,5}, {1,5}, {5,2}, {5,0}}, 6}, + // same as above but missing incoming edges for 0 and 2 (whole graph is a sink, needs correct normalization that guards against zero to not yield nans) + {{ {0,1}, {2,1}, {0,3}, {2,3}, {3,4}, {4,5}, {1,5}}, 6} + }; + for(auto& graph_details : graph_defs) + { + DirectedGraph g(graph_details.first.begin(), graph_details.first.end(), graph_details.second); + std::vector ranks(num_vertices(g)); + auto rank_map = boost::make_iterator_property_map(ranks.begin(), get(boost::vertex_index, g)); + std::vector personalization(num_vertices(g)); + auto personalization_map = boost::make_iterator_property_map(personalization.begin(), get(boost::vertex_index, g)); + personalization[0] = 1; + personalization[1] = 1; + personalization[2] = 1; + personalization[3] = 1; + + std::size_t max_iters(300); // Convergence is just that bad in tested graphs; usually it's much lower.' + auto weight = boost::make_function_property_map([&g,renormalize](DirectedEdge e){ + auto denom = out_degree(source(e, g), g) * out_degree(target(e, g), g); + if(denom==0) return 0.0; + return 1.0 / std::sqrt(renormalize+double(denom)); + }); + auto convergence = custom_rank_convergence(max_iters, 1.E-9); + convergence = boost::graph::personalized_page_rank(g, weight, personalization_map, rank_map, convergence, damping); + + // the following asserts hold for all tested graphs and personalization: 0,1,2,3 plus some other nodes is a mini-cluster with 0 and 2 being structurally symmetric + assert(convergence.get_remaining_iters()0 || (damping<0.0)); // converged (derivatives may not converge) + assert((ranks[0]1.0); // holds for all graphs given low-pass damping + assert((ranks[0]!=ranks[1]) == (damping!=0.0)); // but not the same normally, the same 1.0 personalization if damping is zero + assert(ranks[0]==ranks[2] || (damping<=-1.0)); // equal due to symmetry, even under non-convergence and floating coarseness + assert(std::abs(ranks[0]-ranks[2])<1.E-14); // approximately equal in cases where even addition order matters + assert((ranks[0]<0.0)||(ranks[1]<0.0)||(ranks[num_vertices(g)-1]<=0.0)||(damping>0.0)); // ensure that negative flows are possible for negative damping + } +} + +int main(int, char*[]) +{ + directed_graph_tests(0.9, 0); // normal mode + directed_graph_tests(0.9, 1.0); // with renormalization + directed_graph_tests(0.99, 0); // huge damping (asymptotically exponentially slow convergence as we approach 1.0) + directed_graph_tests(0.0, 0); // just yield the personalization again + directed_graph_tests(-0.8, 1.0); // negative flow = some kind of derivate + directed_graph_tests(-1.0, 1.0); // huge negative flow +}