2: #include <petsc/private/kspimpl.h>   /*I "petscksp.h" I*/

  6: /*@
  7:    KSPInitialResidual - Computes the residual. Either b - A*C*u = b - A*x with right
  8:      preconditioning or C*(b - A*x) with left preconditioning; that later
  9:      residual is often called the "preconditioned residual".

 11:    Collective on KSP

 13:    Input Parameters:
 14: +  vsoln    - solution to use in computing residual
 15: .  vt1, vt2 - temporary work vectors
 16: -  vb       - right-hand-side vector

 18:    Output Parameters:
 19: .  vres     - calculated residual

 21:    Notes:
 22:    This routine assumes that an iterative method, designed for
 23: $     A x = b
 24:    will be used with a preconditioner, C, such that the actual problem is either
 25: $     AC u = b (right preconditioning) or
 26: $     CA x = Cb (left preconditioning).
 27:    This means that the calculated residual will be scaled and/or preconditioned;
 28:    the true residual
 29: $     b-Ax
 30:    is returned in the vt2 temporary.

 32:    Level: developer

 34: .keywords: KSP, residual

 36: .seealso:  KSPMonitor()
 37: @*/

 39: PetscErrorCode  KSPInitialResidual(KSP ksp,Vec vsoln,Vec vt1,Vec vt2,Vec vres,Vec vb)
 40: {
 41:   Mat            Amat,Pmat;

 49:   if (!ksp->pc) {KSPGetPC(ksp,&ksp->pc);}
 50:   PCGetOperators(ksp->pc,&Amat,&Pmat);
 51:   if (!ksp->guess_zero) {
 52:     /* skip right scaling since current guess already has it */
 53:     KSP_MatMult(ksp,Amat,vsoln,vt1);
 54:     VecCopy(vb,vt2);
 55:     VecAXPY(vt2,-1.0,vt1);
 56:     if (ksp->pc_side == PC_RIGHT) {
 57:       PCDiagonalScaleLeft(ksp->pc,vt2,vres);
 58:     } else if (ksp->pc_side == PC_LEFT) {
 59:         KSP_PCApply(ksp,vt2,vres);
 60:         PCDiagonalScaleLeft(ksp->pc,vres,vres);
 61:     } else if (ksp->pc_side == PC_SYMMETRIC) {
 62:       PCApplySymmetricLeft(ksp->pc,vt2,vres);
 63:     } else SETERRQ1(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP, "Invalid preconditioning side %d", (int)ksp->pc_side);
 64:   } else {
 65:     VecCopy(vb,vt2);
 66:     if (ksp->pc_side == PC_RIGHT) {
 67:       PCDiagonalScaleLeft(ksp->pc,vb,vres);
 68:     } else if (ksp->pc_side == PC_LEFT) {
 69:       KSP_PCApply(ksp,vb,vres);
 70:       PCDiagonalScaleLeft(ksp->pc,vres,vres);
 71:     } else if (ksp->pc_side == PC_SYMMETRIC) {
 72:       PCApplySymmetricLeft(ksp->pc, vb, vres);
 73:     } else SETERRQ1(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP, "Invalid preconditioning side %d", (int)ksp->pc_side);
 74:   }
 75:   return(0);
 76: }

 80: /*@
 81:    KSPUnwindPreconditioner - Unwinds the preconditioning in the solution. That is,
 82:      takes solution to the preconditioned problem and gets the solution to the
 83:      original problem from it.

 85:    Collective on KSP

 87:    Input Parameters:
 88: +  ksp  - iterative context
 89: .  vsoln - solution vector
 90: -  vt1   - temporary work vector

 92:    Output Parameter:
 93: .  vsoln - contains solution on output

 95:    Notes:
 96:    If preconditioning either symmetrically or on the right, this routine solves
 97:    for the correction to the unpreconditioned problem.  If preconditioning on
 98:    the left, nothing is done.

100:    Level: advanced

102: .keywords: KSP, unwind, preconditioner

104: .seealso: KSPSetPCSide()
105: @*/
106: PetscErrorCode  KSPUnwindPreconditioner(KSP ksp,Vec vsoln,Vec vt1)
107: {

113:   if (!ksp->pc) {KSPGetPC(ksp,&ksp->pc);}
114:   if (ksp->pc_side == PC_RIGHT) {
115:     KSP_PCApply(ksp,vsoln,vt1);
116:     PCDiagonalScaleRight(ksp->pc,vt1,vsoln);
117:   } else if (ksp->pc_side == PC_SYMMETRIC) {
118:     PCApplySymmetricRight(ksp->pc,vsoln,vt1);
119:     VecCopy(vt1,vsoln);
120:   } else {
121:     PCDiagonalScaleRight(ksp->pc,vsoln,vsoln);
122:   }
123:   return(0);
124: }