1
/* $Id: kspimpl.h,v 1.50 2001/08/06 18:04:29 bsmith Exp $ */
8
typedef struct _KSPOps *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,int*); /* actual solver */
19
int (*setfromoptions)(KSP);
20
int (*publishoptions)(KSP);
21
int (*computeextremesingularvalues)(KSP,PetscReal*,PetscReal*);
22
int (*computeeigenvalues)(KSP,int,PetscReal*,PetscReal*,int *);
24
int (*view)(KSP,PetscViewer);
28
Maximum number of monitors you can run with a single KSP
30
#define MAXKSPMONITORS 5
33
Defines the KSP data structure.
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
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;
52
Vec vec_sol,vec_rhs; /* pointer to where user has stashed
53
the solution and rhs, these are
54
never touched by the code, only
55
passed back to the user */
56
PetscReal *res_hist; /* If !0 stores residual at iterations*/
57
int res_hist_len; /* current size of residual history array */
58
int res_hist_max; /* actual amount of data in residual_history */
59
PetscTruth res_hist_reset; /* reset history to size zero for each new solve */
61
/* --------User (or default) routines (most return -1 on error) --------*/
62
int (*monitor[MAXKSPMONITORS])(KSP,int,PetscReal,void*); /* returns control to user after */
63
int (*monitordestroy[MAXKSPMONITORS])(void*); /* */
64
void *monitorcontext[MAXKSPMONITORS]; /* residual calculation, allows user */
65
int numbermonitors; /* to, for instance, print residual norm, etc. */
67
int (*converged)(KSP,int,PetscReal,KSPConvergedReason*,void*);
72
void *data; /* holder for misc stuff associated
73
with a particular iterative solver */
75
/* ----------------Default work-area management -------------------- */
81
int its; /* number of iterations so far computed */
83
PetscTruth transpose_solve; /* solve transpose system instead */
85
KSPNormType normtype; /* type of norm used for convergence tests */
88
#define KSPLogResidualHistory(ksp,norm) \
89
{if (ksp->res_hist && ksp->res_hist_max > ksp->res_hist_len) \
90
ksp->res_hist[ksp->res_hist_len++] = norm;}
92
#define KSPMonitor(ksp,it,rnorm) \
93
{ int _ierr,_i,_im = ksp->numbermonitors; \
94
for (_i=0; _i<_im; _i++) {\
95
_ierr = (*ksp->monitor[_i])(ksp,it,rnorm,ksp->monitorcontext[_i]);CHKERRQ(_ierr); \
99
EXTERN int KSPDefaultBuildSolution(KSP,Vec,Vec*);
100
EXTERN int KSPDefaultBuildResidual(KSP,Vec,Vec,Vec *);
101
EXTERN int KSPDefaultDestroy(KSP);
102
EXTERN int KSPDefaultGetWork(KSP,int);
103
EXTERN int KSPDefaultFreeWork(KSP);
104
EXTERN int KSPInitialResidual(KSP,Vec,Vec,Vec,Vec,Vec);
105
EXTERN int KSPUnwindPreconditioner(KSP,Vec,Vec);
108
These allow the various Krylov methods to apply to either the linear system
111
#define KSP_MatMult(ksp,A,x,y) (!ksp->transpose_solve) ? MatMult(A,x,y) : MatMultTranspose(A,x,y)
112
#define KSP_MatMultTranspose(ksp,A,x,y) (!ksp->transpose_solve) ? MatMultTranspose(A,x,y) : MatMult(A,x,y)
113
#define KSP_PCApply(ksp,B,x,y) (!ksp->transpose_solve) ? PCApply(B,x,y) : PCApplyTranspose(B,x,y)
114
#define KSP_PCApplyTranspose(ksp,B,x,y) (!ksp->transpose_solve) ? PCApplyTranspose(B,x,y) : PCApply(B,x,y)
115
#define KSP_PCApplyBAorAB(ksp,pc,side,x,y,w) (!ksp->transpose_solve) ? PCApplyBAorAB(pc,side,x,y,w) : PCApplyBAorABTranspose(pc,side,x,y,w)