diff --git a/experimental/algorithm/LAGr_EdgeBetweennessCentrality.c b/experimental/algorithm/LAGr_EdgeBetweennessCentrality.c new file mode 100644 index 0000000000..9fea4ae129 --- /dev/null +++ b/experimental/algorithm/LAGr_EdgeBetweennessCentrality.c @@ -0,0 +1,330 @@ +//------------------------------------------------------------------------------ +// LAGr_EdgeBetweennessCentrality: edge betweenness-centrality +//------------------------------------------------------------------------------ + +// LAGraph, (c) 2019-2022 by The LAGraph Contributors, All Rights Reserved. +// SPDX-Licene-Identifier: BSD-2-Clause +// +// For additional details (including references to third party source code and +// other files) see the LICEnE file or contact permission@sei.cmu.edu. See +// Contributors.txt for a full list of contributors. Created, in part, with +// funding and support from the U.S. Government (see Acknowledgments.txt file). +// DM22-0790 + +// Contributed by Casey Pei and Tim Davis, Texas A&M University; +// Adapted and revised from GraphBLAS C API Spec, Appendix B.4. + +//------------------------------------------------------------------------------ + +// LAGr_EdgeBetweennessCentrality: Exact algorithm for computing +// betweeness centrality. + +// This is an Advanced algorithm (G->AT is required). + + +//------------------------------------------------------------------------------ + +#define LG_FREE_WORK \ +{ \ + GrB_free (&frontier) ; \ + GrB_free (&J_vec) ; \ + GrB_free (&I_vec) ; \ + GrB_free (&J_matrix) ; \ + GrB_free (&I_matrix) ; \ + GrB_free (&Fd1A) ; \ + GrB_free (&paths) ; \ + GrB_free (&bc_vertex_flow) ; \ + GrB_free (&temp_update) ; \ + GrB_free (&Add_One_Divide) ; \ + GrB_free (&Update) ; \ + if (Search != NULL) \ + { \ + for (int64_t i = 0 ; i < n ; i++) \ + { \ + GrB_free (&(Search [i])) ; \ + } \ + LAGraph_Free ((void **) &Search, NULL) ; \ + } \ +} + +#define LG_FREE_ALL \ +{ \ + LG_FREE_WORK ; \ + GrB_free (centrality) ; \ +} + +#include "LG_internal.h" +#include + +#undef LAGRAPH_CATCH +#define LAGRAPH_CATCH(status) \ +{ \ + print ("LAGraph failure (file %s, line %d): status: %d", \ + __FILE__, __LINE__, status) ; \ + LG_ERROR_MSG ("LAGraph failure (file %s, line %d): status: %d", \ + __FILE__, __LINE__, status) ; \ + LG_FREE_ALL ; \ + return (status) ; \ +} + +#undef GRB_CATCH +#define GRB_CATCH(info) \ +{ \ + printf ("GraphBLAS failure (file %s, line %d): info: %d", \ + __FILE__, __LINE__, info) ; \ + LG_ERROR_MSG ("GraphBLAS failure (file %s, line %d): info: %d", \ + __FILE__, __LINE__, info) ; \ + LG_FREE_ALL ; \ + return (info) ; \ +} + +//------------------------------------------------------------------------------ +// (1+x)/y function for double: z = (1 + x) / y +//------------------------------------------------------------------------------ + +void add_one_divide_function (void *z, const void *x, const void *y) +{ + double a = (*((double *) x)) ; + double b = (*((double *) y)) ; + (*((double *) z)) = (1 + a) / b ; +} + +#define ADD_ONE_DIVIDE_FUNCTION_DEFN \ +"void add_one_divide_function (void *z, const void *x, const void *y) \n" \ +"{ \n" \ +" double a = (*((double *) x)) ; \n" \ +" double b = (*((double *) y)) ; \n" \ +" (*((double *) z)) = (1 + a) / b ; \n" \ +"}" + +//------------------------------------------------------------------------------ +// LAGr_EdgeBetweennessCentrality: edge betweenness-centrality +//------------------------------------------------------------------------------ + +int LAGr_EdgeBetweennessCentrality +( + // output: + GrB_Matrix *centrality, // centrality(i): betweeness centrality of i + // input: + LAGraph_Graph G, // input graph + char *msg +) +{ + + //-------------------------------------------------------------------------- + // check inputs + //-------------------------------------------------------------------------- + + LG_CLEAR_MSG ; + + // Array of BFS search matrices. + // Search [i] is a sparse matrix that stores the depth at which each vertex is + // first seen thus far in each BFS at the current depth i. Each column + // corresponds to a BFS traversal starting from a source node. + GrB_Vector *Search = NULL ; + + // Frontier vector, a sparse matrix. + // Stores # of shortest paths to vertices at current BFS depth + GrB_Vector frontier = NULL ; + + // Paths matrix holds the number of shortest paths for each node and + // starting node discovered so far. A dense matrix that is updated with + // sparse updates, and also used as a mask. + GrB_Vector paths = NULL ; + + // The betweenness centrality each vertex. A dense vector. + GrB_Vector bc_vertex_flow = NULL ; + + // Update matrix for betweenness centrality for each edge. A sparse matrix. + GrB_Matrix Update = NULL ; + + GrB_BinaryOp Add_One_Divide = NULL ; + + GrB_Vector J_vec = NULL ; + GrB_Vector I_vec = NULL ; + GrB_Matrix I_matrix = NULL ; + GrB_Matrix J_matrix = NULL ; + GrB_Matrix Fd1A = NULL ; + GrB_Vector temp_update = NULL ; + + GrB_Index n = 0 ; // # nodes in the graph + + LG_ASSERT (centrality != NULL, GrB_NULL_POINTER) ; + (*centrality) = NULL ; + LG_TRY (LAGraph_CheckGraph (G, msg)) ; + + GrB_Matrix A = G->A ; + GrB_Matrix AT ; + if (G->kind == LAGraph_ADJACENCY_UNDIRECTED || + G->is_symmetric_structure == LAGraph_TRUE) + { + // A and A' have the same structure + AT = A ; + } + else + { + // A and A' differ + AT = G->AT ; + LG_ASSERT_MSG (AT != NULL, LAGRAPH_NOT_CACHED, "G->AT is required") ; + } + + // ========================================================================= + // === initialization ===================================================== + // ========================================================================= + + GRB_TRY (GxB_BinaryOp_new (&Add_One_Divide, add_one_divide_function, + GrB_FP64, GrB_FP64, GrB_FP64, + "add_one_divide_function", ADD_ONE_DIVIDE_FUNCTION_DEFN)) ; + + // Initialize the frontier, paths, Update, and bc_vertex_flow + GRB_TRY (GrB_Matrix_nrows (&n, A)) ; + GRB_TRY (GrB_Vector_new (&paths, GrB_FP64, n)) ; + GRB_TRY (GrB_Vector_new (&frontier, GrB_FP64, n)) ; + GRB_TRY (GrB_Matrix_new (&Update, GrB_FP64, n, n)) ; + GRB_TRY (GrB_Vector_new (&bc_vertex_flow, GrB_FP64, n)) ; + + + // Initialize centrality matrix with zeros using A as structural mask + LG_TRY (GrB_Matrix_new(centrality, GrB_FP64, n, n)) ; + GRB_TRY (GrB_assign (*centrality, A, NULL, 0.0, GrB_ALL, n, GrB_ALL, n, GrB_DESC_S)) ; + + // Allocate memory for the array of S vectors + LG_TRY (LAGraph_Calloc ((void **) &Search, n+1, sizeof (GrB_Vector), msg)) ; + + // ========================================================================= + // === Breadth-first search stage ========================================== + // ========================================================================= + + GrB_Index frontier_size, last_frontier_size = 0 ; + GRB_TRY (GrB_Vector_nvals (&frontier_size, frontier)) ; + + int64_t depth, root ; + for (root = 0 ; root < n ; root++) + { + depth = 0 ; + + // root frontier: Search [0](root) = true + GrB_free (&(Search [0])) ; + GRB_TRY (GrB_Vector_new(&(Search [0]), GrB_BOOL, n)) ; + GRB_TRY (GrB_Vector_setElement_BOOL(Search [0], (bool) true, root)) ; + + // clear paths, and then set paths (root) = 1 + GRB_TRY (GrB_Vector_clear (paths)) ; + GRB_TRY (GrB_Vector_setElement (paths, (double) 1.0, root)) ; + + GRB_TRY (GrB_Matrix_clear (Update)) ; + + // Extract row root from A into frontier vector: frontier = AT(root,:) + GRB_TRY (GrB_Col_extract (frontier, NULL, NULL, AT, GrB_ALL, n, root, + NULL)) ; + GRB_TRY (GrB_Vector_nvals (&frontier_size, frontier)) ; + + while (frontier_size != 0) + { + depth++ ; + + //---------------------------------------------------------------------- + // Accumulate path counts: paths += frontier + //---------------------------------------------------------------------- + + GRB_TRY (GrB_assign (paths, NULL, GrB_PLUS_FP64, frontier, GrB_ALL, n, + NULL)) ; + + //---------------------------------------------------------------------- + // Add frontier to S: S(depth, :) = frontier + //---------------------------------------------------------------------- + + GrB_free (&(Search [depth])) ; + LG_TRY (LAGraph_Vector_Structure (&(Search [depth]), frontier, msg)) ; + + //---------------------------------------------------------------------- + // Update frontier: frontier = frontier*A + //---------------------------------------------------------------------- + + GRB_TRY (LG_SET_FORMAT_HINT (frontier, LG_SPARSE)) ; + GRB_TRY (GrB_vxm (frontier, paths, NULL, GxB_PLUS_FIRST_FP64, frontier, + A, GrB_DESC_RSC )) ; + + //---------------------------------------------------------------------- + // Get size of current frontier: frontier_size = nvals(frontier) + //---------------------------------------------------------------------- + + last_frontier_size = frontier_size ; + GRB_TRY (GrB_Vector_nvals (&frontier_size, frontier)) ; + } + + + // ========================================================================= + // === Betweenness centrality computation phase ============================ + // ========================================================================= + + // bc_vertex_flow = ones (n, n) ; a full matrix (and stays full) + GRB_TRY (GrB_Vector_new (&bc_vertex_flow, GrB_FP64, n)) ; + GRB_TRY (GrB_assign(bc_vertex_flow, NULL, NULL, 0.0, GrB_ALL, n, NULL)) ; + + GRB_TRY (GrB_Vector_new(&J_vec, GrB_FP64, n)) ; + GRB_TRY (GrB_Vector_new (&I_vec, GrB_FP64, n)) ; + GRB_TRY (GrB_Matrix_new (&Fd1A, GrB_FP64, n, n)) ; + GRB_TRY (GrB_Vector_new(&temp_update, GrB_FP64, n)) ; // Create a temporary vector + + // Backtrack through the BFS and compute centrality updates for each vertex + // GrB_Index fd1_size; + while (depth >= 1) + { + GrB_Vector f_d = Search [depth] ; + GrB_Vector f_d1 = Search [depth - 1] ; + + // 18 w = S(d, :) ÷ p × v + S(d, :) + // 19 Update = A .× w + // 20 w = S(d − 1, :) × p + // 21 Update = w .× U + + // make J Matrix + GRB_TRY (GrB_eWiseMult(J_vec, f_d, NULL, Add_One_Divide, bc_vertex_flow, paths, GrB_DESC_RS)) ; + // GRB_TRY (GrB_eWiseMult(J_vec, f_d, NULL, GrB_PLUS_FP64, bc_vertex_flow, paths, GrB_DESC_RS)) ; + GRB_TRY (GrB_Matrix_diag(&J_matrix, J_vec, 0)) ; + + // make I matrix + GRB_TRY (GrB_Vector_extract (I_vec, f_d1, NULL, paths, GrB_ALL, n, GrB_DESC_RS)) ; + GRB_TRY (GrB_Matrix_diag(&I_matrix, I_vec, 0)) ; + + // combine + + // intermediate matrix for Fd1 * A + // GRB_TRY (GrB_eWiseMult(Fd1A, NULL, NULL, GrB_TIMES_FP64, I_matrix, AT, NULL)) ; + GRB_TRY(GrB_mxm(Fd1A, NULL, NULL, GrB_PLUS_TIMES_SEMIRING_FP64, + I_matrix, A, NULL)) ; + + // GRB_TRY (GrB_eWiseMult(U, NULL, NULL, GrB_TIMES_FP64, Fd1A, J_matrix, NULL)) ; + GRB_TRY(GrB_mxm(Update, NULL, NULL, GrB_PLUS_TIMES_SEMIRING_FP64, + Fd1A, J_matrix, NULL)) ; + + + // 22 centrality = centrality + Update + // GRB_TRY (GrB_assign(centrality, centrality, GrB_PLUS_FP64, Update, GrB_ALL, n, GrB_ALL, n, + // GrB_DESC_S)) ; + GRB_TRY (GrB_eWiseAdd (*centrality, NULL, NULL, GrB_PLUS_FP64, *centrality, Update, NULL)) ; + + + // 23 v = Update +. + + // Reduce "update" matrix to a vector (sum each column) + GRB_TRY (GrB_reduce(temp_update, NULL, NULL, GrB_PLUS_MONOID_FP64, Update, NULL)) ; + // GRB_TRY (GrB_reduce(temp_update, NULL, NULL, GxB_ANY_FP64_MONOID, Update, NULL)) ; + GRB_TRY (GrB_eWiseAdd(bc_vertex_flow, NULL, NULL, GrB_PLUS_FP64, bc_vertex_flow, temp_update, NULL)) ; + + // 24 d = d − 1 + depth-- ; + } + + + } + + + // ========================================================================= + // === finalize the centrality ============================================= + // ========================================================================= + + LG_FREE_WORK ; + return (GrB_SUCCESS) ; +} diff --git a/experimental/test/LG_check_edgeBetweennessCentrality.c b/experimental/test/LG_check_edgeBetweennessCentrality.c index aa89061d4a..602df6c09a 100644 --- a/experimental/test/LG_check_edgeBetweennessCentrality.c +++ b/experimental/test/LG_check_edgeBetweennessCentrality.c @@ -19,10 +19,10 @@ #define LG_FREE_WORK \ { \ LAGraph_Free ((void **) &queue, NULL) ; \ - LAGraph_Free ((void **) &d, NULL) ; \ - LAGraph_Free ((void **) &delta, NULL) ; \ + LAGraph_Free ((void **) &depth, NULL) ; \ + LAGraph_Free ((void **) &bc_vertex_flow, NULL) ; \ LAGraph_Free ((void **) &S, NULL) ; \ - LAGraph_Free ((void **) &sigma, NULL) ; \ + LAGraph_Free ((void **) &paths, NULL) ; \ LAGraph_Free ((void **) &Pj, NULL) ; \ LAGraph_Free ((void **) &Ptail, NULL) ; \ } @@ -63,10 +63,10 @@ int LG_check_edgeBetweennessCentrality double* result ; // Holds the distances (depth levels) from the source vertex. - int64_t *d = NULL ; + int64_t *depth = NULL ; - // Stores dependency scores for each vertex. - double *delta = NULL ; + // Stores the betweenness centrality for each vertex. + double *bc_vertex_flow = NULL ; // Stack used for backtracking phase int64_t *S = NULL ; @@ -78,7 +78,8 @@ int LG_check_edgeBetweennessCentrality GrB_Index *Ptail = NULL ; GrB_Index *Phead = NULL ; - double *sigma = NULL ; + // Holds the number of shortest paths for current node + double *paths = NULL ; //-------------------------------------------------------------------------- // check inputs @@ -114,25 +115,7 @@ int LG_check_edgeBetweennessCentrality LG_ASSERT_MSG (AT != NULL, LAGRAPH_NOT_CACHED, "G->AT is required") ; } - // better: as basic algo - /* - GrB_Matrix A_temp = NULL - if nself_edges is unknown - compute it - if any self edges - copy G->A into A_temp - remove self edges from A_temp - A = A_temp - now A has no self edges - when done: free A_temp - */ - - // hack: - // G->nself_edges = LAGRAPH_UNKNOWN ; <=== overkill - - // GRB_TRY (LAGraph_Graph_Print (G, 500, stdout, msg)) ; - // GxB_print (A, 5) ; - // GxB_print (AT, 5) ; + //-------------------------------------------------------------------------- LG_CLEAR_MSG ; @@ -141,9 +124,9 @@ int LG_check_edgeBetweennessCentrality // allocate workspace //-------------------------------------------------------------------------- - LG_TRY(LAGraph_Malloc((void **)&d, n, sizeof(int64_t), msg)); + LG_TRY(LAGraph_Malloc((void **)&depth, n, sizeof(int64_t), msg)); - LG_TRY(LAGraph_Calloc((void **)&delta, n, sizeof(double), msg)); + LG_TRY(LAGraph_Calloc((void **)&bc_vertex_flow, n, sizeof(double), msg)); LG_TRY(LAGraph_Malloc((void **)&S, n, sizeof(int64_t), msg)); @@ -162,7 +145,6 @@ int LG_check_edgeBetweennessCentrality // Initialize centrality matrix result to 0 // 1. result [(v, w)] ← 0, ∀(v, w) ∈ E - // TODO make this a copy of A except with 1 = 0 // A temporary result centrality matrix initialized to 0 for all vertice, // -- further changes would need to be made to make it a dictionary of edges. GrB_Index result_size = n * n ; @@ -186,17 +168,12 @@ int LG_check_edgeBetweennessCentrality Phead = ATp ; - // for (size_t i = 0; i < Aj_size; i++) { - // printf("%ld ", Aj[i]); - // } - // printf("\n"); - //-------------------------------------------------------------------------- LG_TRY(LAGraph_Malloc((void **)&Pj, nvals, sizeof(GrB_Index), msg)); LG_TRY(LAGraph_Malloc((void **)&Ptail, n, sizeof(GrB_Index), msg)); // might need to be + 1 - LAGraph_Calloc ((void **) &sigma, n, sizeof (double), msg) ; + LAGraph_Calloc ((void **) &paths, n, sizeof (double), msg) ; // 2. for ∀s ∈ V for (int64_t s = 0; s < n; s++) { @@ -204,23 +181,22 @@ int LG_check_edgeBetweennessCentrality size_t sp = 0; // Initialize predecessors list P[w] to empty - // TODO // 5. P [w] ← empty queue, ∀w ∈ V memcpy (Ptail, ATp, n * sizeof (GrB_Index)) ; - // Initialize sigma[t], d[t] for all t + // Initialize paths[t], d[t] for all t // 6. σ[t] ← 0, ∀t ∈ V , σ[s] ← 1 // Keeps track of the number of shortest paths for each vertex. for (int64_t i = 0; i < n; i++) { - sigma [i] = 0 ; + paths [i] = 0 ; } - sigma [s] = 1 ; + paths [s] = 1 ; // 7. d[t] ← −1, ∀t ∈ V , d[s] ← 0 for (size_t t = 0; t < n; t++) { - d [t] = -1; + depth [t] = -1; } - d [s] = 0; + depth [s] = 0; // Initialize queue and enqueue starting node s // 8. Q ← empty queue @@ -237,35 +213,25 @@ int LG_check_edgeBetweennessCentrality // 13. push(S, v) S [sp++] = v; - // if (s < 10) printf("v: %ld\n", v); - // TODO // traverse all entries in A(v,:) for (int64_t p = Ap [v] ; p < Ap [v+1] ; p++) { int64_t w = Aj [p] ; - // if (s < 10) printf(" w: %ld\n", w); - - // printf(" %ld - d[w] (%ld) < 0?\n", w, d [w]); // 16. if d[w] < 0 - if (d [w] < 0) { + if (depth [w] < 0) { // Update depth and enqueue // 18. enqueue(Q, w) queue [qt++] = w ; // 19. d[w] ← d[v] + 1 - d [w] = d [v] + 1 ; - // printf(" changed d[w] (%ld) = d[v] (%ld) + 1\n", d [w], d [v] ); - + depth [w] = depth [v] + 1 ; } - // printf(" %ld - d[w] (%ld) == d[v] (%ld) + 1?\n", w, d [w], d [v]); - // 20. if d[w] = d[v] + 1 - if (d [w] == d [v] + 1) { + if (depth [w] == depth [v] + 1) { // Update shortest path count and add predecessor // 22. σ[w] ← σ[w] + σ[v] - // printf(" sigma[w] (%g) = sigma[w] (%g) + sigma[v] (%g)\n", sigma[w]+sigma[v], sigma[w], sigma[v]); - sigma [w] = sigma [w] + sigma [v] ; + paths [w] = paths [w] + paths [v] ; // 23. append(P [w], v) if (Ptail [w] >= Phead [w+1] || Ptail [w] < Phead [w]) { @@ -274,131 +240,40 @@ int LG_check_edgeBetweennessCentrality fflush (stdout) ; abort ( ) ; } Pj [Ptail [w]++] = v ; - // printf(" added\n Pj: "); - - // for (int64_t p = Phead [w] ; p < Ptail [w] ; p++) { - // printf("%ld ", Pj [p]) ; - // } - // printf("\n"); - } } - // printf("\n"); - } - - if (s < 10) { - - printf("==========================================\n"); - - for (int64_t w = 0; w < n; w++) { - printf("%ld's: ", w); - for (int64_t p = Phead [w] ; p < Ptail [w] ; p++) { - printf("%ld ", Pj [p]) ; - } - printf("\n"); - } - printf("\n"); - } // Set dependency score δ[v] ← 0 // 24. δ[v] ← 0, ∀v ∈ V for (size_t v = 0; v < n; v++) { - delta [v] = 0 ; + bc_vertex_flow [v] = 0 ; } - - // PRINT OUT STUFF - // printf("==========================================\n"); - // printf("d:\n"); - // for (size_t i = 0; i < n; i++) { - // printf("(%ld, %ld) ", i, d[i]); - // } - // printf("\n"); - - // printf("delta:\n"); - // for (size_t i = 0; i < n; i++) { - // printf("(%ld, %g) ", i, delta[i]); - // } - // printf("\n"); - - // printf("S:\n"); - // for (size_t i = 0; i < n; i++) { - // printf("(%ld, %ld) ", i, S[i]); - // } - // printf("\n"); - - // printf("queue:\n"); - // for (size_t i = 0; i < n; i++) { - // printf("(%ld, %ld) ", i, queue[i]); - // } - // printf("\n"); - - // printf("sigma:\n"); - // for (size_t i = 0; i < n; i++) { - // printf("(%ld, %g) ", i, sigma[i]); - // } - // printf("\n==========================================\n"); - // printf("\n"); - - // Process stack S // 25. while ¬empty(S) while (sp > 0) { // 27. w ← pop(S) int64_t w = S [--sp] ; - if (s < 10) printf("w: %ld\n", w); - // 28. for v ∈ P [w] for (int64_t p = Phead [w] ; p < Ptail [w] ; p++) { int64_t v = Pj [p] ; - if (s < 10) printf(" v: %ld\n", v); // Update dependency and centrality values // 30. δ[v] ← δ[v] + σ[v] × ( δ[w]/σ[w] + 1) - if (s < 10) printf(" %g = %g * (%g + 1)/%g\n", sigma [v] * ((delta [w] + 1) / sigma [w]), sigma [v], delta [w], sigma [w]) ; - - // if (v == w) { printf ("Ack!!\n") ; fflush (stdout) ; abort ( ) ; } if (v == w) { printf (" Ack!!\n") ; fflush (stdout) ; abort ( ) ; continue; } - double centrality = sigma [v] * ((delta [w] + 1) / sigma [w]) ; - delta [v] += centrality ; + double centrality = paths [v] * ((bc_vertex_flow [w] + 1) / paths [w]) ; + bc_vertex_flow [v] += centrality ; // 31. result [(v, w)] ← result [(v, w)] + σ[v] × ( δ[w]/σ[w] + 1) result [INDEX (v,w)] += centrality; - // result [INDEX (w,v)] += centrality; - if (s < 10) printf(" result: %g\n", result[INDEX(v,w)]); - - // if (result [INDEX (v,w)] == 0) { - // result [INDEX (w,v)] += centrality; - // if (s < 10) printf(" result: %g\n", result[INDEX(w,v)]); - - // } - // else { - // result [INDEX (v,w)] += centrality; - // if (s < 10) printf(" result: %g\n", result[INDEX(v,w)]); - // } - - // if (s < 10) { - // printf(" delta:\n "); - // for (size_t i = 0; i < n; i++) { - // printf(" (%ld, %g)", i, delta[i]); - // } - // printf("\n"); - // printf(" sigma:\n "); - // for (size_t i = 0; i < n; i++) { - // printf(" (%ld, %g)", i, sigma[i]); - // } - // printf("\n"); - // } - } - // if (s < 10) printf("\n"); } @@ -422,20 +297,6 @@ int LG_check_edgeBetweennessCentrality &ATp, &ATj, &ATx, ATp_size, ATj_size, ATx_size, AT_iso, false, NULL)) ; #endif - // printf("result: \n") ; - // for (int64_t i = 0 ; i < n ; i++) - // { - // printf ("row: %ld\n", i) ; - // for (int64_t j = 0 ; j < n ; j++) - // { - // int64_t p = INDEX (i,j) ; - // double aij = result [p] ; - // printf(" C(%ld,%ld) = %g\n", i, j, aij) ; - // // numerical value of A(i,j) - // } - // printf("\n") ; - // } - #if 0 GrB_Info GxB_Matrix_pack_FullR // pack a full matrix, held by row ( @@ -454,8 +315,6 @@ GrB_Info GxB_Matrix_pack_FullR // pack a full matrix, held by row LG_TRY (GrB_assign(C_temp, A, NULL, C_temp, GrB_ALL, n, GrB_ALL, n, GrB_DESC_RS)) ; - GxB_print(C_temp, GxB_COMPLETE) ; - *C = C_temp; //-------------------------------------------------------------------------- diff --git a/experimental/test/test_edgeBetweennessCentrality.c b/experimental/test/test_edgeBetweennessCentrality.c index 2a11dcba0b..acec3cced8 100644 --- a/experimental/test/test_edgeBetweennessCentrality.c +++ b/experimental/test/test_edgeBetweennessCentrality.c @@ -56,7 +56,7 @@ double difference(GrB_Matrix bc, double* gap_result, GrB_Index rows, GrB_Index c // GxB_print (diff, 5) ; OK(GrB_apply(diff, NULL, NULL, GrB_ABS_FP64, diff, NULL)); - float err = 0; + double err = 0; OK(GrB_reduce(&err, NULL, GrB_MAX_MONOID_FP64, diff, NULL)); OK(GrB_free(&diff)); @@ -139,8 +139,10 @@ void test_diamonds_ebc (void) { LAGraph_Init (msg) ; GrB_Matrix A = NULL ; + GrB_Matrix AT = NULL ; GrB_Matrix centrality = NULL ; int niters = 0 ; + LAGraph_Kind kind = LAGraph_ADJACENCY_DIRECTED ; // create the karate graph snprintf (filename, LEN, LG_DATA_DIR "%s", "diamonds.mtx") ; @@ -148,21 +150,31 @@ void test_diamonds_ebc (void) TEST_CHECK (f != NULL) ; OK (LAGraph_MMRead (&A, f, msg)) ; OK (fclose (f)) ; - OK (LAGraph_New (&G, &A, LAGraph_ADJACENCY_DIRECTED, msg)) ; + OK (LAGraph_New (&G, &A, kind, msg)) ; TEST_CHECK (A == NULL) ; // A has been moved into G->A - // compute its betweenness centrality - OK (LG_check_edgeBetweennessCentrality (¢rality, G, msg)) ; + // check that AT is cached + int ok_result = (kind == LAGraph_ADJACENCY_UNDIRECTED) ? + LAGRAPH_CACHE_NOT_NEEDED : GrB_SUCCESS ; + int result = LAGraph_Cached_AT (G, msg) ; + TEST_CHECK (result == ok_result) ; - // compare with GAP: + // compute its betweenness centrality with C version + OK (LG_check_edgeBetweennessCentrality (¢rality, G, msg)) ; double err = difference(centrality, &diamonds_ebc[0][0], 8, 8) ; - printf ("diamonds: err: %e\n", err) ; + printf ("diamonds: err: %e (C version)\n", err) ; TEST_CHECK (err < 1e-4) ; OK (GrB_free (¢rality)) ; - OK (LAGraph_Delete (&G, msg)) ; - LAGraph_Finalize (msg) ; + // compute its betweenness centrality with GraphBLAS version + OK (LAGr_EdgeBetweennessCentrality (¢rality, G, msg)) ; + err = difference(centrality, &diamonds_ebc[0][0], 8, 8) ; + printf ("diamonds: err: %e (pure GraphBLAS)\n", err) ; + TEST_CHECK (err < 1e-4) ; + OK (GrB_free (¢rality)) ; + OK (LAGraph_Delete (&G, msg)) ; + LAGraph_Finalize (msg) ; } //------------------------------------------------------------------------------ @@ -185,16 +197,21 @@ void test_karate_ebc (void) OK (LAGraph_New (&G, &A, LAGraph_ADJACENCY_UNDIRECTED, msg)) ; TEST_CHECK (A == NULL) ; // A has been moved into G->A - // compute its betweenness centrality + // compute its betweenness centrality (C version) OK (LG_check_edgeBetweennessCentrality (¢rality, G, msg)) ; + double err = difference(centrality, &karate_ebc[0][0], 34, 34) ; + printf ("karate: err: %e (C version)\n", err) ; + TEST_CHECK (err < 1e-4) ; + OK (GrB_free (¢rality)) ; - // compare with GAP: - float err = difference(centrality, &karate_ebc[0][0], 34, 34) ; - printf ("karate: err: %e\n", err) ; + // compute its betweenness centrality (GraphBLAS version) + OK (LAGr_EdgeBetweennessCentrality (¢rality, G, msg)) ; + err = difference(centrality, &karate_ebc[0][0], 34, 34) ; + printf ("karate: err: %e (GraphBLAS version)\n", err) ; TEST_CHECK (err < 1e-4) ; OK (GrB_free (¢rality)) ; - OK (LAGraph_Delete (&G, msg)) ; + OK (LAGraph_Delete (&G, msg)) ; LAGraph_Finalize (msg) ; } @@ -203,8 +220,11 @@ void test_karate_ebc (void) // list of tests //------------------------------------------------------------------------------ +// FIXME: add more matrices + TEST_LIST = { {"test_diamonds_ebc", test_diamonds_ebc}, {"test_karate_ebc", test_karate_ebc}, +// {"test_many", test_many}, FIXME ADD THIS {NULL, NULL} }; diff --git a/include/LAGraphX.h b/include/LAGraphX.h index e02214dffa..148c47007c 100644 --- a/include/LAGraphX.h +++ b/include/LAGraphX.h @@ -1276,6 +1276,20 @@ int LAGr_HITS char *msg ) ; +//------------------------------------------------------------------------------ +// edge betweenness centrality +//------------------------------------------------------------------------------ + +LAGRAPHX_PUBLIC +int LAGr_EdgeBetweennessCentrality +( + // output: + GrB_Matrix *centrality, // centrality(i): betweeness centrality of i + // input: + LAGraph_Graph G, // input graph + char *msg +); + //------------------------------------------------------------------------------ // graph clustering with quality metrics //------------------------------------------------------------------------------