Any BV type will do. The default BVSVEC is generally best. Jose
> El 30 oct 2017, a las 17:18, Franck Houssen <[email protected]> > escribió: > > It was not clear to me when I read the doc. That's OK now: got it to work, > thanks Jose ! > Just to make sure, to make it work, I had to set a BV type: I chose BVMAT as > I use BVCreateFromMat. Is that the good type ? (BVVECS works too) > > Franck > > ----- Mail original ----- >> De: "Jose E. Roman" <[email protected]> >> À: "Franck Houssen" <[email protected]> >> Cc: "For users of the development version of PETSc" <[email protected]> >> Envoyé: Samedi 28 Octobre 2017 16:56:22 >> Objet: Re: [petsc-dev] QR factorization of dense matrix >> >> Matrix R must be mxm. >> BVOrthogonalize computes Z=Q*R, where Q overwrites Z. >> Jose >> >>> El 28 oct 2017, a las 13:11, Franck Houssen <[email protected]> >>> escribió: >>> >>> I've seen that !... But can't get BVOrthogonalize to work. >>> >>> I tried: >>> Mat Z; MatCreateSeqDense(PETSC_COMM_SELF, n, m, NULL, &Z); >>> ...; // MatSetValues(Z, ...) >>> BVCreate(PETSC_COMM_SELF, &bv); >>> BVCreateFromMat(Z, &bv); // Z is tall-skinny >>> Mat R; MatCreateSeqDense(PETSC_COMM_SELF, n, m, NULL, &R); // Same n, m >>> than Z. >>> BVOrthogonalize(bv, R); >>> >>> But BVOrthogonalize fails with : >>>> [0]PETSC ERROR: Nonconforming object sizes >>>> [0]PETSC ERROR: Mat argument is not square, it has 1 rows and 3 columns >>> >>> So, as I didn't get what's wrong, I was looking for another way to do this. >>> >>> Franck >>> >>> ----- Mail original ----- >>>> De: "Jose E. Roman" <[email protected]> >>>> À: "Franck Houssen" <[email protected]> >>>> Cc: "For users of the development version of PETSc" >>>> <[email protected]> >>>> Envoyé: Vendredi 27 Octobre 2017 19:03:37 >>>> Objet: Re: [petsc-dev] QR factorization of dense matrix >>>> >>>> Franck, >>>> >>>> SLEPc has some support for this, but it is intended only for tall-skinny >>>> matrices, that is, when the number of columns is much smaller than rows. >>>> For >>>> an almost square matrix you should not use it. >>>> >>>> Have a look at this >>>> http://slepc.upv.es/documentation/current/docs/manualpages/BV/BVOrthogonalize.html >>>> http://slepc.upv.es/documentation/current/docs/manualpages/BV/BVOrthogBlockType.html >>>> >>>> You can see there are three methods. All of them have drawbacks: >>>> GS: This is a Gram-Schmidt QR, computed column by column, so it is slower >>>> than the other two. However, it is robust. >>>> CHOL: Cholesky QR, it is not numerically stable. In the future we will add >>>> Cholesky QR2. >>>> TSQR: Unfortunately this is not implemented in parallel. I wanted to add >>>> the >>>> parallel version for 3.8, but didn't have time. It will be added soon. >>>> >>>> You can use BVCreateFromMat() to create a BV object from a Mat. >>>> >>>> Jose >>>> >>>> >>>>> El 27 oct 2017, a las 18:39, Franck Houssen <[email protected]> >>>>> escribió: >>>>> >>>>> I am looking for QR factorization of (sequential) dense matrix: is this >>>>> available in PETSc ? I "just" need the diagonal of R (I do not need >>>>> neither the full content of R, nor Q) >>>>> >>>>> I found that (old !) thread >>>>> https://lists.mcs.anl.gov/pipermail/petsc-users/2013-November/019577.html >>>>> that says it could be implemented: has it been done ? >>>>> As for a direct solve, the way to go is "KSPSetType(ksp, KSPPREONLY); >>>>> PCSetType(pc, PCLU);", I was expecting something like "KSPSetType(ksp, >>>>> KSPPREONLY); PCSetType(pc, PCQR);"... But it seems there is no PCQR >>>>> available. Or is it possible to do that using "an iterative way" with a >>>>> specific kind of KSP that triggers a Gram Schmidt orthogonalization in >>>>> back-end ? (I have seen a KSPLSQR but could I get Q and R back ? As I >>>>> understand this, I would say no: I would say the user can only get the >>>>> solution) >>>>> >>>>> Is it possible to QR a (sequential) dense matrix in PETSc ? If yes, what >>>>> are the steps to follow ? >>>>> >>>>> Franck >>>>> >>>>> My understanding is that DGEQRF from lapack can do "more" than what I >>>>> need, >>>>> but, no sure to get if I can use it from PETSc through a KSP: >>>>>>> git grep DGEQRF >>>>> include/petscblaslapack_stdcall.h:# define LAPACKgeqrf_ DGEQRF >>>>>>> git grep LAPACKgeqrf_ >>>>> include/petscblaslapack.h:PETSC_EXTERN void >>>>> LAPACKgeqrf_(PetscBLASInt*,PetscBLASInt*,PetscScalar*,PetscBLASInt*,PetscScalar*,PetscScalar*,PetscBLASInt*,PetscBLASInt*); >>>>> include/petscblaslapack_mangle.h:#define LAPACKgeqrf_ >>>>> PETSCBLAS(geqrf,GEQRF) >>>>> include/petscblaslapack_stdcall.h:# define LAPACKgeqrf_ SGEQRF >>>>> include/petscblaslapack_stdcall.h:# define LAPACKgeqrf_ DGEQRF >>>>> include/petscblaslapack_stdcall.h:# define LAPACKgeqrf_ CGEQRF >>>>> include/petscblaslapack_stdcall.h:# define LAPACKgeqrf_ ZGEQRF >>>>> include/petscblaslapack_stdcall.h:PETSC_EXTERN void PETSC_STDCALL >>>>> LAPACKgeqrf_(PetscBLASInt*,PetscBLASInt*,PetscScalar*,PetscBLASInt*,PetscScalar*,PetscScalar*,PetscBLASInt*,PetscBLASInt*); >>>>> src/dm/dt/interface/dt.c: >>>>> PetscStackCallBLAS("LAPACKgeqrf",LAPACKgeqrf_(&M,&N,A,&lda,tau,work,&ldwork,&info)); >>>>> src/dm/dt/interface/dtfv.c: >>>>> LAPACKgeqrf_(&M,&N,A,&lda,tau,work,&ldwork,&info); >>>>> src/ksp/ksp/impls/gmres/agmres/agmres.c: >>>>> PetscStackCallBLAS("LAPACKgeqrf",LAPACKgeqrf_(&lC, &KspSize, >>>>> agmres->hh_origin, &ldH, agmres->tau, agmres->work, &lwork, &info)); >>>>> src/ksp/pc/impls/bddc/bddcprivate.c: >>>>> PetscStackCallBLAS("LAPACKgeqrf",LAPACKgeqrf_(&Blas_M,&Blas_N,qr_basis,&Blas_LDA,qr_tau,&lqr_work_t,&lqr_work,&lierr)); >>>>> src/ksp/pc/impls/bddc/bddcprivate.c: >>>>> PetscStackCallBLAS("LAPACKgeqrf",LAPACKgeqrf_(&Blas_M,&Blas_N,qr_basis,&Blas_LDA,qr_tau,qr_work,&lqr_work,&lierr)); >>>>> src/ksp/pc/impls/gamg/agg.c: >>>>> PetscStackCallBLAS("LAPACKgeqrf",LAPACKgeqrf_(&Mdata, &N, qqc, &LDA, >>>>> TAU, WORK, &LWORK, &INFO)); >>>>> src/tao/leastsquares/impls/pounders/pounders.c: >>>>> PetscStackCallBLAS("LAPACKgeqrf",LAPACKgeqrf_(&blasnp,&blasnplus1,mfqP->Q_tmp,&blasnpmax,mfqP->tau_tmp,mfqP->mwork,&blasmaxmn,&info)); >>>>> src/tao/leastsquares/impls/pounders/pounders.c: >>>>> PetscStackCallBLAS("LAPACKgeqrf",LAPACKgeqrf_(&blasn,&blask,mfqP->Q,&blasnpmax,mfqP->tau,mfqP->mwork,&blasmaxmn,&info)); >>>> >>>> >> >>
