1 //------------------------------------------------------------------------------
2 // irowscale: scale the rows of an adjacency matrix by out-degree
3 //------------------------------------------------------------------------------
4
5 // SuiteSparse:GraphBLAS, Timothy A. Davis, (c) 2017-2021, All Rights Reserved.
6 // SPDX-License-Identifier: Apache-2.0
7
8 //------------------------------------------------------------------------------
9
10 // on input, A is a square unsymmetric binary matrix of size n-by-n, of any
11 // built-in type. On output, C is a rowscaled version of A, of type
12 // GrB_UINT64, with C = D*A + I. The diagonal matrix D has diagonal entries
13 // D(i,i)=(2^30)//sum(A(i,:)), or D(i,i) is not present if A(i,:) is empty.
14 // The matrix I has entries only on the diagonal, all equal to zero. This
15 // optional step ensures that C has no empty columns, which speeds up the
16 // subsequent PageRank computation.
17
18 /* MATLAB equivalent (excluding the addition of I):
19
20 function C = rowscale (A)
21 %ROWSCALE row scale an adjacency matrix by out-degree
22 % C = rowscale (A)
23
24 % scale the adjacency matrix by out-degree
25 dout = sum (A,2) ; % dout(i) is the out-degree of node i
26 is_nonempty = (dout > 0) ; % find vertices with outgoing edges
27 nonempty = find (is_nonempty) ; % list of vertices with outgoing edges
28 n = size (A,1) ;
29
30 % divide each non-empty row by its out-degree
31 dinv = 1 ./ dout (is_nonempty) ;
32 D = sparse (nonempty, nonempty, dinv, n, n) ;
33 C = D*A ;
34
35 C = floor ((2^30) * C) ; % scale the result to integer
36 */
37
38 #include "GraphBLAS.h"
39
40 //------------------------------------------------------------------------------
41 // helper macros
42 //------------------------------------------------------------------------------
43
44 // free all workspace
45 #define FREEWORK \
46 { \
47 GrB_Vector_free (&dout) ; \
48 GrB_Matrix_free (&D) ; \
49 if (I != NULL) free (I) ; \
50 if (X != NULL) free (X) ; \
51 }
52
53 // error handler: free workspace and the output matrix C
54 #define FREE_ALL \
55 { \
56 GrB_Matrix_free (&C) ; \
57 FREEWORK ; \
58 }
59
60 #undef GB_PUBLIC
61 #define GB_LIBRARY
62 #include "graphblas_demos.h"
63
64 //------------------------------------------------------------------------------
65 // irowscale: C = D*A + I*0 where D(i,i) = ZSCALE/sum(A(i,:)
66 //------------------------------------------------------------------------------
67
68 GB_PUBLIC
irowscale(GrB_Matrix * Chandle,GrB_Matrix A)69 GrB_Info irowscale // GrB_SUCCESS or error condition
70 (
71 GrB_Matrix *Chandle, // output matrix C = rowscale (A)
72 GrB_Matrix A // input matrix, not modified
73 )
74 {
75
76 //--------------------------------------------------------------------------
77 // intializations
78 //--------------------------------------------------------------------------
79
80 GrB_Info info ;
81 GrB_Vector dout = NULL ;
82 GrB_Matrix D = NULL, C = NULL ;
83 GrB_Index *I = NULL ;
84 uint64_t *X = NULL ;
85
86 // n = size (A,1) ;
87 GrB_Index n ;
88 OK (GrB_Matrix_nrows (&n, A)) ;
89
90 //--------------------------------------------------------------------------
91 // dout = sum (A,2) ; // dout(i) is the out-degree of node i
92 //--------------------------------------------------------------------------
93
94 OK (GrB_Vector_new (&dout, GrB_UINT64, n)) ;
95 OK (GrB_Matrix_reduce_Monoid (dout, NULL, NULL, GrB_PLUS_MONOID_UINT64,
96 A, NULL)) ;
97
98 //--------------------------------------------------------------------------
99 // construct scaling matrix D
100 //--------------------------------------------------------------------------
101
102 // D is diagonal with D(i,i) = ZSCALE/sum(A(i,:)), or D(i,i) is not present
103 // if row i of A has no entries.
104
105 // [I,~,X] = find (dout) ;
106 I = (GrB_Index *) malloc ((n+1) * sizeof (GrB_Index)) ;
107 X = (uint64_t *) malloc ((n+1) * sizeof (uint64_t)) ;
108 CHECK (I != NULL && X != NULL, GrB_OUT_OF_MEMORY) ;
109 GrB_Index nvals = n ;
110 OK (GrB_Vector_extractTuples_UINT64 (I, X, &nvals, dout)) ;
111
112 // I and X exclude empty columns of A. This condition is always true.
113 CHECK (nvals <= n, GrB_INVALID_VALUE) ;
114
115 // D = diag (ZSCALE./dout) ;
116 OK (GrB_Matrix_new (&D, GrB_UINT64, n, n)) ;
117 for (int64_t i = 0 ; i < nvals ; i++)
118 {
119 // X (i) = ZSCALE / X (i), but make sure it doesn't underflow to zero.
120 // Underflow would only occur if a node has degree higher than 2^30.
121 uint64_t y = ZSCALE / X [i] ;
122 if (y == 0) y = 1 ;
123 X [i] = y ;
124 }
125 OK (GrB_Matrix_build_UINT64 (D, I, I, X, nvals, GrB_PLUS_UINT64)) ;
126
127 //--------------------------------------------------------------------------
128 // C = diagonal matrix with explicit zeros on diagonal
129 //--------------------------------------------------------------------------
130
131 // This step is optional, but it ensures that C has no non-empty columns,
132 // which speeds up the pagerank computation by ensuring r*C remains a dense
133 // vector.
134 for (int64_t i = 0 ; i < n ; i++)
135 {
136 I [i] = i ;
137 X [i] = 0 ;
138 }
139 OK (GrB_Matrix_new (&C, GrB_UINT64, n, n)) ;
140 OK (GrB_Matrix_build_UINT64 (C, I, I, X, n, GrB_PLUS_UINT64)) ;
141
142 //--------------------------------------------------------------------------
143 // C += D*A
144 //--------------------------------------------------------------------------
145
146 OK (GrB_mxm (C, NULL, GrB_PLUS_UINT64, GxB_PLUS_TIMES_UINT64, D, A, NULL)) ;
147
148 //--------------------------------------------------------------------------
149 // free workspace and return result
150 //--------------------------------------------------------------------------
151
152 FREEWORK ;
153 (*Chandle) = C ;
154 return (GrB_SUCCESS) ;
155 }
156
157