Actual source code: snesngmres.h
1: #ifndef _SNESNGMRES_H
2: #define _SNESNGMRES_H
4: #include <petsc/private/snesimpl.h>
6: /* Data structure for the Nonlinear GMRES method. */
7: typedef struct {
8: /* Solver parameters and counters */
9: PetscInt msize; /* maximum size of krylov space */
10: PetscInt restart_it; /* number of iterations the restart conditions persist before restart */
11: PetscViewer monitor; /* debugging output for NGMRES */
12: PetscInt restart_periodic; /* number of iterations to restart after */
14: SNESNGMRESRestartType restart_type;
15: SNESNGMRESSelectType select_type;
17: /* History and subspace data */
18: Vec *Fdot; /* residual history -- length msize */
19: Vec *Xdot; /* solution history -- length msize */
20: PetscReal *fnorms; /* the residual norm history */
21: PetscReal *xnorms; /* the solution norm history */
23: /* General minimization problem context */
24: PetscScalar *h; /* the constraint matrix */
25: PetscScalar *beta; /* rhs for the minimization problem */
26: PetscScalar *xi; /* the dot-product of the current and previous res. */
28: /* Line searches */
29: SNESLineSearch additive_linesearch; /* Line search for the additive variant */
31: /* Selection constants */
32: PetscBool candidate; /* use candidate storage approach */
33: PetscBool approxfunc; /* approximate the function rather than recomputing it */
34: PetscBool singlereduction; /* use a single reduction (with more local work) for tolerance selection */
35: PetscReal gammaA; /* Criterion A residual tolerance */
36: PetscReal epsilonB; /* Criterion B difference tolerance */
37: PetscReal deltaB; /* Criterion B residual tolerance */
38: PetscReal gammaC; /* Restart residual tolerance */
39: PetscBool restart_fm_rise; /* Restart on F_M residual increase */
41: PetscReal andersonBeta; /* Relaxation parameter for Anderson Mixing */
43: /* Least squares minimization solve context */
44: PetscScalar *q; /* the matrix formed as q_ij = (rdot_i, rdot_j) */
45: PetscBLASInt m; /* matrix dimension */
46: PetscBLASInt n; /* matrix dimension */
47: PetscBLASInt nrhs; /* the number of right hand sides */
48: PetscBLASInt lda; /* the padded matrix dimension */
49: PetscBLASInt ldb; /* the padded vector dimension */
50: PetscReal *s; /* the singular values */
51: PetscReal rcond; /* the exit condition */
52: PetscBLASInt rank; /* the effective rank */
53: PetscScalar *work; /* the work vector */
54: PetscReal *rwork; /* the real work vector used for complex */
55: PetscBLASInt lwork; /* the size of the work vector */
56: PetscBLASInt info; /* the output condition */
58: PetscBool setup_called; /* indicates whether SNESSetUp_NGMRES() has been called */
59: } SNES_NGMRES;
61: #define H(i, j) ngmres->h[i * ngmres->msize + j]
62: #define Q(i, j) ngmres->q[i * ngmres->msize + j]
64: /* private functions that are shared components of the methods */
65: PETSC_INTERN PetscErrorCode SNESNGMRESUpdateSubspace_Private(SNES, PetscInt, PetscInt, Vec, PetscReal, Vec);
66: PETSC_INTERN PetscErrorCode SNESNGMRESFormCombinedSolution_Private(SNES, PetscInt, PetscInt, Vec, Vec, PetscReal, Vec, Vec, Vec);
67: PETSC_INTERN PetscErrorCode SNESNGMRESNorms_Private(SNES, PetscInt, Vec, Vec, Vec, Vec, Vec, Vec, Vec, PetscReal *, PetscReal *, PetscReal *, PetscReal *, PetscReal *, PetscReal *, PetscReal *, PetscReal *);
68: PETSC_INTERN PetscErrorCode SNESNGMRESSelect_Private(SNES, PetscInt, Vec, Vec, PetscReal, PetscReal, PetscReal, Vec, Vec, PetscReal, PetscReal, PetscReal, PetscReal, PetscReal, PetscReal, Vec, Vec, Vec, PetscReal *, PetscReal *, PetscReal *);
69: PETSC_INTERN PetscErrorCode SNESNGMRESSelectRestart_Private(SNES, PetscInt, PetscReal, PetscReal, PetscReal, PetscReal, PetscReal, PetscBool *);
71: PETSC_INTERN PetscErrorCode SNESDestroy_NGMRES(SNES);
72: PETSC_INTERN PetscErrorCode SNESReset_NGMRES(SNES);
73: PETSC_INTERN PetscErrorCode SNESSetUp_NGMRES(SNES);
74: PETSC_INTERN PetscErrorCode SNESView_NGMRES(SNES, PetscViewer);
76: #endif