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