I’m sorry, solving overdetermined systems, alongside (overlapping) domain decomposition preconditioners and solving systems with multiple right-hand sides, is one of the topic for which I need to stop pushing new features and update the users manual instead… The very basic documentation of PCQR is here: https://petsc.org/main/docs/manualpages/PC/PCQR <https://petsc.org/main/docs/manualpages/PC/PCQR> (I’m guessing you are using the release documentation in which it’s indeed not present). Some comments about your problem of solving Ax=b with a rectangular matrix A. 1) if you switch to KSPLSQR, it is wrong to use KSPSetOperators(ksp, A, A). You can get away with murder if your PCType is PCNONE, but otherwise, you should always supply the normal equations as the Pmat (you will get runtime errors otherwise). To get the normal equations, you can use https://petsc.org/main/docs/manualpages/Mat/MatCreateNormal/ <https://petsc.org/main/docs/manualpages/Mat/MatCreateNormal/> The following two points only applies if your Pmat is sparse (or sparse with some dense rows). 2) there are a couple of PC that handle MATNORMAL: PCNONE, PCQR, PCJACOBI, PCBJACOBI, PCASM, PCHPDDM Currently, PCQR needs SuiteSparse, and thus runs only if the Pmat is distributed on a single MPI process (if your Pmat is distributed on multiple processes, you should first use PCREDUNDANT). 3) if you intend to make your problem scale in sizes, most of these preconditioners will not be very robust. In that case, if your problem does not have any dense rows, you should either use PCHPDDM or MatConvert(Pmat, MATAIJ, PmatAIJ) and then use PCCHOLESKY, PCHYPRE or PCGAMG (you can have a look at https://epubs.siam.org/doi/abs/10.1137/21M1434891 <https://epubs.siam.org/doi/abs/10.1137/21M1434891> for a comparison) If your problem has dense rows, I have somewhere the code to recast it into an augmented system then solved using PCFIELDSPLIT (see https://link.springer.com/article/10.1007/s11075-018-0478-2 <https://link.springer.com/article/10.1007/s11075-018-0478-2>). I can send it to you if needed. Let me know if I can be of further assistance or if something is not clear to you.
Thanks, Pierre > On 26 Sep 2022, at 10:56 AM, fujisan <[email protected]> wrote: > > OK thank you. > > On Mon, Sep 26, 2022 at 10:53 AM Jose E. Roman <[email protected] > <mailto:[email protected]>> wrote: > The QR factorization from SuiteSparse is sequential only, cannot be used in > parallel. > In parallel you can try PCBJACOBI with a PCQR local preconditioner. > Pierre may have additional suggestions. > > Jose > > > > El 26 sept 2022, a las 10:47, fujisan <[email protected] > > <mailto:[email protected]>> escribió: > > > > I did configure Petsc with the option --download-suitesparse. > > > > The error is more like this: > > PETSC ERROR: Could not locate a solver type for factorization type QR and > > matrix type mpiaij. > > > > Fuji > > > > On Mon, Sep 26, 2022 at 10:25 AM Jose E. Roman <[email protected] > > <mailto:[email protected]>> wrote: > > If the error message is "Could not locate a solver type for factorization > > type QR" then you should configure PETSc with --download-suitesparse > > > > Jose > > > > > > > El 26 sept 2022, a las 9:06, fujisan <[email protected] > > > <mailto:[email protected]>> escribió: > > > > > > Thank you Pierre, > > > > > > I used PCNONE along with KSPLSQR and it worked. > > > But as for PCQR, it cannot be found. There is nothing about it in the > > > documentation as well. > > > > > > Fuji > > > > > > On Wed, Sep 21, 2022 at 12:20 PM Pierre Jolivet <[email protected] > > > <mailto:[email protected]>> wrote: > > > Yes, but you need to use a KSP that handles rectangular Mat, such as > > > KSPLSQR (-ksp_type lsqr). > > > PCLU does not handle rectangular Pmat. The only PC that handle > > > rectangular Pmat are PCQR, PCNONE. > > > If you supply the normal equations as the Pmat for LSQR, then you can use > > > “standard” PC. > > > You can have a look at > > > https://petsc.org/main/src/ksp/ksp/tutorials/ex27.c.html > > > <https://petsc.org/main/src/ksp/ksp/tutorials/ex27.c.html> that covers > > > most of these cases. > > > > > > Thanks, > > > Pierre > > > > > > (sorry for the earlier answer sent wrongfully to petsc-maint, please > > > discard the previous email) > > > > > >> On 21 Sep 2022, at 10:03 AM, fujisan <[email protected] > > >> <mailto:[email protected]>> wrote: > > >> > > >> I'm trying to solve Ax=b with a sparse rectangular matrix A (of size > > >> 33x17 in my test) using > > >> options '-ksp_type stcg -pc_type lu' on 1 or 2 cpus. > > >> > > >> And I always get an error saying "Incompatible vector local lengths" > > >> (see below). > > >> > > >> Here is the relevant lines of my code: > > >> > > >> program test > > >> ... > > >> ! Variable declarations > > >> > > >> PetscCallA(PetscInitialize(PETSC_NULL_CHARACTER,ierr)) > > >> > > >> PetscCall(MatCreate(PETSC_COMM_WORLD,A,ierr)) > > >> PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,m,n,ierr)) > > >> PetscCall(MatSetType(A,MATMPIAIJ,ierr)) > > >> PetscCall(MatSetFromOptions(A,ierr)) > > >> PetscCall(MatSetUp(A,ierr)) > > >> PetscCall(MatGetOwnershipRange(A,istart,iend,ierr)) > > >> > > >> do irow=istart,iend-1 > > >> ... Reading from file ... > > >> PetscCall(MatSetValues(A,1,irow,nzv,col,val,ADD_VALUES,ierr)) > > >> ... > > >> enddo > > >> > > >> PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)) > > >> PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)) > > >> > > >> ! Creating vectors x and b > > >> PetscCallA(MatCreateVecs(A,x,b,ierr)) > > >> > > >> ! Duplicating x in u. > > >> PetscCallA(VecDuplicate(x,u,ierr)) > > >> > > >> ! u is used to calculate b > > >> PetscCallA(VecSet(u,1.0,ierr)) > > >> > > >> PetscCallA(VecAssemblyBegin(u,ierr)) > > >> PetscCallA(VecAssemblyEnd(u,ierr)) > > >> > > >> ! Calculating Au = b > > >> PetscCallA(MatMult(A,u,b,ierr)) ! A.u = b > > >> > > >> PetscCallA(KSPSetType(ksp,KSPCG,ierr)) > > >> > > >> PetscCallA(KSPSetOperators(ksp,A,A,ierr)) > > >> > > >> PetscCallA(KSPSetFromOptions(ksp,ierr)) > > >> > > >> ! Solving Ax = b, x unknown > > >> PetscCallA(KSPSolve(ksp,b,x,ierr)) > > >> > > >> PetscCallA(VecDestroy(x,ierr)) > > >> PetscCallA(VecDestroy(u,ierr)) > > >> PetscCallA(VecDestroy(b,ierr)) > > >> PetscCallA(MatDestroy(A,ierr)) > > >> PetscCallA(KSPDestroy(ksp,ierr)) > > >> > > >> call PetscFinalize(ierr) > > >> end program > > >> > > >> The code reads a sparse matrix from a binary file. > > >> I also output the sizes of matrix A and vectors b, x, u. > > >> They all seem consistent. > > >> > > >> What am I doing wrong? > > >> Is it possible to solve Ax=b with A rectangular? > > >> > > >> Thank you in advance for your help. > > >> Have a nice day. > > >> > > >> Fuji > > >> > > >> Matrix size : m= 33 n= 17 cpu size: 1 > > >> Size of matrix A : 33 17 > > >> Size of vector b : 33 > > >> Size of vector x : 17 > > >> Size of vector u : 17 > > >> [0]PETSC ERROR: --------------------- Error Message > > >> -------------------------------------------------------------- > > >> [0]PETSC ERROR: Arguments are incompatible > > >> [0]PETSC ERROR: Incompatible vector local lengths parameter # 1 local > > >> size 33 != parameter # 2 local size 17 > > >> [0]PETSC ERROR: See https://petsc.org/release/faq/ > > >> <https://petsc.org/release/faq/> for trouble shooting. > > >> [0]PETSC ERROR: Petsc Development GIT revision: v3.17.4-1341-g91b2b62a00 > > >> GIT Date: 2022-09-15 19:26:07 +0000 > > >> [0]PETSC ERROR: ./bin/solve on a x86_64 named master by fujisan Tue Sep > > >> 20 16:56:37 2022 > > >> [0]PETSC ERROR: Configure options --with-petsc-arch=x86_64 > > >> --COPTFLAGS="-g -O3" --FOPTFLAGS="-g -O3" --CXXOPTFLAGS="-g -O3" > > >> --with-debugging=0 --with-cc=mpiicc --with-cxx=mpiicpc > > >> --with-fc=mpiifort --with-single-library=1 --with-mpiexec=mpiexec > > >> --with-precision=double --with-fortran-interfaces=1 --with-make=1 > > >> --with-mpi=1 --with-mpi-compilers=1 --download-fblaslapack=0 > > >> --download-hypre=1 --download-cmake=0 --with-cmake=1 --download-metis=1 > > >> --download-parmetis=1 --download-ptscotch=0 --download-suitesparse=1 > > >> --download-triangle=1 --download-superlu=1 --download-superlu_dist=1 > > >> --download-scalapack=1 --download-mumps=1 --download-elemental=1 > > >> --download-spai=0 --download-parms=1 --download-moab=1 > > >> --download-chaco=0 --download-fftw=1 --with-petsc4py=1 > > >> --download-mpi4py=1 --download-saws --download-concurrencykit=1 > > >> --download-revolve=1 --download-cams=1 --download-p4est=0 --with-zlib=1 > > >> --download-mfem=1 --download-glvis=0 --with-opengl=0 --download-libpng=1 > > >> --download-libjpeg=1 --download-slepc=1 --download-hpddm=1 > > >> --download-bamg=1 --download-mmg=0 --download-parmmg=0 > > >> --download-htool=1 --download-egads=0 --download-opencascade=0 > > >> PETSC_ARCH=x86_64 > > >> [0]PETSC ERROR: #1 VecCopy() at > > >> /data/softs/petsc/src/vec/vec/interface/vector.c:1607 > > >> [0]PETSC ERROR: #2 KSPSolve_BiCG() at > > >> /data/softs/petsc/src/ksp/ksp/impls/bicg/bicg.c:40 > > >> [0]PETSC ERROR: #3 KSPSolve_Private() at > > >> /data/softs/petsc/src/ksp/ksp/interface/itfunc.c:877 > > >> [0]PETSC ERROR: #4 KSPSolve() at > > >> /data/softs/petsc/src/ksp/ksp/interface/itfunc.c:1048 > > >> [0]PETSC ERROR: #5 solve.F90:218 > > >> Abort(75) on node 0 (rank 0 in comm 16): application called > > >> MPI_Abort(MPI_COMM_SELF, 75) - process 0 > > >> > > > > > >
