Actual source code: gltrimpl.h
1: /*****************************************************************************/
2: /* Context for using preconditioned conjugate gradient method to minimized a */
3: /* quadratic function subject to a trust region constraint. If the matrix */
4: /* is indefinite, a direction of negative curvature may be encountered. If */
5: /* a direction of negative curvature is found, we continue to build the */
6: /* tridiagonal Lanczos matrix for a fixed number of iterations. After this */
7: /* matrix is computed, we compute a global solution to solve the trust- */
8: /* region problem with the tridiagonal approximation by using a variant of */
9: /* the More'-Sorenson algorithm. The direction is then constructed from */
10: /* this solution. */
11: /* */
12: /* This method is described in: */
13: /* N. Gould, S. Lucidi, M. Roma, and Ph. Toint, "Solving the Trust-Region */
14: /* Subproblem using the Lanczos Method", SIAM Journal on Optimization, */
15: /* 9, pages 504-525, 1999. */
16: /*****************************************************************************/
18: #ifndef __GLTR
21: typedef struct {
22: PetscReal *diag; /* Diagonal part of Lanczos matrix */
23: PetscReal *offd; /* Off-diagonal part of Lanczos matrix */
24: PetscReal *alpha; /* Record of alpha values from CG */
25: PetscReal *beta; /* Record of beta values from CG */
26: PetscReal *norm_r; /* Record of residual values from CG */
28: PetscReal *rwork; /* Real workspace for solver computations */
29: PetscBLASInt *iwork; /* Integer workspace for solver computations */
31: PetscReal radius;
32: PetscReal norm_d;
33: PetscReal e_min;
34: PetscReal o_fcn;
35: PetscReal lambda;
37: PetscReal init_pert; /* Initial perturbation for solve */
38: PetscReal eigen_tol; /* Tolerance used when computing eigenvalue */
39: PetscReal newton_tol; /* Tolerance used for newton method */
41: PetscInt alloced; /* Size of workspace vectors allocated */
42: PetscInt init_alloc; /* Initial size for workspace vectors */
44: PetscInt max_lanczos_its; /* Maximum lanczos iterations */
45: PetscInt max_newton_its; /* Maximum newton iterations */
46: PetscInt dtype; /* Method used to measure the norm of step */
47: } KSP_GLTR;
49: #endif