[petsc-users] Least Square solver LSQR

[Natacha BEREUX]

Dear  PETSc users,
I am trying to solve an over determined linear system of equations Ax =b
with the least square solver LSQR.
The matrix A is rectangular : it has m rows and n columns with m>n.
I set the preconditioner to "PCNONE"  (following for instance post
https://lists.mcs.anl.gov/mailman/htdig/petsc-users/2012-May/013591.html)

Matrix A has the form :
A(1:n,1:n)        = Identity(n)
A(n+1:m, 1:n) = 0.

I define a vector u of size n and set u(1:n) to 1.0

The right hand side b= Au.

The (fortran) test program is attached (it's
petsc-3.5.2/src/ksp/ksp/examples/tutorials/ex1f.F slightly modified)

The size of matrix A is defined at runtime

The solver ends OK with m=10, n=2 but goes wrong  with  m=10, n=3 (KSP
Converged reason = -5)

The output is shown below:

First run with n=2 columns
====================
./ex1f -ksp_type lsqr -vec_type seq -ksp_monitor -m 10 -n 2

  0 KSP Residual norm 1.414213562373e+00
  1 KSP Residual norm 3.790370795009e-16
KSP Object: 1 MPI processes
  type: lsqr
  maximum iterations=10000, initial guess is zero
  tolerances:  relative=1e-07, absolute=1e-50, divergence=10000
  left preconditioning
  using UNPRECONDITIONED norm type for convergence test
PC Object: 1 MPI processes
  type: none
  linear system matrix = precond matrix:
  Mat Object:   1 MPI processes
    type: seqaij
    rows=10, cols=2
    total: nonzeros=2, allocated nonzeros=50
    total number of mallocs used during MatSetValues calls =0
      using I-node routines: found 4 nodes, limit used is 5
 KSP Converged Reason  2
Norm of error < 1.e-12,Iterations =     1

Second run, with n=3
================
./ex1f -ksp_type lsqr -vec_type seq -ksp_monitor -m 10 -n 3

 0 KSP Residual norm 1.732050807569e+00
KSP Object: 1 MPI processes
  type: lsqr
  maximum iterations=10000, initial guess is zero
  tolerances:  relative=1e-07, absolute=1e-50, divergence=10000
  left preconditioning
  using UNPRECONDITIONED norm type for convergence test
PC Object: 1 MPI processes
  type: none
  linear system matrix = precond matrix:
  Mat Object:   1 MPI processes
    type: seqaij
    rows=10, cols=3
    total: nonzeros=3, allocated nonzeros=50
    total number of mallocs used during MatSetValues calls =0
      using I-node routines: found 5 nodes, limit used is 5
 KSP Converged Reason -5
Norm of error =  0.1732E+01,  Iterations =     0


I am using version 3.5.2 of PETSc library

I am probably doing something wrong, but I don't understand what the
problem is.
Does anyone have an idea of what is going on ?
Best regards
Natacha
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20150106/c241d5e7/attachment.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: ex1f.F
Type: text/x-fortran
Size: 7436 bytes
Desc: not available
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20150106/c241d5e7/attachment.bin>