Actual source code: itres.c
2: #include <petsc/private/kspimpl.h>
4: /*@
5: KSPInitialResidual - Computes the residual. Either b - A*C*u = b - A*x with right
6: preconditioning or C*(b - A*x) with left preconditioning; the latter
7: residual is often called the "preconditioned residual".
9: Collective
11: Input Parameters:
12: + vsoln - solution to use in computing residual
13: . vt1, vt2 - temporary work vectors
14: - vb - right-hand-side vector
16: Output Parameters:
17: . vres - calculated residual
19: Level: developer
21: Note:
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: .seealso: [](chapter_ksp), `KSP`, `KSPSolve()`, `KSPMonitor()`
33: @*/
35: PetscErrorCode KSPInitialResidual(KSP ksp, Vec vsoln, Vec vt1, Vec vt2, Vec vres, Vec vb)
36: {
37: Mat Amat, Pmat;
43: if (!ksp->pc) KSPGetPC(ksp, &ksp->pc);
44: PCGetOperators(ksp->pc, &Amat, &Pmat);
45: if (!ksp->guess_zero) {
46: /* skip right scaling since current guess already has it */
47: KSP_MatMult(ksp, Amat, vsoln, vt1);
48: VecCopy(vb, vt2);
49: VecAXPY(vt2, -1.0, vt1);
50: if (ksp->pc_side == PC_RIGHT) {
51: PCDiagonalScaleLeft(ksp->pc, vt2, vres);
52: } else if (ksp->pc_side == PC_LEFT) {
53: KSP_PCApply(ksp, vt2, vres);
54: PCDiagonalScaleLeft(ksp->pc, vres, vres);
55: } else if (ksp->pc_side == PC_SYMMETRIC) {
56: PCApplySymmetricLeft(ksp->pc, vt2, vres);
57: } else SETERRQ(PetscObjectComm((PetscObject)ksp), PETSC_ERR_SUP, "Invalid preconditioning side %d", (int)ksp->pc_side);
58: } else {
59: VecCopy(vb, vt2);
60: if (ksp->pc_side == PC_RIGHT) {
61: PCDiagonalScaleLeft(ksp->pc, vb, vres);
62: } else if (ksp->pc_side == PC_LEFT) {
63: KSP_PCApply(ksp, vb, vres);
64: PCDiagonalScaleLeft(ksp->pc, vres, vres);
65: } else if (ksp->pc_side == PC_SYMMETRIC) {
66: PCApplySymmetricLeft(ksp->pc, vb, vres);
67: } else SETERRQ(PetscObjectComm((PetscObject)ksp), PETSC_ERR_SUP, "Invalid preconditioning side %d", (int)ksp->pc_side);
68: }
69: /* This may be true only on a subset of MPI ranks; setting it here so it will be detected by the first norm computation in the Krylov method */
70: if (ksp->reason == KSP_DIVERGED_PC_FAILED) VecSetInf(vres);
71: return 0;
72: }
74: /*@
75: KSPUnwindPreconditioner - Unwinds the preconditioning in the solution. That is,
76: takes solution to the preconditioned problem and gets the solution to the
77: original problem from it.
79: Collective
81: Input Parameters:
82: + ksp - iterative context
83: . vsoln - solution vector
84: - vt1 - temporary work vector
86: Output Parameter:
87: . vsoln - contains solution on output
89: Note:
90: If preconditioning either symmetrically or on the right, this routine solves
91: for the correction to the unpreconditioned problem. If preconditioning on
92: the left, nothing is done.
94: Level: advanced
96: .seealso: [](chapter_ksp), `KSP`, `KSPSetPCSide()`
97: @*/
98: PetscErrorCode KSPUnwindPreconditioner(KSP ksp, Vec vsoln, Vec vt1)
99: {
102: if (!ksp->pc) KSPGetPC(ksp, &ksp->pc);
103: if (ksp->pc_side == PC_RIGHT) {
104: KSP_PCApply(ksp, vsoln, vt1);
105: PCDiagonalScaleRight(ksp->pc, vt1, vsoln);
106: } else if (ksp->pc_side == PC_SYMMETRIC) {
107: PCApplySymmetricRight(ksp->pc, vsoln, vt1);
108: VecCopy(vt1, vsoln);
109: } else {
110: PCDiagonalScaleRight(ksp->pc, vsoln, vsoln);
111: }
112: return 0;
113: }