Actual source code: nlsimpl.h
1: /*
2: Context for a Newton line search method (unconstrained minimization)
3: */
5: #pragma once
6: #include <petsc/private/taoimpl.h>
8: typedef struct {
9: Mat M;
10: PC bfgs_pre;
12: Vec D;
13: Vec W;
15: Vec Xold;
16: Vec Gold;
18: /* Parameters when updating the perturbation added to the Hessian matrix
19: according to the following scheme:
21: pert = sval;
23: do until convergence
24: shift Hessian by pert
25: solve Newton system
27: if (linear solver failed or did not compute a descent direction)
28: use steepest descent direction and increase perturbation
30: if (0 == pert)
31: initialize perturbation
32: pert = min(imax, max(imin, imfac * norm(G)))
33: else
34: increase perturbation
35: pert = min(pmax, max(pgfac * pert, pmgfac * norm(G)))
36: fi
37: else
38: use linear solver direction and decrease perturbation
40: pert = min(psfac * pert, pmsfac * norm(G))
41: if (pert < pmin)
42: pert = 0
43: fi
44: fi
46: perform line search
47: function and gradient evaluation
48: check convergence
49: od
50: */
51: PetscReal sval; /* Starting perturbation value, default zero */
53: PetscReal imin; /* Minimum perturbation added during initialization */
54: PetscReal imax; /* Maximum perturbation added during initialization */
55: PetscReal imfac; /* Merit function factor during initialization */
57: PetscReal pmin; /* Minimum perturbation value */
58: PetscReal pmax; /* Maximum perturbation value */
59: PetscReal pgfac; /* Perturbation growth factor */
60: PetscReal psfac; /* Perturbation shrink factor */
61: PetscReal pmgfac; /* Merit function growth factor */
62: PetscReal pmsfac; /* Merit function shrink factor */
64: /* Parameters when updating the trust-region radius based on steplength
65: if step < nu1 (very bad step)
66: radius = omega1 * min(norm(d), radius)
67: elif step < nu2 (bad step)
68: radius = omega2 * min(norm(d), radius)
69: elif step < nu3 (okay step)
70: radius = omega3 * radius;
71: elif step < nu4 (good step)
72: radius = max(omega4 * norm(d), radius)
73: else (very good step)
74: radius = max(omega5 * norm(d), radius)
75: fi
76: */
77: PetscReal nu1; /* used to compute trust-region radius */
78: PetscReal nu2; /* used to compute trust-region radius */
79: PetscReal nu3; /* used to compute trust-region radius */
80: PetscReal nu4; /* used to compute trust-region radius */
82: PetscReal omega1; /* factor used for trust-region update */
83: PetscReal omega2; /* factor used for trust-region update */
84: PetscReal omega3; /* factor used for trust-region update */
85: PetscReal omega4; /* factor used for trust-region update */
86: PetscReal omega5; /* factor used for trust-region update */
88: /* Parameters when updating the trust-region radius based on reduction
90: kappa = ared / pred
91: if kappa < eta1 (very bad step)
92: radius = alpha1 * min(norm(d), radius)
93: elif kappa < eta2 (bad step)
94: radius = alpha2 * min(norm(d), radius)
95: elif kappa < eta3 (okay step)
96: radius = alpha3 * radius;
97: elif kappa < eta4 (good step)
98: radius = max(alpha4 * norm(d), radius)
99: else (very good step)
100: radius = max(alpha5 * norm(d), radius)
101: fi
102: */
103: PetscReal eta1; /* used to compute trust-region radius */
104: PetscReal eta2; /* used to compute trust-region radius */
105: PetscReal eta3; /* used to compute trust-region radius */
106: PetscReal eta4; /* used to compute trust-region radius */
108: PetscReal alpha1; /* factor used for trust-region update */
109: PetscReal alpha2; /* factor used for trust-region update */
110: PetscReal alpha3; /* factor used for trust-region update */
111: PetscReal alpha4; /* factor used for trust-region update */
112: PetscReal alpha5; /* factor used for trust-region update */
114: /* Parameters when updating the trust-region radius based on interpolation
116: kappa = ared / pred
117: if kappa >= 1.0 - mu1 (very good step)
118: choose tau in [gamma3, gamma4]
119: radius = max(tau * norm(d), radius)
120: elif kappa >= 1.0 - mu2 (good step)
121: choose tau in [gamma2, gamma3]
122: if (tau >= 1.0)
123: radius = max(tau * norm(d), radius)
124: else
125: radius = tau * min(norm(d), radius)
126: fi
127: else (bad step)
128: choose tau in [gamma1, 1.0]
129: radius = tau * min(norm(d), radius)
130: fi
131: */
132: PetscReal mu1; /* used for model agreement in interpolation */
133: PetscReal mu2; /* used for model agreement in interpolation */
135: PetscReal gamma1; /* factor used for interpolation */
136: PetscReal gamma2; /* factor used for interpolation */
137: PetscReal gamma3; /* factor used for interpolation */
138: PetscReal gamma4; /* factor used for interpolation */
140: PetscReal theta; /* factor used for interpolation */
142: /* Parameters when initializing trust-region radius based on interpolation */
143: PetscReal mu1_i; /* used for model agreement in interpolation */
144: PetscReal mu2_i; /* used for model agreement in interpolation */
146: PetscReal gamma1_i; /* factor used for interpolation */
147: PetscReal gamma2_i; /* factor used for interpolation */
148: PetscReal gamma3_i; /* factor used for interpolation */
149: PetscReal gamma4_i; /* factor used for interpolation */
151: PetscReal theta_i; /* factor used for interpolation */
153: /* Other parameters */
154: PetscReal min_radius; /* lower bound on initial radius value */
155: PetscReal max_radius; /* upper bound on trust region radius */
156: PetscReal epsilon; /* tolerance used when computing ared/pred */
158: PetscInt newt; /* Newton directions attempted */
159: PetscInt bfgs; /* BFGS directions attempted */
160: PetscInt grad; /* Gradient directions attempted */
162: PetscInt init_type; /* Trust-region initialization method */
163: PetscInt update_type; /* Trust-region update method */
165: PetscInt ksp_atol;
166: PetscInt ksp_rtol;
167: PetscInt ksp_ctol;
168: PetscInt ksp_negc;
169: PetscInt ksp_dtol;
170: PetscInt ksp_iter;
171: PetscInt ksp_othr;
172: } TAO_NLS;