[petsc-users] Suggestion for KSPRichardson

[Gaetan Kenway]

Hi Everyone

I have a small suggestion for the KSPRichardson implementation.  I often
use a KSPRichardson object to perform a small number (2-4) applications of
a preconditioner (Additive Schwartz for example). Since I'm using a fixed
number of iterations, the object uses KSP_NORM_NONE to save the vecNorm.

Near the end of the function KSPSolve_Richardson, (line 138 according to
the current docs, link below), the new residual is evaluated. However, for
the last iteration, when using KSP_NORM_NONE, the last matMult and VecAYPX
is not necessary. This results in additional unnecessary matMults(). In one
system I'm solving, this adds approximately 20% to the computation cost.
 (30 secs -> 36 secs)

I think the following code should work:
 if (ksp->its < maxit  || ksp->normtype != KSP_NORM_NONE){
        ierr = KSP_MatMult(ksp,Amat,x,r);CHKERRQ(ierr);      /*   r  <- b -
Ax      */
        ierr = VecAYPX(r,-1.0,b);CHKERRQ(ierr);
      }

Also, it makes the degenerate case of 1 iteration have same cost as if the
KSPRichardson wasn't there.

Hopefully this is useful.

Gaetan Kenway

Link to source code:
http://www.mcs.anl.gov/petsc/petsc-current/src/ksp/ksp/impls/rich/rich.c.html#KSPRICHARDSON
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20130102/fd0ca9c1/attachment.html>