Actual source code: itres.c

  1: /*$Id: itres.c,v 1.54 2001/08/07 03:03:45 balay Exp $*/

 3:  #include src/sles/ksp/kspimpl.h

  5: /*@C
  6:    KSPInitialResidual - Computes the residual.

  8:    Collective on KSP

 10:    Input Parameters:
 11: +  vsoln    - solution to use in computing residual
 12: .  vt1, vt2 - temporary work vectors
 13: .  vres     - calculated residual
 14: .  vbinvf   - the result of binv^{-1} b.  If null, don't do it.
 15: -  vb       - right-hand-side vector

 17:    Notes:
 18:    This routine assumes that an iterative method, designed for
 19: $     A x = b
 20:    will be used with a preconditioner, C, such that the actual problem is either
 21: $     AC u = f (right preconditioning) or
 22: $     CA x = Cf (left preconditioning).

 24:    Level: developer

 26: .keywords: KSP, residual

 28: .seealso:  KSPMonitor()
 29: @*/
 30: int KSPInitialResidual(KSP ksp,Vec vsoln,Vec vt1,Vec vt2,Vec vres,Vec vbinvf,Vec vb)
 31: {
 32:   PetscScalar   mone = -1.0;
 33:   MatStructure  pflag;
 34:   Mat           Amat,Pmat;
 35:   int           ierr;

 39:   PCGetOperators(ksp->B,&Amat,&Pmat,&pflag);
 40:   if (ksp->pc_side == PC_RIGHT) {
 41:     PCDiagonalScaleLeft(ksp->B,vb,vbinvf);
 42:   } else if (ksp->pc_side == PC_LEFT) {
 43:     KSP_PCApply(ksp,ksp->B,vb,vbinvf);
 44:     PCDiagonalScaleLeft(ksp->B,vbinvf,vbinvf);
 45:   } else if (ksp->pc_side == PC_SYMMETRIC) {
 46:     PCApplySymmetricLeft(ksp->B, vb, vbinvf);
 47:   } else {
 48:     SETERRQ1(PETSC_ERR_SUP, "Invalid preconditioning side %d", ksp->pc_side);
 49:   }
 50:   if (!ksp->guess_zero) {
 51:     /* skip right scaling since current guess already has it */
 52:     KSP_MatMult(ksp,Amat,vsoln,vt1);
 53:     VecCopy(vb,vt2);
 54:     VecAXPY(&mone,vt1,vt2);
 55:     (ksp->pc_side == PC_RIGHT)?(VecCopy(vt2,vres)):(KSP_PCApply(ksp,ksp->B,vt2,vres));
 56:     PCDiagonalScaleLeft(ksp->B,vres,vres);
 57:   } else {
 58:     VecCopy(vbinvf,vres);
 59:   }
 60:   return(0);
 61: }

 63: /*@
 64:    KSPUnwindPreconditioner - Unwinds the preconditioning in the solution. That is,
 65:      takes solution to the preconditioned problem and gets the solution to the 
 66:      original problem from it.

 68:    Collective on KSP

 70:    Input Parameters:
 71: +  ksp  - iterative context
 72: .  vsoln - solution vector 
 73: -  vt1   - temporary work vector

 75:    Output Parameter:
 76: .  vsoln - contains solution on output  

 78:    Notes:
 79:    If preconditioning either symmetrically or on the right, this routine solves 
 80:    for the correction to the unpreconditioned problem.  If preconditioning on 
 81:    the left, nothing is done.

 83:    Level: advanced

 85: .keywords: KSP, unwind, preconditioner

 87: .seealso: KSPSetPreconditionerSide()
 88: @*/
 89: int KSPUnwindPreconditioner(KSP ksp,Vec vsoln,Vec vt1)
 90: {

 95:   if (ksp->pc_side == PC_RIGHT) {
 96:     KSP_PCApply(ksp,ksp->B,vsoln,vt1);
 97:     PCDiagonalScaleRight(ksp->B,vt1,vsoln);
 98:   } else if (ksp->pc_side == PC_SYMMETRIC) {
 99:     PCApplySymmetricRight(ksp->B,vsoln,vt1);
100:     VecCopy(vt1,vsoln);
101:   } else {
102:     PCDiagonalScaleRight(ksp->B,vsoln,vsoln);
103:   }
104:   return(0);
105: }