Actual source code: itres.c
1: #define PETSCKSP_DLL
3: #include private/kspimpl.h
7: /*@
8: KSPInitialResidual - Computes the residual. Either b - A*C*x with right
9: preconditioning or C*b - C*A*x with left preconditioning; that later
10: residual is often called the "preconditioned residual".
12: Collective on KSP
14: Input Parameters:
15: + vsoln - solution to use in computing residual
16: . vt1, vt2 - temporary work vectors
17: - vb - right-hand-side vector
19: Output Parameters:
20: . vres - calculated residual
22: Notes:
23: This routine assumes that an iterative method, designed for
24: $ A x = b
25: will be used with a preconditioner, C, such that the actual problem is either
26: $ AC u = b (right preconditioning) or
27: $ CA x = Cb (left preconditioning).
28: This means that the calculated residual will be scaled and/or preconditioned;
29: the true residual
30: $ b-Ax
31: is returned in the vt2 temporary.
33: Level: developer
35: .keywords: KSP, residual
37: .seealso: KSPMonitor()
38: @*/
40: PetscErrorCode KSPInitialResidual(KSP ksp,Vec vsoln,Vec vt1,Vec vt2,Vec vres,Vec vb)
41: {
42: MatStructure pflag;
43: Mat Amat,Pmat;
51: if (!ksp->pc) {KSPGetPC(ksp,&ksp->pc);}
52: PCGetOperators(ksp->pc,&Amat,&Pmat,&pflag);
53: if (!ksp->guess_zero) {
54: /* skip right scaling since current guess already has it */
55: KSP_MatMult(ksp,Amat,vsoln,vt1);
56: VecCopy(vb,vt2);
57: VecAXPY(vt2,-1.0,vt1);
58: (ksp->pc_side == PC_RIGHT)?(VecCopy(vt2,vres)):(KSP_PCApply(ksp,vt2,vres));
59: PCDiagonalScaleLeft(ksp->pc,vres,vres);
60: } else {
61: VecCopy(vb,vt2);
62: if (ksp->pc_side == PC_RIGHT) {
63: PCDiagonalScaleLeft(ksp->pc,vb,vres);
64: } else if (ksp->pc_side == PC_LEFT) {
65: KSP_PCApply(ksp,vb,vres);
66: PCDiagonalScaleLeft(ksp->pc,vres,vres);
67: } else if (ksp->pc_side == PC_SYMMETRIC) {
68: PCApplySymmetricLeft(ksp->pc, vb, vres);
69: } else {
70: SETERRQ1(PETSC_ERR_SUP, "Invalid preconditioning side %d", (int)ksp->pc_side);
71: }
72: }
73: return(0);
74: }
78: /*@
79: KSPUnwindPreconditioner - Unwinds the preconditioning in the solution. That is,
80: takes solution to the preconditioned problem and gets the solution to the
81: original problem from it.
83: Collective on KSP
85: Input Parameters:
86: + ksp - iterative context
87: . vsoln - solution vector
88: - vt1 - temporary work vector
90: Output Parameter:
91: . vsoln - contains solution on output
93: Notes:
94: If preconditioning either symmetrically or on the right, this routine solves
95: for the correction to the unpreconditioned problem. If preconditioning on
96: the left, nothing is done.
98: Level: advanced
100: .keywords: KSP, unwind, preconditioner
102: .seealso: KSPSetPreconditionerSide()
103: @*/
104: PetscErrorCode KSPUnwindPreconditioner(KSP ksp,Vec vsoln,Vec vt1)
105: {
111: if (!ksp->pc) {KSPGetPC(ksp,&ksp->pc);}
112: if (ksp->pc_side == PC_RIGHT) {
113: KSP_PCApply(ksp,vsoln,vt1);
114: PCDiagonalScaleRight(ksp->pc,vt1,vsoln);
115: } else if (ksp->pc_side == PC_SYMMETRIC) {
116: PCApplySymmetricRight(ksp->pc,vsoln,vt1);
117: VecCopy(vt1,vsoln);
118: } else {
119: PCDiagonalScaleRight(ksp->pc,vsoln,vsoln);
120: }
121: return(0);
122: }