20 #include "common/frobeniusupdate.h"
26 #ifndef DOXYGEN_SHOULD_SKIP_THIS
27 static pastix_complex64_t zone = 1.0;
28 static pastix_complex64_t zzero = 0.0;
131 pastix_complex64_t *A,
133 pastix_complex64_t *tau,
134 pastix_complex64_t *B,
136 pastix_complex64_t *tau_b,
137 pastix_complex64_t *work,
141 int SEED[4] = {26, 67, 52, 197};
149 pastix_complex64_t *omega = work;
150 pastix_complex64_t *subw = work;
153 sublw = pastix_imax( sublw, size_O );
156 #if defined(PRECISION_c) || defined(PRECISION_z)
164 work[0] = (pastix_complex64_t)lwkopt;
174 if (lda < pastix_imax(1, m)) {
177 if( lwork < lwkopt ) {
182 minMN = pastix_imin(m, n);
186 maxrank = pastix_imin( maxrank, minMN );
190 normA = LAPACKE_zlange_work( LAPACK_COL_MAJOR,
'f', m, n,
199 if ( maxrank == 0 ) {
214 #if defined(PASTIX_DEBUG_LR)
215 omega = malloc( size_O *
sizeof(pastix_complex64_t) );
216 subw = malloc( sublw *
sizeof(pastix_complex64_t) );
220 ret = LAPACKE_zlarnv_work(3, SEED, size_O, omega);
224 while ( (rk < maxrank) && loop )
231 ib = pastix_imin( bp, maxrank-rk );
235 cblas_zgemm( CblasColMajor, CblasConjTrans, CblasNoTrans,
237 CBLAS_SADDR(zone), A + rk*lda + rk, lda,
239 CBLAS_SADDR(zzero), B + rk*ldb + rk, ldb );
247 cblas_zgemm( CblasColMajor, CblasNoTrans, CblasNoTrans,
249 CBLAS_SADDR(zone), A + rk*lda + rk, lda,
250 B + rk*ldb + rk, ldb,
251 CBLAS_SADDR(zzero), omega, ldo );
253 cblas_zgemm( CblasColMajor, CblasConjTrans, CblasNoTrans,
255 CBLAS_SADDR(zone), A + rk*lda + rk, lda,
257 CBLAS_SADDR(zzero), B + rk*ldb + rk, ldb );
261 ret = LAPACKE_zlarnv_work(3, SEED, size_O, omega);
269 ret = LAPACKE_zgeqrf_work( LAPACK_COL_MAJOR, n-rk, d,
270 B + rk*ldb + rk, ldb, tau_b+rk,
277 ret = LAPACKE_zunmqr_work( LAPACK_COL_MAJOR,
'R',
'N',
279 B + rk*ldb + rk, ldb, tau_b+rk,
280 A + rk*lda + rk, lda,
287 ret = LAPACKE_zgeqrf_work( LAPACK_COL_MAJOR, m-rk, d,
288 A + rk*lda + rk, lda, tau + rk,
296 ret = LAPACKE_zunmqr_work( LAPACK_COL_MAJOR,
'L', trans,
298 A + rk *lda + rk, lda, tau + rk,
299 A + (rk+d)*lda + rk, lda,
305 normR = LAPACKE_zlange_work( LAPACK_COL_MAJOR,
'f', m-rk-d, n-rk-d, A + (rk+d) * (lda+1), lda, NULL );
313 for( i=d-1; i>=0; i-- ) {
314 double normRk = cblas_dznrm2( n-rk-i, A + (rk+i) * lda + (rk+i), lda );
316 frobenius_update( 1., &scl, &ssq, &normRk );
317 normRk = scl * sqrt( ssq );
319 if ( normRk > tol ) {
334 #if defined(PASTIX_DEBUG_LR)
340 if ( (loop && (rk < minMN)) || (rk > maxrank) ) {
BEGIN_C_DECLS typedef int pastix_int_t
int core_zrqrrt(double tol, pastix_int_t maxrank, pastix_int_t nb, pastix_int_t m, pastix_int_t n, pastix_complex64_t *A, pastix_int_t lda, pastix_complex64_t *tau, pastix_complex64_t *B, pastix_int_t ldb, pastix_complex64_t *tau_b, pastix_complex64_t *work, pastix_int_t lwork, double normA)
Compute a randomized QR factorization with rotation technique.
The block low-rank structure to hold a matrix in low-rank form.
pastix_fixdbl_t core_zge2lr_rqrrt(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 RQRRT.
pastix_fixdbl_t core_zge2lr_qrrt(core_zrrqr_rt_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.
Manage nancheck for lowrank kernels. This header describes all the LAPACKE functions used for low-ran...