Actual source code: kspimpl.h
2: #ifndef _KSPIMPL_H
3: #define _KSPIMPL_H
5: #include petscksp.h
7: typedef struct _KSPOps *KSPOps;
9: struct _KSPOps {
10: PetscErrorCode (*buildsolution)(KSP,Vec,Vec*); /* Returns a pointer to the solution, or
11: calculates the solution in a
12: user-provided area. */
13: PetscErrorCode (*buildresidual)(KSP,Vec,Vec,Vec*); /* Returns a pointer to the residual, or
14: calculates the residual in a
15: user-provided area. */
16: PetscErrorCode (*solve)(KSP); /* actual solver */
17: PetscErrorCode (*setup)(KSP);
18: PetscErrorCode (*setfromoptions)(KSP);
19: PetscErrorCode (*publishoptions)(KSP);
20: PetscErrorCode (*computeextremesingularvalues)(KSP,PetscReal*,PetscReal*);
21: PetscErrorCode (*computeeigenvalues)(KSP,PetscInt,PetscReal*,PetscReal*,PetscInt *);
22: PetscErrorCode (*destroy)(KSP);
23: PetscErrorCode (*view)(KSP,PetscViewer);
24: };
26: /*
27: Maximum number of monitors you can run with a single KSP
28: */
29: #define MAXKSPMONITORS 5
31: /*
32: Defines the KSP data structure.
33: */
34: struct _p_KSP {
35: PETSCHEADER(struct _KSPOps);
36: /*------------------------- User parameters--------------------------*/
37: PetscInt max_it; /* maximum number of iterations */
38: PetscTruth guess_zero, /* flag for whether initial guess is 0 */
39: calc_sings, /* calculate extreme Singular Values */
40: guess_knoll; /* use initial guess of PCApply(ksp->B,b */
41: PCSide pc_side; /* flag for left, right, or symmetric
42: preconditioning */
43: PetscReal rtol, /* relative tolerance */
44: abstol, /* absolute tolerance */
45: ttol, /* (not set by user) */
46: divtol; /* divergence tolerance */
47: PetscTruth defaultconvergedinitialrtol; /* default relative residual decrease is computing from initial residual, not rhs */
48: PetscTruth defaultconvergedmininitialrtol; /* default relative residual decrease is computing from min of initial residual and rhs */
49: PetscReal rnorm0; /* initial residual norm (used for divergence testing) */
50: PetscReal rnorm; /* current residual norm */
51: KSPConvergedReason reason;
52: PetscTruth printreason; /* prints converged reason after solve */
54: Vec vec_sol,vec_rhs; /* pointer to where user has stashed
55: the solution and rhs, these are
56: never touched by the code, only
57: passed back to the user */
58: PetscReal *res_hist; /* If !0 stores residual at iterations*/
59: PetscReal *res_hist_alloc; /* If !0 means user did not provide buffer, needs deallocation */
60: PetscInt res_hist_len; /* current size of residual history array */
61: PetscInt res_hist_max; /* actual amount of data in residual_history */
62: PetscTruth res_hist_reset; /* reset history to size zero for each new solve */
64: /* --------User (or default) routines (most return -1 on error) --------*/
65: PetscErrorCode (*monitor[MAXKSPMONITORS])(KSP,PetscInt,PetscReal,void*); /* returns control to user after */
66: PetscErrorCode (*monitordestroy[MAXKSPMONITORS])(void*); /* */
67: void *monitorcontext[MAXKSPMONITORS]; /* residual calculation, allows user */
68: PetscInt numbermonitors; /* to, for instance, print residual norm, etc. */
70: PetscErrorCode (*converged)(KSP,PetscInt,PetscReal,KSPConvergedReason*,void*);
71: void *cnvP;
73: PC pc;
75: void *data; /* holder for misc stuff associated
76: with a particular iterative solver */
78: /* ----------------Default work-area management -------------------- */
79: PetscInt nwork;
80: Vec *work;
82: PetscInt setupcalled;
84: PetscInt its; /* number of iterations so far computed */
86: PetscTruth transpose_solve; /* solve transpose system instead */
88: KSPNormType normtype; /* type of norm used for convergence tests */
90: /* Allow diagonally scaling the matrix before computing the preconditioner or using
91: the Krylov method. Note this is NOT just Jacobi preconditioning */
93: PetscTruth dscale; /* diagonal scale system; used with KSPSetDiagonalScale() */
94: PetscTruth dscalefix; /* unscale system after solve */
95: PetscTruth dscalefix2; /* system has been unscaled */
96: Vec diagonal; /* 1/sqrt(diag of matrix) */
97: Vec truediagonal;
99: MatNullSpace nullsp; /* Null space of the operator, removed from Krylov space */
100: };
102: #define KSPLogResidualHistory(ksp,norm) \
103: {if (ksp->res_hist && ksp->res_hist_max > ksp->res_hist_len) \
104: ksp->res_hist[ksp->res_hist_len++] = norm;}
106: #define KSPMonitor(ksp,it,rnorm) \
107: { PetscErrorCode _ierr; PetscInt _i,_im = ksp->numbermonitors; \
108: for (_i=0; _i<_im; _i++) {\
109: _(*ksp->monitor[_i])(ksp,it,rnorm,ksp->monitorcontext[_i]);CHKERRQ(_ierr); \
110: } \
111: }
113: EXTERN PetscErrorCode KSPDefaultDestroy(KSP);
114: EXTERN PetscErrorCode KSPGetVecs(KSP,PetscInt,Vec**,PetscInt,Vec**);
115: EXTERN PetscErrorCode KSPDefaultGetWork(KSP,PetscInt);
116: EXTERN PetscErrorCode KSPDefaultFreeWork(KSP);
118: /*
119: These allow the various Krylov methods to apply to either the linear system or its transpose.
120: */
121: #define KSP_RemoveNullSpace(ksp,y) ((ksp->nullsp && ksp->pc_side == PC_LEFT) ? MatNullSpaceRemove(ksp->nullsp,y,PETSC_NULL) : 0)
123: #define KSP_MatMult(ksp,A,x,y) (!ksp->transpose_solve) ? MatMult(A,x,y) : MatMultTranspose(A,x,y)
124: #define KSP_MatMultTranspose(ksp,A,x,y) (!ksp->transpose_solve) ? MatMultTranspose(A,x,y) : MatMult(A,x,y)
125: #define KSP_PCApply(ksp,x,y) (!ksp->transpose_solve) ? (PCApply(ksp->pc,x,y) || KSP_RemoveNullSpace(ksp,y)) : PCApplyTranspose(ksp->pc,x,y)
126: #define KSP_PCApplyTranspose(ksp,x,y) (!ksp->transpose_solve) ? PCApplyTranspose(ksp->pc,x,y) : (PCApply(ksp->pc,x,y) || KSP_RemoveNullSpace(ksp,y))
127: #define KSP_PCApplyBAorAB(ksp,x,y,w) (!ksp->transpose_solve) ? (PCApplyBAorAB(ksp->pc,ksp->pc_side,x,y,w) || KSP_RemoveNullSpace(ksp,y)) : PCApplyBAorABTranspose(ksp->pc,ksp->pc_side,x,y,w)
131: #endif