Actual source code: itres.c
2: #include src/ksp/ksp/kspimpl.h
6: /*@C
7: KSPInitialResidual - Computes the residual.
9: Collective on KSP
11: Input Parameters:
12: + vsoln - solution to use in computing residual
13: . vt1, vt2 - temporary work vectors
14: . vres - calculated residual
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: PetscErrorCode KSPInitialResidual(KSP ksp,Vec vsoln,Vec vt1,Vec vt2,Vec vres,Vec vb)
31: {
32: PetscScalar mone = -1.0;
33: MatStructure pflag;
34: Mat Amat,Pmat;
42: PCGetOperators(ksp->pc,&Amat,&Pmat,&pflag);
43: if (!ksp->guess_zero) {
44: /* skip right scaling since current guess already has it */
45: KSP_MatMult(ksp,Amat,vsoln,vt1);
46: VecCopy(vb,vt2);
47: VecAXPY(&mone,vt1,vt2);
48: (ksp->pc_side == PC_RIGHT)?(VecCopy(vt2,vres)):(KSP_PCApply(ksp,vt2,vres));
49: PCDiagonalScaleLeft(ksp->pc,vres,vres);
50: } else {
51: if (ksp->pc_side == PC_RIGHT) {
52: PCDiagonalScaleLeft(ksp->pc,vb,vres);
53: } else if (ksp->pc_side == PC_LEFT) {
54: KSP_PCApply(ksp,vb,vres);
55: PCDiagonalScaleLeft(ksp->pc,vres,vres);
56: } else if (ksp->pc_side == PC_SYMMETRIC) {
57: PCApplySymmetricLeft(ksp->pc, vb, vres);
58: } else {
59: SETERRQ1(PETSC_ERR_SUP, "Invalid preconditioning side %d", (int)ksp->pc_side);
60: }
62: }
63: return(0);
64: }
68: /*@
69: KSPUnwindPreconditioner - Unwinds the preconditioning in the solution. That is,
70: takes solution to the preconditioned problem and gets the solution to the
71: original problem from it.
73: Collective on KSP
75: Input Parameters:
76: + ksp - iterative context
77: . vsoln - solution vector
78: - vt1 - temporary work vector
80: Output Parameter:
81: . vsoln - contains solution on output
83: Notes:
84: If preconditioning either symmetrically or on the right, this routine solves
85: for the correction to the unpreconditioned problem. If preconditioning on
86: the left, nothing is done.
88: Level: advanced
90: .keywords: KSP, unwind, preconditioner
92: .seealso: KSPSetPreconditionerSide()
93: @*/
94: PetscErrorCode KSPUnwindPreconditioner(KSP ksp,Vec vsoln,Vec vt1)
95: {
101: if (ksp->pc_side == PC_RIGHT) {
102: KSP_PCApply(ksp,vsoln,vt1);
103: PCDiagonalScaleRight(ksp->pc,vt1,vsoln);
104: } else if (ksp->pc_side == PC_SYMMETRIC) {
105: PCApplySymmetricRight(ksp->pc,vsoln,vt1);
106: VecCopy(vt1,vsoln);
107: } else {
108: PCDiagonalScaleRight(ksp->pc,vsoln,vsoln);
109: }
110: return(0);
111: }