PaStiX Handbook  6.2.1
core_clrothu.c
Go to the documentation of this file.
1 /**
2  *
3  * @file core_clrothu.c
4  *
5  * PaStiX low-rank kernel routines to othogonalize the U matrix with QR approximations.
6  *
7  * @copyright 2016-2021 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria,
8  * Univ. Bordeaux. All rights reserved.
9  *
10  * @version 6.0.3
11  * @author Alfredo Buttari
12  * @author Gregoire Pichon
13  * @author Esragul Korkmaz
14  * @author Mathieu Faverge
15  * @date 2019-11-12
16  * @generated from /builds/solverstack/pastix/kernels/core_zlrothu.c, normal z -> c, Tue Apr 12 09:38:40 2022
17  *
18  **/
19 #include "common.h"
20 #include <cblas.h>
21 #include <lapacke.h>
22 #include "flops.h"
23 #include "kernels_trace.h"
24 #include "blend/solver.h"
25 #include "pastix_ccores.h"
26 #include "pastix_clrcores.h"
27 #include "c_nan_check.h"
28 #include "pastix_lowrank.h"
29 
30 #ifndef DOXYGEN_SHOULD_SKIP_THIS
31 static pastix_complex32_t mcone = -1.0;
32 static pastix_complex32_t cone = 1.0;
33 static pastix_complex32_t czero = 0.0;
34 #endif /* DOXYGEN_SHOULD_SKIP_THIS */
35 
36 /**
37  *******************************************************************************
38  *
39  * @brief Try to orthognalize the u part of the low-rank form, and update the v
40  * part accordingly using full QR.
41  *
42  * This function considers a low-rank matrix resulting from the addition of two
43  * matrices B += A, with A of smaller or equal size to B.
44  * The product has the form: U * V^t
45  *
46  * The U part of the low-rank form must be orthognalized to get the smaller
47  * possible rank during the rradd operation. This function perfoms this by
48  * applying a full QR factorization on the U part.
49  *
50  * U = Q R, then U' = Q, and V' = R * V
51  *
52  *******************************************************************************
53  *
54  * @param[in] M
55  * The number of rows of the u1u2 matrix.
56  *
57  * @param[in] N
58  * The number of columns of the v1v2 matrix.
59  *
60  * @param[in] rank
61  * The number of columns of the U matrix, and the number of rows of the
62  * V part in the v1v2 matrix.
63  *
64  * @param[inout] U
65  * The U matrix of size ldu -by- rank. On exit, Q from U = Q R.
66  *
67  * @param[in] ldu
68  * The leading dimension of the U matrix. ldu >= max(1, M)
69  *
70  * @param[inout] V
71  * The V matrix of size ldv -by- N.
72  * On exit, R * V, with R from U = Q R.
73  *
74  * @param[in] ldv
75  * The leading dimension of the V matrix. ldv >= max(1, rank)
76  *
77  *******************************************************************************
78  *
79  * @return The number of flops required to perform the operation.
80  *
81  *******************************************************************************/
82 pastix_fixdbl_t
83 core_clrorthu_fullqr( pastix_int_t M, pastix_int_t N, pastix_int_t rank,
84  pastix_complex32_t *U, pastix_int_t ldu,
85  pastix_complex32_t *V, pastix_int_t ldv )
86 {
87  pastix_int_t minMK = pastix_imin( M, rank );
88  pastix_int_t lwork = M * 32 + minMK;
89  pastix_int_t ret;
90  pastix_complex32_t *W = malloc( lwork * sizeof(pastix_complex32_t) );
91  pastix_complex32_t *tau, *work;
92  pastix_fixdbl_t flops = 0.;
93 
94  tau = W;
95  work = W + minMK;
96  lwork -= minMK;
97 
98  assert( M >= rank );
99 
100  /* Compute U = Q * R */
101  ret = LAPACKE_cgeqrf_work( LAPACK_COL_MAJOR, M, rank,
102  U, ldu, tau, work, lwork );
103  assert( ret == 0 );
104  flops += FLOPS_CGEQRF( M, rank );
105 
106  /* Compute V' = R * V' */
107  cblas_ctrmm( CblasColMajor,
108  CblasLeft, CblasUpper,
109  CblasNoTrans, CblasNonUnit,
110  rank, N, CBLAS_SADDR(cone),
111  U, ldu, V, ldv );
112  flops += FLOPS_CTRMM( PastixLeft, rank, N );
113 
114  /* Generate the Q */
115  ret = LAPACKE_cungqr_work( LAPACK_COL_MAJOR, M, rank, rank,
116  U, ldu, tau, work, lwork );
117  assert( ret == 0 );
118  flops += FLOPS_CUNGQR( M, rank, rank );
119 
120  free(W);
121 
122  (void)ret;
123  return flops;
124 }
125 
126 /**
127  *******************************************************************************
128  *
129  * @brief Try to orthognalize the U part of the low-rank form, and update the V
130  * part accordingly using partial QR.
131  *
132  * This function considers a low-rank matrix resulting from the addition of two
133  * matrices B += A, with A of smaller or equal size to B.
134  * The product has the form: U * V^t
135  *
136  * The U part of the low-rank form must be orthognalized to get the smaller
137  * possible rank during the rradd operation. This function perfoms this by
138  * applying a full QR factorization on the U part.
139  *
140  * In that case, it takes benefit from the fact that U = [ u1, u2 ], and V = [
141  * v1, v2 ] with u2 and v2 wich are matrices of respective size M2-by-r2, and
142  * r2-by-N2, offset by offx and offy
143  *
144  * The steps are:
145  * - Scaling of u2 with removal of the null columns
146  * - Orthogonalization of u2 relatively to u1
147  * - Application of the update to v2
148  * - orthogonalization through QR of u2
149  * - Update of V
150  *
151  *******************************************************************************
152  *
153  * @param[in] M
154  * The number of rows of the u1u2 matrix.
155  *
156  * @param[in] N
157  * The number of columns of the v1v2 matrix.
158  *
159  * @param[in] r1
160  * The number of columns of the U matrix in the u1 part, and the number
161  * of rows of the V part in the v1 part.
162  *
163  * @param[inout] r2ptr
164  * The number of columns of the U matrix in the u2 part, and the number
165  * of rows of the V part in the v2 part. On exit, this rank is reduced
166  * y the number of null columns found in U.
167  *
168  * @param[in] offx
169  * The row offset of the matrix u2 in U.
170  *
171  * @param[in] offy
172  * The column offset of the matrix v2 in V.
173  *
174  * @param[inout] U
175  * The U matrix of size ldu -by- rank. On exit, the orthogonalized U.
176  *
177  * @param[in] ldu
178  * The leading dimension of the U matrix. ldu >= max(1, M)
179  *
180  * @param[inout] V
181  * The V matrix of size ldv -by- N.
182  * On exit, the updated V matrix.
183  *
184  * @param[in] ldv
185  * The leading dimension of the V matrix. ldv >= max(1, rank)
186  *
187  *******************************************************************************
188  *
189  * @return The number of flops required to perform the operation.
190  *
191  *******************************************************************************/
192 pastix_fixdbl_t
193 core_clrorthu_partialqr( pastix_int_t M, pastix_int_t N,
194  pastix_int_t r1, pastix_int_t *r2ptr,
195  pastix_int_t offx, pastix_int_t offy,
196  pastix_complex32_t *U, pastix_int_t ldu,
197  pastix_complex32_t *V, pastix_int_t ldv )
198 {
199  pastix_int_t r2 = *r2ptr;
200  pastix_int_t minMN = pastix_imin( M, r2 );
201  pastix_int_t ldwork = pastix_imax( r1 * r2, M * 32 + minMN );
202  pastix_int_t ret, i;
203  pastix_complex32_t *u1 = U;
204  pastix_complex32_t *u2 = U + r1 * ldu;
205  pastix_complex32_t *v1 = V;
206  pastix_complex32_t *v2 = V + r1;
207  pastix_complex32_t *W = malloc( ldwork * sizeof(pastix_complex32_t) );
208  pastix_complex32_t *tau, *work;
209  pastix_fixdbl_t flops = 0.;
210  float norm, eps;
211 
212  tau = W;
213  work = W + minMN;
214  ldwork -= minMN;
215 
216  eps = LAPACKE_slamch_work('e');
217 
218  /* Scaling */
219  for (i=0; i<r2; i++, u2 += ldu, v2++) {
220  norm = cblas_scnrm2( M, u2, 1 );
221  if ( norm > (M * eps) ) {
222  cblas_csscal( M, 1. / norm, u2, 1 );
223  cblas_csscal( N, norm, v2, ldv );
224  }
225  else {
226  if ( i < (r2-1) ) {
227  cblas_cswap( M, u2, 1, U + (r1+r2-1) * ldu, 1 );
228  memset( U + (r1+r2-1) * ldu, 0, M * sizeof(pastix_complex32_t) );
229 
230  cblas_cswap( N, v2, ldv, V + (r1+r2-1), ldv );
231  LAPACKE_claset_work( LAPACK_COL_MAJOR, 'A', 1, N,
232  0., 0., V + (r1+r2-1), ldv );
233  r2--;
234  i--;
235  u2-= ldu;
236  v2--;
237  }
238  else {
239  memset( u2, 0, M * sizeof(pastix_complex32_t) );
240  LAPACKE_claset_work( LAPACK_COL_MAJOR, 'A', 1, N,
241  0., 0., v2, ldv );
242  r2--;
243  }
244  }
245  }
246  u2 = U + r1 * ldu;
247  v2 = V + r1;
248 
249  *r2ptr = r2;
250 
251  if ( r2 == 0 ) {
252  free( W );
253  return 0.;
254  }
255 
256  /* Compute W = u1^t u2 */
257  cblas_cgemm( CblasColMajor, CblasConjTrans, CblasNoTrans,
258  r1, r2, M,
259  CBLAS_SADDR(cone), u1, ldu,
260  u2, ldu,
261  CBLAS_SADDR(czero), W, r1 );
262  flops += FLOPS_CGEMM( r1, r2, M );
263 
264  /* Compute u2 = u2 - u1 ( u1^t u2 ) = u2 - u1 * W */
265  cblas_cgemm( CblasColMajor, CblasNoTrans, CblasNoTrans,
266  M, r2, r1,
267  CBLAS_SADDR(mcone), u1, ldu,
268  W, r1,
269  CBLAS_SADDR(cone), u2, ldu );
270  flops += FLOPS_CGEMM( M, r2, r1 );
271 
272  /* Update v1 = v1 + ( u1^t u2 ) v2 = v1 + W * v2 */
273  cblas_cgemm( CblasColMajor, CblasNoTrans, CblasNoTrans,
274  r1, N, r2,
275  CBLAS_SADDR(cone), W, r1,
276  v2, ldv,
277  CBLAS_SADDR(cone), v1, ldv );
278  flops += FLOPS_CGEMM( r1, N, r2 );
279 
280 #if !defined(PASTIX_LR_CGS1)
281  /* Compute W = u1^t u2 */
282  cblas_cgemm( CblasColMajor, CblasConjTrans, CblasNoTrans,
283  r1, r2, M,
284  CBLAS_SADDR(cone), u1, ldu,
285  u2, ldu,
286  CBLAS_SADDR(czero), W, r1 );
287  flops += FLOPS_CGEMM( r1, r2, M );
288 
289  /* Compute u2 = u2 - u1 ( u1^t u2 ) = u2 - u1 * W */
290  cblas_cgemm( CblasColMajor, CblasNoTrans, CblasNoTrans,
291  M, r2, r1,
292  CBLAS_SADDR(mcone), u1, ldu,
293  W, r1,
294  CBLAS_SADDR(cone), u2, ldu );
295  flops += FLOPS_CGEMM( M, r2, r1 );
296 
297  /* Update v1 = v1 + ( u1^t u2 ) v2 = v1 + W * v2 */
298  cblas_cgemm( CblasColMajor, CblasNoTrans, CblasNoTrans,
299  r1, N, r2,
300  CBLAS_SADDR(cone), W, r1,
301  v2, ldv,
302  CBLAS_SADDR(cone), v1, ldv );
303  flops += FLOPS_CGEMM( r1, N, r2 );
304 #endif
305 
306 #if defined(PASTIX_DEBUG_LR)
307  if ( core_clrdbg_check_orthogonality_AB( M, r1, r2, u1, ldu, u2, ldu ) != 0 ) {
308  fprintf(stderr, "partialQR: u2 not correctly projected with u1\n" );
309  }
310 #endif
311 
312  /* Compute u2 = Q * R */
313  ret = LAPACKE_cgeqrf_work( LAPACK_COL_MAJOR, M, r2,
314  u2, ldu, tau, work, ldwork );
315  assert( ret == 0 );
316  flops += FLOPS_CGEQRF( M, r2 );
317 
318  /* Compute v2' = R * v2 */
319  cblas_ctrmm( CblasColMajor,
320  CblasLeft, CblasUpper,
321  CblasNoTrans, CblasNonUnit,
322  r2, N, CBLAS_SADDR(cone),
323  u2, ldu, v2, ldv);
324  flops += FLOPS_CTRMM( PastixLeft, r2, N );
325 
326  /* Generate the Q */
327  ret = LAPACKE_cungqr_work( LAPACK_COL_MAJOR, M, r2, r2,
328  u2, ldu, tau, work, ldwork );
329  assert( ret == 0 );
330  flops += FLOPS_CUNGQR( M, r2, r2 );
331 
332 #if defined(PASTIX_DEBUG_LR)
333  if ( core_clrdbg_check_orthogonality_AB( M, r1, r2, u1, ldu, u2, ldu ) != 0 ) {
334  fprintf(stderr, "partialQR: Final u2 not orthogonal to u1\n" );
335  }
336 #endif
337 
338  free( W );
339 
340  (void)ret;
341  (void)offx;
342  (void)offy;
343 
344  return flops;
345 }
346 
347 /**
348  *******************************************************************************
349  *
350  * @brief Try to orthognalize the U part of the low-rank form, and update the V
351  * part accordingly using CGS.
352  *
353  * This function considers a low-rank matrix resulting from the addition of two
354  * matrices B += A, with A of smaller or equal size to B.
355  * The product has the form: U * V^t
356  *
357  * The U part of the low-rank form must be orthognalized to get the smaller
358  * possible rank during the rradd operation. This function perfoms this by
359  * applying a full QR factorization on the U part.
360  *
361  * In that case, it takes benefit from the fact that U = [ u1, u2 ], and V = [
362  * v1, v2 ] with u2 and v2 wich are matrices of respective size M2-by-r2, and
363  * r2-by-N2, offset by offx and offy
364  *
365  * The steps are:
366  * - for each column of u2
367  * - Scaling of u2 with removal of the null columns
368  * - Orthogonalization of u2 relatively to u1
369  * - Remove the column if null
370  *
371  *******************************************************************************
372  *
373  * @param[in] M1
374  * The number of rows of the U matrix.
375  *
376  * @param[in] N1
377  * The number of columns of the U matrix.
378  *
379  * @param[in] M2
380  * The number of rows of the u2 part of the U matrix.
381  *
382  * @param[in] N2
383  * The number of columns of the v2 part of the V matrix.
384  *
385  * @param[in] r1
386  * The number of columns of the U matrix in the u1 part, and the number
387  * of rows of the V part in the v1 part.
388  *
389  * @param[inout] r2ptr
390  * The number of columns of the U matrix in the u2 part, and the number
391  * of rows of the V part in the v2 part. On exit, this rank is reduced
392  * y the number of null columns found in U.
393  *
394  * @param[in] offx
395  * The row offset of the matrix u2 in U.
396  *
397  * @param[in] offy
398  * The column offset of the matrix v2 in V.
399  *
400  * @param[inout] U
401  * The U matrix of size ldu -by- rank. On exit, the orthogonalized U.
402  *
403  * @param[in] ldu
404  * The leading dimension of the U matrix. ldu >= max(1, M)
405  *
406  * @param[inout] V
407  * The V matrix of size ldv -by- N.
408  * On exit, the updated V matrix.
409  *
410  * @param[in] ldv
411  * The leading dimension of the V matrix. ldv >= max(1, rank)
412  *
413  *******************************************************************************
414  *
415  * @return The number of flops required to perform the operation.
416  *
417  *******************************************************************************/
418 pastix_fixdbl_t
419 core_clrorthu_cgs( pastix_int_t M1, pastix_int_t N1,
420  pastix_int_t M2, pastix_int_t N2,
421  pastix_int_t r1, pastix_int_t *r2ptr,
422  pastix_int_t offx, pastix_int_t offy,
423  pastix_complex32_t *U, pastix_int_t ldu,
424  pastix_complex32_t *V, pastix_int_t ldv )
425 {
426  pastix_int_t r2 = *r2ptr;
427  pastix_complex32_t *u1 = U;
428  pastix_complex32_t *u2 = U + r1 * ldu;
429  pastix_complex32_t *v1 = V;
430  pastix_complex32_t *v2 = V + r1;
431  pastix_complex32_t *W;
432  pastix_fixdbl_t flops = 0.0;
433  pastix_int_t i, rank = r1 + r2;
434  pastix_int_t ldwork = rank;
435  float eps, norm;
436  float norm_before, alpha;
437 
438  assert( M1 >= (M2 + offx) );
439  assert( N1 >= (N2 + offy) );
440 
441  W = malloc(ldwork * sizeof(pastix_complex32_t));
442  eps = LAPACKE_slamch_work( 'e' );
443  alpha = 1. / sqrtf(2);
444 
445  /* Classical Gram-Schmidt */
446  for (i=r1; i<rank; i++, u2 += ldu, v2++) {
447 
448  norm = cblas_scnrm2( M2, u2 + offx, 1 );
449  if ( norm > ( M2 * eps ) ) {
450  cblas_csscal( M2, 1. / norm, u2 + offx, 1 );
451  cblas_csscal( N2, norm, v2 + offy * ldv, ldv );
452  }
453  else {
454  rank--; r2--;
455  if ( i < rank ) {
456  cblas_cswap( M2, u2 + offx, 1, U + rank * ldu + offx, 1 );
457 #if !defined(NDEBUG)
458  memset( U + rank * ldu, 0, M1 * sizeof(pastix_complex32_t) );
459 #endif
460 
461  cblas_cswap( N2, v2 + offy * ldv, ldv, V + offy * ldv + rank, ldv );
462 
463 #if !defined(NDEBUG)
464  LAPACKE_claset_work( LAPACK_COL_MAJOR, 'A', 1, N1,
465  0., 0., V + rank, ldv );
466 #endif
467  i--;
468  u2-= ldu;
469  v2--;
470  }
471 #if !defined(NDEBUG)
472  else {
473  memset( u2, 0, M1 * sizeof(pastix_complex32_t) );
474  LAPACKE_claset_work( LAPACK_COL_MAJOR, 'A', 1, N1,
475  0., 0., v2, ldv );
476  }
477 #endif
478  continue;
479  }
480 
481  /* Compute W = u1^t u2 */
482  cblas_cgemv( CblasColMajor, CblasConjTrans,
483  M2, i,
484  CBLAS_SADDR(cone), u1+offx, ldu,
485  u2+offx, 1,
486  CBLAS_SADDR(czero), W, 1 );
487  flops += FLOPS_CGEMM( M2, i, 1 );
488 
489  /* Compute u2 = u2 - u1 ( u1^t u2 ) = u2 - u1 * W */
490  cblas_cgemv( CblasColMajor, CblasNoTrans,
491  M1, i,
492  CBLAS_SADDR(mcone), u1, ldu,
493  W, 1,
494  CBLAS_SADDR(cone), u2, 1 );
495  flops += FLOPS_CGEMM( M1, i, 1 );
496 
497  /* Update v1 = v1 + ( u1^t u2 ) v2 = v1 + W * v2 */
498  cblas_cgemm( CblasColMajor, CblasNoTrans, CblasNoTrans,
499  i, N1, 1,
500  CBLAS_SADDR(cone), W, i,
501  v2, ldv,
502  CBLAS_SADDR(cone), v1, ldv );
503  flops += FLOPS_CGEMM( i, N1, 1 );
504 
505  norm_before = cblas_scnrm2( i, W, 1 );
506  norm = cblas_scnrm2( M1, u2, 1 );
507 
508 #if !defined(PASTIX_LR_CGS1)
509  if ( norm <= (alpha * norm_before) ){
510  /* Compute W = u1^t u2 */
511  cblas_cgemv( CblasColMajor, CblasConjTrans,
512  M1, i,
513  CBLAS_SADDR(cone), u1, ldu,
514  u2, 1,
515  CBLAS_SADDR(czero), W, 1 );
516  flops += FLOPS_CGEMM( M1, i, 1 );
517 
518  /* Compute u2 = u2 - u1 ( u1^t u2 ) = u2 - u1 * W */
519  cblas_cgemv( CblasColMajor, CblasNoTrans,
520  M1, i,
521  CBLAS_SADDR(mcone), u1, ldu,
522  W, 1,
523  CBLAS_SADDR(cone), u2, 1 );
524  flops += FLOPS_CGEMM( M1, i, 1 );
525 
526  /* Update v1 = v1 + ( u1^t u2 ) v2 = v1 + W * v2 */
527  cblas_cgemm( CblasColMajor, CblasNoTrans, CblasNoTrans,
528  i, N1, 1,
529  CBLAS_SADDR(cone), W, i,
530  v2, ldv,
531  CBLAS_SADDR(cone), v1, ldv );
532  flops += FLOPS_CGEMM( i, N1, 1 );
533 
534  norm = cblas_scnrm2( M1, u2, 1 );
535  }
536 #endif
537 
538  if ( norm > M1 * eps ) {
539  cblas_csscal( M1, 1. / norm, u2, 1 );
540  cblas_csscal( N1, norm, v2, ldv );
541  }
542  else {
543  rank--; r2--;
544  if ( i < rank ) {
545  cblas_cswap( M1, u2, 1, U + rank * ldu, 1 );
546  memset( U + rank * ldu, 0, M1 * sizeof(pastix_complex32_t) );
547 
548  cblas_cswap( N1, v2, ldv, V + rank, ldv );
549  LAPACKE_claset_work( LAPACK_COL_MAJOR, 'A', 1, N1,
550  0., 0., V + rank, ldv );
551  i--;
552  u2-= ldu;
553  v2--;
554  }
555  else {
556  memset( u2, 0, M1 * sizeof(pastix_complex32_t) );
557  LAPACKE_claset_work( LAPACK_COL_MAJOR, 'A', 1, N1,
558  0., 0., v2, ldv );
559  }
560  }
561  }
562  free(W);
563 
564 #if defined(PASTIX_DEBUG_LR)
565  {
566  u2 = U + r1 * ldu;
567  if ( core_clrdbg_check_orthogonality_AB( M1, r1, r2, u1, ldu, u2, ldu ) != 0 ) {
568  fprintf(stderr, "cgs: Final u2 not orthogonal to u1\n" );
569  }
570  }
571 #endif
572 
573  *r2ptr = r2;
574 
575  (void)offy;
576  (void)N2;
577  return flops;
578 }
solver.h
pastix_lowrank.h
core_clrorthu_partialqr
pastix_fixdbl_t core_clrorthu_partialqr(pastix_int_t M, pastix_int_t N, pastix_int_t r1, pastix_int_t *r2ptr, pastix_int_t offx, pastix_int_t offy, pastix_complex32_t *U, pastix_int_t ldu, pastix_complex32_t *V, pastix_int_t ldv)
Try to orthognalize the U part of the low-rank form, and update the V part accordingly using partial ...
Definition: core_clrothu.c:193
c_nan_check.h
Manage nancheck for lowrank kernels. This header describes all the LAPACKE functions used for low-ran...
core_clrorthu_fullqr
pastix_fixdbl_t core_clrorthu_fullqr(pastix_int_t M, pastix_int_t N, pastix_int_t rank, pastix_complex32_t *U, pastix_int_t ldu, pastix_complex32_t *V, pastix_int_t ldv)
Try to orthognalize the u part of the low-rank form, and update the v part accordingly using full QR.
Definition: core_clrothu.c:83
pastix_ccores.h
core_clrorthu_cgs
pastix_fixdbl_t core_clrorthu_cgs(pastix_int_t M1, pastix_int_t N1, pastix_int_t M2, pastix_int_t N2, pastix_int_t r1, pastix_int_t *r2ptr, pastix_int_t offx, pastix_int_t offy, pastix_complex32_t *U, pastix_int_t ldu, pastix_complex32_t *V, pastix_int_t ldv)
Try to orthognalize the U part of the low-rank form, and update the V part accordingly using CGS.
Definition: core_clrothu.c:419
core_clrdbg_check_orthogonality_AB
int core_clrdbg_check_orthogonality_AB(pastix_int_t M, pastix_int_t NA, pastix_int_t NB, const pastix_complex32_t *A, pastix_int_t lda, const pastix_complex32_t *B, pastix_int_t ldb)
Check the orthogonality of the matrix A relatively to the matrix B.
Definition: core_clrdbg.c:181
pastix_clrcores.h
PastixLeft
@ PastixLeft
Definition: api.h:473