Actual source code: itres.c
petsc-3.7.7 2017-09-25
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: }