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