Line data Source code
1 : //------------------------------------------------------------------------------
2 : // LAGraph/experimental/test/test_Fiedler.c: test the HDIP Fiedler method
3 : //------------------------------------------------------------------------------
4 :
5 : // LAGraph, (c) 2022 by The LAGraph Contributors, All Rights Reserved.
6 : // SPDX-License-Identifier: BSD-2-Clause
7 : // See additional acknowledgments in the LICENSE file,
8 : // or contact permission@sei.cmu.edu for the full terms.
9 :
10 : // Contributed by Timothy A Davis, Texas A&M University
11 :
12 : //------------------------------------------------------------------------------
13 :
14 : #include <acutest.h>
15 : #include <LAGraph_test.h>
16 : #include "LG_Xtest.h"
17 : #include "LG_internal.h"
18 :
19 : typedef struct
20 : {
21 : const char *name ;
22 : double lambda ;
23 : double *fiedler ;
24 : }
25 : matrix_info ;
26 :
27 : float difference (GrB_Vector centrality, double *matlab_result) ;
28 :
29 3 : float difference (GrB_Vector centrality, double *matlab_result)
30 : {
31 3 : GrB_Vector diff = NULL, cmatlab = NULL ;
32 3 : GrB_Index n = 0 ;
33 3 : GrB_Vector_size (&n, centrality) ;
34 3 : GrB_Vector_new (&cmatlab, GrB_FP32, n) ;
35 164 : for (int i = 0 ; i < n ; i++)
36 : {
37 161 : GrB_Vector_setElement_FP64 (cmatlab, matlab_result [i], i) ;
38 : }
39 : // diff = max (abs (cmatlab - centrality))
40 3 : GrB_Vector_new (&diff, GrB_FP32, n) ;
41 3 : GrB_eWiseAdd (diff, NULL, NULL, GrB_MINUS_FP32, cmatlab, centrality,
42 : NULL) ;
43 3 : GrB_apply (diff, NULL, NULL, GrB_ABS_FP32, diff, NULL) ;
44 3 : float err = 0 ;
45 3 : GrB_reduce (&err, NULL, GrB_MAX_MONOID_FP32, diff, NULL) ;
46 3 : GrB_free (&diff) ;
47 3 : GrB_free (&cmatlab) ;
48 3 : return (err) ;
49 : }
50 :
51 : // Bucky test data
52 : double bucky_fiedler [60] = {
53 : -0.2236,
54 : -0.2071,
55 : -0.1804,
56 : -0.1804,
57 : -0.2071,
58 : -0.2022,
59 : -0.1669,
60 : -0.1098,
61 : -0.1098,
62 : -0.1669,
63 : -0.1669,
64 : -0.1481,
65 : -0.0744,
66 : -0.0477,
67 : -0.1049,
68 : -0.1098,
69 : -0.0744,
70 : 0.0094,
71 : 0.0259,
72 : -0.0477,
73 : -0.1098,
74 : -0.0477,
75 : 0.0259,
76 : 0.0094,
77 : -0.0744,
78 : -0.1669,
79 : -0.1049,
80 : -0.0477,
81 : -0.0744,
82 : -0.1481,
83 : 0.1481,
84 : 0.0745,
85 : 0.0477,
86 : 0.1049,
87 : 0.1669,
88 : 0.0745,
89 : -0.0094,
90 : -0.0259,
91 : 0.0477,
92 : 0.1098,
93 : 0.0477,
94 : -0.0259,
95 : -0.0094,
96 : 0.0745,
97 : 0.1098,
98 : 0.1049,
99 : 0.0477,
100 : 0.0745,
101 : 0.1481,
102 : 0.1669,
103 : 0.1669,
104 : 0.1098,
105 : 0.1098,
106 : 0.1669,
107 : 0.2022,
108 : 0.2071,
109 : 0.1804,
110 : 0.1804,
111 : 0.2071,
112 : 0.2236} ;
113 :
114 : // Karate test data
115 : double karate_fiedler [34] = {
116 : -0.3561,
117 : -0.1036,
118 : -0.0156,
119 : -0.1243,
120 : -0.2280,
121 : -0.2097,
122 : -0.2097,
123 : -0.1224,
124 : 0.0163,
125 : 0.1108,
126 : -0.2280,
127 : -0.2463,
128 : -0.1853,
129 : -0.0725,
130 : 0.1900,
131 : 0.1900,
132 : -0.1548,
133 : -0.1749,
134 : 0.1900,
135 : -0.0741,
136 : 0.1900,
137 : -0.1749,
138 : 0.1900,
139 : 0.1792,
140 : 0.1703,
141 : 0.1794,
142 : 0.2155,
143 : 0.1428,
144 : 0.1002,
145 : 0.1937,
146 : 0.0732,
147 : 0.0790,
148 : 0.1427,
149 : 0.1274} ;
150 :
151 : double west0067_fiedler [67] = {
152 : -0.7918,
153 : -0.0506,
154 : -0.0329,
155 : -0.0366,
156 : -0.1569,
157 : -0.1608,
158 : -0.1776,
159 : -0.1747,
160 : -0.1529,
161 : -0.0391,
162 : -0.0320,
163 : -0.0046,
164 : -0.1138,
165 : -0.0140,
166 : -0.0314,
167 : -0.0114,
168 : -0.0066,
169 : -0.0862,
170 : 0.0245,
171 : -0.0117,
172 : 0.0232,
173 : 0.0338,
174 : 0.0052,
175 : 0.0185,
176 : -0.0731,
177 : -0.0520,
178 : -0.0602,
179 : -0.0711,
180 : -0.0623,
181 : 0.0528,
182 : -0.0016,
183 : 0.0447,
184 : 0.0566,
185 : 0.0444,
186 : 0.0610,
187 : 0.0220,
188 : -0.0008,
189 : 0.0171,
190 : 0.0305,
191 : 0.0519,
192 : 0.0414,
193 : 0.0491,
194 : 0.0482,
195 : 0.0912,
196 : 0.0660,
197 : 0.1074,
198 : 0.1016,
199 : 0.1078,
200 : 0.0683,
201 : 0.0871,
202 : 0.0777,
203 : 0.0839,
204 : 0.0901,
205 : 0.1092,
206 : 0.0850,
207 : 0.0752,
208 : -0.0019,
209 : 0.0239,
210 : 0.0442,
211 : 0.0772,
212 : -0.0179,
213 : 0.0770,
214 : 0.1072,
215 : 0.0342,
216 : 0.0762,
217 : 0.1115,
218 : 0.1000} ;
219 :
220 : const matrix_info files [ ] =
221 : {
222 : { "bucky.mtx", 0.2434, bucky_fiedler },
223 : { "karate.mtx", 1.3297, karate_fiedler },
224 : { "west0067.mtx", 6.5586, west0067_fiedler },
225 : { "", 0, NULL },
226 : } ;
227 :
228 1 : void test_fiedler (void)
229 : {
230 : #if LAGRAPH_SUITESPARSE
231 :
232 : //--------------------------------------------------------------------------
233 : // startup LAGraph and GraphBLAS
234 : //--------------------------------------------------------------------------
235 :
236 : char msg[LAGRAPH_MSG_LEN]; // for error messages from LAGraph
237 1 : LAGraph_Graph G = NULL;
238 1 : GrB_Matrix Y = NULL;
239 1 : GrB_Matrix A = NULL;
240 1 : GrB_Matrix indiag = NULL;
241 1 : GrB_Vector x = NULL;
242 : // Additional variables and modifications needed to test MYPCG2
243 1 : GrB_Vector steper = NULL;
244 1 : GrB_Vector u = NULL; // a vector of size nrowsLap, filled with 1.
245 : // set u[0] = 1+sqrt(nrowsLap)
246 : // Additional variables needed to test Hdip
247 1 : GrB_Vector iters = NULL;
248 1 : float lambda_result = 0;
249 1 : GrB_Vector fiedler_vector = NULL;
250 1 : GrB_Vector kmax = NULL;
251 : #define LEN 512
252 : char filename [LEN+1] ;
253 : GrB_Index n ;
254 :
255 : // start GraphBLAS and LAGraph
256 1 : bool burble = false; // set true for diagnostic outputs
257 1 : LAGraph_Init (msg) ;
258 :
259 : //--------------------------------------------------------------------------
260 : // read in the graphs and test them
261 : //--------------------------------------------------------------------------
262 :
263 1 : for (int k = 0 ; ; k++)
264 3 : {
265 :
266 : // load the matrix as A
267 4 : const char *aname = files [k].name ;
268 4 : if (strlen (aname) == 0) break;
269 3 : printf ("\n %s: ==================================\n", aname) ;
270 3 : TEST_CASE (aname) ;
271 3 : snprintf (filename, LEN, LG_DATA_DIR "%s", aname) ;
272 3 : FILE *f = fopen (filename, "r") ;
273 3 : TEST_CHECK (f != NULL) ;
274 3 : OK (LAGraph_MMRead (&A, f, msg)) ;
275 3 : TEST_MSG ("Loading of adjacency matrix failed") ;
276 3 : fclose (f) ;
277 :
278 : // set all entries to 1
279 3 : OK (GrB_Matrix_nrows (&n, A)) ;
280 3 : OK (GrB_assign(A, A, NULL, (double)1,
281 : GrB_ALL, n, GrB_ALL, n, GrB_DESC_S));
282 :
283 : // ensure A is symmetric and remove self edges
284 3 : OK (GrB_eWiseAdd (A, NULL, NULL, GrB_ONEB_FP64, A, A, GrB_DESC_T1)) ;
285 3 : OK (GrB_select (A, NULL, NULL, GrB_OFFDIAG, A, 0, NULL)) ;
286 :
287 : //----------------------------------------------------------------------
288 : // try the LAGraph_Laplacian algorithm
289 : //----------------------------------------------------------------------
290 :
291 : // Variables needed to test Laplacian
292 : float infnorm;
293 3 : OK (LAGraph_Laplacian(&Y, &infnorm, A, msg));
294 :
295 3 : printf("\n===========================The laplacian matrix: \n");
296 3 : OK (LAGraph_Matrix_Print(Y, LAGraph_SHORT, stdout, msg)) ;
297 :
298 : //----------------------------------------------------------------------
299 : // try the LAGraph_mypcg2 algorithm
300 : //----------------------------------------------------------------------
301 :
302 : GrB_Index kk;
303 : GrB_Index nrows;
304 : float nrowsLap; // number of rows of laplacian matrix
305 : float alpha;
306 3 : OK (GrB_Matrix_nrows(&nrows, Y));
307 :
308 3 : nrowsLap = (float)n;
309 :
310 3 : OK (GrB_Vector_new(&u, GrB_FP32, n));
311 : // u = all ones vector
312 3 : OK (GrB_assign(u, NULL, NULL, 1, GrB_ALL, n, NULL));
313 : // u [0] = 1+sqrt(n)
314 3 : OK (GrB_Vector_setElement_FP32(u, 1 + sqrt(nrowsLap), 0));
315 :
316 3 : alpha = nrowsLap + sqrt(nrowsLap);
317 :
318 3 : OK (GrB_Matrix_new(&indiag, GrB_FP32, n, n));
319 3 : OK (GrB_select(indiag, NULL, NULL, GrB_DIAG, Y, 0, NULL));
320 3 : OK (GrB_apply(indiag, NULL, NULL, GrB_MINV_FP32, indiag, NULL));
321 :
322 3 : OK (GrB_Vector_new(&x, GrB_FP32, n));
323 3 : OK (GrB_assign(x, NULL, NULL, 1, GrB_ALL, n, NULL));
324 3 : OK (GrB_Vector_setElement_FP32(x, 0, 0));
325 :
326 3 : OK (LAGraph_mypcg2 (&steper, &kk, Y, u, alpha, indiag, x, .000001, 50,
327 : msg)) ;
328 :
329 : //--------------------------------------------------------------------------
330 : // try the LAGraph_Hdip_Fiedler algorithm
331 : //--------------------------------------------------------------------------
332 :
333 : //Set kmax = [20,50]
334 3 : OK (GrB_Vector_new(&kmax, GrB_FP32, 2));
335 3 : OK (GrB_Vector_setElement_FP32(kmax, 20, 0));
336 3 : OK (GrB_Vector_setElement_FP32(kmax, 50, 1));
337 :
338 3 : OK (LAGraph_Hdip_Fiedler (&iters, &lambda_result,
339 : &fiedler_vector, Y, infnorm, kmax, 0.000001, 0.000001, msg)) ;
340 :
341 : //--------------------------------------------------------------------------
342 : // check the results
343 : //--------------------------------------------------------------------------
344 :
345 3 : float err = difference (fiedler_vector, files [k].fiedler) ;
346 3 : TEST_CHECK (err < 1e-4) ;
347 3 : err = fabs(lambda_result - files [k].lambda) ;
348 3 : TEST_CHECK (err < 1e-4) ;
349 :
350 : //--------------------------------------------------------------------------
351 : // print the results
352 : //--------------------------------------------------------------------------
353 :
354 3 : printf("\n===============================The result vector x:\n");
355 3 : LAGraph_Vector_Print (fiedler_vector, 3, stdout, msg) ;
356 3 : printf("\n===============================The lambda: %f\n", lambda_result);
357 3 : printf("\n===============================The iters: \n");
358 3 : LAGraph_Vector_Print (iters, 3, stdout, msg) ;
359 :
360 : //--------------------------------------------------------------------------
361 : // free everyting and finish
362 : //--------------------------------------------------------------------------
363 :
364 3 : GrB_free (&A) ;
365 3 : GrB_free (&Y) ;
366 3 : GrB_free (&x) ;
367 3 : GrB_free (&fiedler_vector) ;
368 3 : GrB_free (&u) ;
369 3 : GrB_free (&steper) ;
370 3 : GrB_free (&indiag) ;
371 3 : GrB_free (&iters) ;
372 3 : GrB_free (&kmax) ;
373 : }
374 :
375 1 : OK (LAGraph_Finalize (msg)) ;
376 : #endif
377 1 : }
378 :
379 : TEST_LIST = {
380 : {"Fieder", test_fiedler},
381 : {NULL, NULL}
382 : } ;
383 :
|