@@ -249,54 +249,50 @@ static GrB_Info LAGraph_RPQMatrixKleene(RPQMatrixPlan *plan, char *msg)
249249 // Creating identity matrix.
250250 GrB_Index n ;
251251 GRB_TRY (GrB_Matrix_nrows (& n , B ));
252- GrB_Matrix E ;
253- GRB_TRY (GrB_Matrix_new (& E , GrB_BOOL , n , n ));
252+ GrB_Matrix I ;
253+ GRB_TRY (GrB_Matrix_new (& I , GrB_BOOL , n , n ));
254254
255255 GrB_Vector v ;
256256 GRB_TRY (GrB_Vector_new (& v , GrB_BOOL , n ));
257257 GRB_TRY (GrB_Vector_assign_BOOL (v , NULL , NULL , true, GrB_ALL , n , NULL ));
258258
259- GRB_TRY (GrB_Matrix_diag (& E , v , 0 ));
259+ GRB_TRY (GrB_Matrix_diag (& I , v , 0 ));
260260
261261 GRB_TRY (GrB_Vector_free (& v ));
262262
263- // B + E
264- GrB_Matrix BPE ;
265- GRB_TRY (GrB_Matrix_new (& BPE , GrB_BOOL , n , n ));
266- GRB_TRY (GrB_eWiseAdd (BPE , GrB_NULL , GrB_NULL ,
267- op , B , E , GrB_DESC_R ));
268- // S <- S x (B + E)
269- GrB_Matrix S , T ;
270- GRB_TRY (GrB_Matrix_dup (& S , E ));
263+ // B + I
264+ GrB_Matrix BPI ;
265+ GRB_TRY (GrB_Matrix_new (& BPI , GrB_BOOL , n , n ));
266+ GRB_TRY (GrB_eWiseAdd (BPI , GrB_NULL , GrB_NULL ,
267+ op , B , I , GrB_DESC_R ));
268+ // S <- I
269+ GrB_Matrix S ;
270+ GRB_TRY (GrB_Matrix_dup (& S , I ));
271271
272272 bool changed = true;
273- GrB_Index nnz_S = n , nnz_T = 0 ;
273+ GrB_Index nnz_S = n , nnz_Sold = 0 ;
274274
275275 while (changed )
276276 {
277- // T = S * (B + E)
278- GRB_TRY (GrB_Matrix_new (& T , GrB_BOOL , n , n ));
279- GRB_TRY (GrB_mxm (T , GrB_NULL , GrB_NULL ,
280- sr , S , BPE , GrB_DESC_R ));
277+ // S <- S x (B + I)
278+ GRB_TRY (GrB_mxm (S , GrB_NULL , GrB_NULL ,
279+ sr , S , BPI , GrB_DESC_R ));
281280
282- GRB_TRY (GrB_Matrix_nvals (& nnz_T , T ));
283- if (nnz_T != nnz_S )
281+ GRB_TRY (GrB_Matrix_nvals (& nnz_S , S ));
282+ if (nnz_S != nnz_Sold )
284283 {
285284 changed = true;
286- nnz_S = nnz_T ;
287- GRB_TRY (GrB_Matrix_free (& S ));
288- S = T ;
285+ nnz_Sold = nnz_S ;
289286 }
290287 else
291288 {
292289 changed = false;
293- GRB_TRY (GrB_Matrix_free (& T ));
294290 }
295291 }
296292 plan -> res_mat = S ;
297293
298- GRB_TRY (GrB_Matrix_free (& E ));
299- GRB_TRY (GrB_Matrix_free (& BPE ));
294+ GRB_TRY (GrB_Matrix_free (& I ));
295+ GRB_TRY (GrB_Matrix_free (& BPI ));
300296 return (GrB_SUCCESS );
301297}
302298
@@ -310,8 +306,8 @@ static GrB_Info LAGraph_RPQMatrixKleene(RPQMatrixPlan *plan, char *msg)
310306// ┌─┬─┴─┬─┐
311307// │*│ │b│
312308// └┬┘ └─┘
313- // ┌┴┐
314- // │a│
309+ // ┌┴┐
310+ // │a│
315311// └─┘
316312// If matrix B is sparse and A is dense, then instead of naive
317313// way:
@@ -328,7 +324,7 @@ static GrB_Info LAGraph_RPQMatrixKleene(RPQMatrixPlan *plan, char *msg)
328324// └─┬─┘
329325// ┌─┬─┴─┬─┐
330326// │a│ │b│
331- // └─┘ └─┘
327+ // └─┘ └─┘
332328static GrB_Info LAGraph_RPQMatrixKleene_L (RPQMatrixPlan * plan , char * msg )
333329{
334330 LG_ASSERT (plan != NULL , GrB_NULL_POINTER );
@@ -347,44 +343,56 @@ static GrB_Info LAGraph_RPQMatrixKleene_L(RPQMatrixPlan *plan, char *msg)
347343 GrB_Matrix A = lhs -> res_mat ;
348344 GrB_Matrix B = rhs -> res_mat ;
349345
346+ // Creating identity matrix.
350347 GrB_Index n ;
351348 GRB_TRY (GrB_Matrix_nrows (& n , A ));
349+ GrB_Matrix I ;
350+ GRB_TRY (GrB_Matrix_new (& I , GrB_BOOL , n , n ));
351+
352+ GrB_Vector v ;
353+ GRB_TRY (GrB_Vector_new (& v , GrB_BOOL , n ));
354+ GRB_TRY (GrB_Vector_assign_BOOL (v , NULL , NULL , true, GrB_ALL , n , NULL ));
355+
356+ GRB_TRY (GrB_Matrix_diag (& I , v , 0 ));
357+
358+ GRB_TRY (GrB_Vector_free (& v ));
359+
360+ // B + I
361+ GrB_Matrix API ;
362+ GRB_TRY (GrB_Matrix_new (& API , GrB_BOOL , n , n ));
363+ GRB_TRY (GrB_eWiseAdd (API , GrB_NULL , GrB_NULL ,
364+ op , A , I , GrB_DESC_R ));
352365
353366 // S <- B
354- GrB_Matrix S , T ;
367+ GrB_Matrix S ;
355368 GRB_TRY (GrB_Matrix_dup (& S , B ));
356- GRB_TRY (GrB_Matrix_new (& T , GrB_BOOL , n , n ));
357369
358370 bool changed = true;
359- GrB_Index nnz_S = 0 , nnz_T = 0 ;
371+ GrB_Index nnz_S = 0 , nnz_Sold = 0 ;
360372
361373 while (changed )
362374 {
363- // T <- A x S
364- GRB_TRY (GrB_mxm (T , NULL , NULL , sr , A , S , NULL ));
365-
366- // S <- S + T
367- GRB_TRY (GrB_eWiseAdd (S , NULL , NULL , op , S , T , NULL ));
375+ // T <- (A + I) x S
376+ GRB_TRY (GrB_mxm (S , NULL , NULL , sr , API , S , NULL ));
368377
369378 GRB_TRY (GrB_Matrix_nvals (& nnz_S , S ));
370- if (nnz_S == nnz_T )
379+ if (nnz_S != nnz_Sold )
371380 {
372- changed = false;
381+ changed = true;
382+ nnz_Sold = nnz_S ;
373383 }
374384 else
375385 {
376- changed = true;
377- nnz_T = nnz_S ;
378-
386+ changed = false;
379387 }
380388 }
381389
382390 plan -> res_mat = S ;
383- GRB_TRY (GrB_Matrix_free (& T ));
391+ GRB_TRY (GrB_Matrix_free (& I ));
392+ GRB_TRY (GrB_Matrix_free (& API ));
384393 return GrB_SUCCESS ;
385394}
386395
387-
388396// this function need to handle special case where some optimization
389397// are available.
390398// consider following AST:
@@ -412,7 +420,7 @@ static GrB_Info LAGraph_RPQMatrixKleene_L(RPQMatrixPlan *plan, char *msg)
412420// └─┬─┘
413421// ┌─┬─┴─┬─┐
414422// │a│ │b│
415- // └─┘ └─┘
423+ // └─┘ └─┘
416424
417425static GrB_Info LAGraph_RPQMatrixKleene_R (RPQMatrixPlan * plan , char * msg )
418426{
@@ -432,39 +440,53 @@ static GrB_Info LAGraph_RPQMatrixKleene_R(RPQMatrixPlan *plan, char *msg)
432440 GrB_Matrix A = lhs -> res_mat ;
433441 GrB_Matrix B = rhs -> res_mat ;
434442
443+ // Creating identity matrix.
435444 GrB_Index n ;
436- GRB_TRY (GrB_Matrix_nrows (& n , A ));
445+ GRB_TRY (GrB_Matrix_nrows (& n , B ));
446+ GrB_Matrix I ;
447+ GRB_TRY (GrB_Matrix_new (& I , GrB_BOOL , n , n ));
448+
449+ GrB_Vector v ;
450+ GRB_TRY (GrB_Vector_new (& v , GrB_BOOL , n ));
451+ GRB_TRY (GrB_Vector_assign_BOOL (v , NULL , NULL , true, GrB_ALL , n , NULL ));
452+
453+ GRB_TRY (GrB_Matrix_diag (& I , v , 0 ));
454+
455+ GRB_TRY (GrB_Vector_free (& v ));
456+
457+ // B + I
458+ GrB_Matrix BPI ;
459+ GRB_TRY (GrB_Matrix_new (& BPI , GrB_BOOL , n , n ));
460+ GRB_TRY (GrB_eWiseAdd (BPI , GrB_NULL , GrB_NULL ,
461+ op , B , I , GrB_DESC_R ));
437462
438463 // S <- A
439- GrB_Matrix S , T ;
464+ GrB_Matrix S ;
440465 GRB_TRY (GrB_Matrix_dup (& S , A ));
441- GRB_TRY (GrB_Matrix_new (& T , GrB_BOOL , n , n ));
442466
443467 bool changed = true;
444- GrB_Index nnz_S = 0 , nnz_T = 0 ;
468+ GrB_Index nnz_S = 0 , nnz_Sold = 0 ;
445469
446470 while (changed )
447471 {
448- // T <- S × B
449- GRB_TRY (GrB_mxm (T , NULL , NULL , sr , S , B , NULL ));
450-
451- // S <- S + T
452- GRB_TRY (GrB_eWiseAdd (S , NULL , NULL , op , S , T , NULL ));
472+ // S <- S x (B + I)
473+ GRB_TRY (GrB_mxm (S , NULL , NULL , sr , S , BPI , NULL ));
453474
454475 GRB_TRY (GrB_Matrix_nvals (& nnz_S , S ));
455- if (nnz_S == nnz_T )
476+ if (nnz_S != nnz_Sold )
456477 {
457- changed = false;
478+ changed = true;
479+ nnz_Sold = nnz_S ;
458480 }
459481 else
460482 {
461- changed = true;
462- nnz_T = nnz_S ;
483+ changed = false;
463484 }
464485 }
465486
466487 plan -> res_mat = S ;
467- GRB_TRY (GrB_Matrix_free (& T ));
488+ GRB_TRY (GrB_Matrix_free (& I ));
489+ GRB_TRY (GrB_Matrix_free (& BPI ));
468490 return GrB_SUCCESS ;
469491}
470492
@@ -479,7 +501,7 @@ GrB_Info LAGraph_RPQMatrix_solver(RPQMatrixPlan *plan, char *msg)
479501 {
480502 case RPQ_MATRIX_OP_LABEL :
481503 LG_ASSERT_MSG (plan -> lhs == NULL && plan -> rhs == NULL ,
482- GrB_INVALID_VALUE , "label node should not have any children nodes" );
504+ GrB_INVALID_VALUE , "label node should not have any children nodes" );
483505 plan -> res_mat = plan -> mat ;
484506 return (GrB_SUCCESS );
485507 case RPQ_MATRIX_OP_LOR :
0 commit comments