27 #ifndef DOXYGEN_SHOULD_SKIP_THIS
28 static pastix_complex64_t mzone = -1.0;
29 static pastix_complex64_t zone = 1.0;
30 static pastix_complex64_t zzero = 0.0;
120 pastix_complex64_t *A,
123 pastix_complex64_t *tau,
124 pastix_complex64_t *work,
128 int SEED[4] = {26, 67, 52, 197};
140 double tolB = sqrt( (
double)(bp) ) * tol;
142 pastix_complex64_t *B = work;
143 pastix_complex64_t *tau_b = B + size_B;
144 pastix_complex64_t *omega = tau_b + n;
145 pastix_complex64_t *subw = tau_b + n;
147 sublw = pastix_imax( sublw, size_O );
150 #if defined(PRECISION_c) || defined(PRECISION_z)
156 lwkopt = size_B + n + sublw;
158 work[0] = (pastix_complex64_t)lwkopt;
168 if (lda < pastix_imax(1, m)) {
171 if( lwork < lwkopt ) {
176 minMN = pastix_imin(m, n);
180 maxrank = pastix_imin( maxrank, minMN );
187 if ( maxrank == 0 ) {
192 norm = LAPACKE_zlange_work( LAPACK_COL_MAJOR,
'f', m, n,
200 #if defined(PASTIX_DEBUG_LR)
201 B = malloc( size_B *
sizeof(pastix_complex64_t) );
202 tau_b = malloc( n *
sizeof(pastix_complex64_t) );
203 omega = malloc( size_O *
sizeof(pastix_complex64_t) );
204 subw = malloc( sublw *
sizeof(pastix_complex64_t) );
208 for (j=0; j<n; j++) jpvt[j] = j;
211 ret = LAPACKE_zlarnv_work(3, SEED, size_O, omega);
215 cblas_zgemm( CblasColMajor, CblasNoTrans, CblasNoTrans,
217 CBLAS_SADDR(zone), omega, ldo,
219 CBLAS_SADDR(zzero), B, ldb );
224 ib = pastix_imin( b, minMN-rk );
245 if ( (rk + d) > maxrank ) {
251 for (j = rk; j < rk + d; j++) {
252 if (jpvt_b[j] >= 0) {
257 jpvt_b[k] = - jpvt_b[k] - 1;
259 while( jpvt_b[in] >= 0 ) {
262 cblas_zswap( m, A + k * lda, 1,
270 jpvt_b[in] = - jpvt_b[in] - 1;
280 ret = LAPACKE_zgeqrf_work( LAPACK_COL_MAJOR, m-rk, d,
281 A + rk*lda + rk, lda, tau + rk,
290 ret = LAPACKE_zunmqr_work( LAPACK_COL_MAJOR,
'L', trans,
292 A + rk *lda + rk, lda, tau + rk,
293 A + (rk+d)*lda + rk, lda,
301 ret = LAPACKE_zlaset_work( LAPACK_COL_MAJOR,
'L', d-1, d-1,
302 0, 0, B + rk*ldb + 1, ldb );
308 cblas_ztrsm( CblasColMajor, CblasRight, CblasUpper,
309 CblasNoTrans, CblasNonUnit,
311 CBLAS_SADDR(zone), A + rk*lda + rk, lda,
314 cblas_zgemm( CblasColMajor, CblasNoTrans, CblasNoTrans,
316 CBLAS_SADDR(mzone), B + rk *ldb, ldb,
317 A + (rk+d)*lda + rk, lda,
318 CBLAS_SADDR(zone), B + (rk+d)*ldb, ldb );
325 if ( refine && !loop && (rk < maxrank) ) {
332 A + rk * lda + rk, lda,
334 work, lwork, rwork );
337 for (j=0; j<d; j++) {
338 if (jpvt_b[j] >= 0) {
343 jpvt_b[k] = - jpvt_b[k] - 1;
345 while( jpvt_b[in] >= 0 ) {
348 cblas_zswap( rk, A + (rk + k ) * lda, 1,
349 A + (rk + in) * lda, 1 );
352 jpvt[rk + k] = jpvt[rk + in];
353 jpvt[rk + in] = itmp;
356 jpvt_b[in] = - jpvt_b[in] - 1;
364 if ( rk > maxrank ) {
371 #if defined(PASTIX_DEBUG_LR)
493 const void *alphaptr,
504 M1, N1, A, M2, N2, B, offx, offy );
BEGIN_C_DECLS typedef int pastix_int_t
int core_zrqrcp(double tol, pastix_int_t maxrank, int refine, pastix_int_t nb, pastix_int_t m, pastix_int_t n, pastix_complex64_t *A, pastix_int_t lda, pastix_int_t *jpvt, pastix_complex64_t *tau, pastix_complex64_t *work, pastix_int_t lwork, double *rwork)
Compute a randomized QR factorization.
int core_zpqrcp(double tol, pastix_int_t maxrank, int full_update, pastix_int_t nb, pastix_int_t m, pastix_int_t n, pastix_complex64_t *A, pastix_int_t lda, pastix_int_t *jpvt, pastix_complex64_t *tau, pastix_complex64_t *work, pastix_int_t lwork, double *rwork)
Compute a rank-reavealing QR factorization.
Structure to define the type of function to use for the low-rank kernels and their parameters.
The block low-rank structure to hold a matrix in low-rank form.
pastix_fixdbl_t core_zge2lr_rqrcp(int use_reltol, pastix_fixdbl_t tol, pastix_int_t rklimit, pastix_int_t m, pastix_int_t n, const void *A, pastix_int_t lda, pastix_lrblock_t *Alr)
Convert a full rank matrix in a low rank matrix, using RQRCP.
pastix_fixdbl_t core_zge2lr_qrcp(core_zrrqr_cp_t rrqrfct, int use_reltol, pastix_fixdbl_t tol, pastix_int_t rklimit, pastix_int_t m, pastix_int_t n, const void *Avoid, pastix_int_t lda, pastix_lrblock_t *Alr)
Template to convert a full rank matrix into a low rank matrix through QR decompositions.
pastix_fixdbl_t core_zrradd_qr(core_zrrqr_cp_t rrqrfct, const pastix_lr_t *lowrank, pastix_trans_t transA1, const void *alphaptr, pastix_int_t M1, pastix_int_t N1, const pastix_lrblock_t *A, pastix_int_t M2, pastix_int_t N2, pastix_lrblock_t *B, pastix_int_t offx, pastix_int_t offy)
Template to perform the addition of two low-rank structures with compression kernel based on QR decom...
pastix_fixdbl_t core_zrradd_rqrcp(const pastix_lr_t *lowrank, pastix_trans_t transA1, const void *alphaptr, pastix_int_t M1, pastix_int_t N1, const pastix_lrblock_t *A, pastix_int_t M2, pastix_int_t N2, pastix_lrblock_t *B, pastix_int_t offx, pastix_int_t offy)
Add two LR structures A=(-u1) v1^T and B=u2 v2^T into u2 v2^T.
enum pastix_trans_e pastix_trans_t
Transpostion.
Manage nancheck for lowrank kernels. This header describes all the LAPACKE functions used for low-ran...