Actual source code: itres.c
1: /*$Id: itres.c,v 1.54 2001/08/07 03:03:45 balay Exp $*/
3: #include src/sles/ksp/kspimpl.h
5: /*@C
6: KSPInitialResidual - Computes the residual.
8: Collective on KSP
10: Input Parameters:
11: + vsoln - solution to use in computing residual
12: . vt1, vt2 - temporary work vectors
13: . vres - calculated residual
14: - vb - right-hand-side vector
16: Notes:
17: This routine assumes that an iterative method, designed for
18: $ A x = b
19: will be used with a preconditioner, C, such that the actual problem is either
20: $ AC u = f (right preconditioning) or
21: $ CA x = Cf (left preconditioning).
23: Level: developer
25: .keywords: KSP, residual
27: .seealso: KSPMonitor()
28: @*/
29: int KSPInitialResidual(KSP ksp,Vec vsoln,Vec vt1,Vec vt2,Vec vres,Vec vb)
30: {
31: PetscScalar mone = -1.0;
32: MatStructure pflag;
33: Mat Amat,Pmat;
34: int ierr;
38: PCGetOperators(ksp->B,&Amat,&Pmat,&pflag);
39: if (!ksp->guess_zero) {
40: /* skip right scaling since current guess already has it */
41: KSP_MatMult(ksp,Amat,vsoln,vt1);
42: VecCopy(vb,vt2);
43: VecAXPY(&mone,vt1,vt2);
44: (ksp->pc_side == PC_RIGHT)?(VecCopy(vt2,vres)):(KSP_PCApply(ksp,ksp->B,vt2,vres));
45: PCDiagonalScaleLeft(ksp->B,vres,vres);
46: } else {
47: if (ksp->pc_side == PC_RIGHT) {
48: PCDiagonalScaleLeft(ksp->B,vb,vres);
49: } else if (ksp->pc_side == PC_LEFT) {
50: KSP_PCApply(ksp,ksp->B,vb,vres);
51: PCDiagonalScaleLeft(ksp->B,vres,vres);
52: } else if (ksp->pc_side == PC_SYMMETRIC) {
53: PCApplySymmetricLeft(ksp->B, vb, vres);
54: } else {
55: SETERRQ1(PETSC_ERR_SUP, "Invalid preconditioning side %d", ksp->pc_side);
56: }
58: }
59: return(0);
60: }
62: /*@
63: KSPUnwindPreconditioner - Unwinds the preconditioning in the solution. That is,
64: takes solution to the preconditioned problem and gets the solution to the
65: original problem from it.
67: Collective on KSP
69: Input Parameters:
70: + ksp - iterative context
71: . vsoln - solution vector
72: - vt1 - temporary work vector
74: Output Parameter:
75: . vsoln - contains solution on output
77: Notes:
78: If preconditioning either symmetrically or on the right, this routine solves
79: for the correction to the unpreconditioned problem. If preconditioning on
80: the left, nothing is done.
82: Level: advanced
84: .keywords: KSP, unwind, preconditioner
86: .seealso: KSPSetPreconditionerSide()
87: @*/
88: int KSPUnwindPreconditioner(KSP ksp,Vec vsoln,Vec vt1)
89: {
94: if (ksp->pc_side == PC_RIGHT) {
95: KSP_PCApply(ksp,ksp->B,vsoln,vt1);
96: PCDiagonalScaleRight(ksp->B,vt1,vsoln);
97: } else if (ksp->pc_side == PC_SYMMETRIC) {
98: PCApplySymmetricRight(ksp->B,vsoln,vt1);
99: VecCopy(vt1,vsoln);
100: } else {
101: PCDiagonalScaleRight(ksp->B,vsoln,vsoln);
102: }
103: return(0);
104: }