Actual source code: bncg.h
1: /*
2: Context for bound-constrained nonlinear conjugate gradient method
3: */
5: #ifndef __TAO_BNCG_H
8: #include <petsc/private/taoimpl.h>
10: typedef struct {
11: Mat B;
12: Mat pc;
13: Vec G_old, X_old, W, work;
14: Vec g_work, y_work, d_work;
15: Vec sk, yk;
16: Vec unprojected_gradient, unprojected_gradient_old;
17: Vec inactive_grad, inactive_step;
19: IS active_lower, active_upper, active_fixed, active_idx, inactive_idx, inactive_old, new_inactives;
21: PetscReal alpha; /* convex ratio in the scalar scaling */
22: PetscReal hz_theta;
23: PetscReal xi; /* Parameter for KD, DK, and HZ methods. */
24: PetscReal theta; /* The convex combination parameter in the SSML_Broyden method. */
25: PetscReal zeta; /* Another parameter, exclusive to Kou-Dai, modifying the y_k scalar contribution */
26: PetscReal hz_eta, dk_eta;
27: PetscReal bfgs_scale, dfp_scale; /* Scaling of the bfgs/dfp tau parameter in the bfgs and broyden methods. Default 1. */
28: PetscReal tau_bfgs, tau_dfp;
29: PetscReal as_step, as_tol, yts, yty, sts;
30: PetscReal f, delta_min, delta_max;
31: PetscReal epsilon; /* Machine precision unless changed by user */
32: PetscReal eps_23; /* Two-thirds power of machine precision */
34: PetscInt cg_type; /* Formula to use */
35: PetscInt min_restart_num; /* Restarts every x*n iterations, where n is the dimension */
36: PetscInt ls_fails, resets, descent_error, skipped_updates, pure_gd_steps;
37: PetscInt iter_quad, min_quad; /* Dynamic restart variables in Dai-Kou, SIAM J. Optim. Vol 23, pp. 296-320, Algorithm 4.1 */
38: PetscInt as_type;
40: PetscBool inv_sig;
41: PetscReal tol_quad; /* tolerance for Dai-Kou dynamic restart */
42: PetscBool dynamic_restart; /* Keeps track of whether or not to do a dynamic (KD) restart */
43: PetscBool spaced_restart; /* If true, restarts the CG method every x*n iterations */
44: PetscBool use_dynamic_restart;
45: PetscBool neg_xi;
46: PetscBool unscaled_restart; /* Gradient descent restarts are done without rescaling*/
47: PetscBool diag_scaling;
48: PetscBool no_scaling;
50: } TAO_BNCG;
52: #endif /* ifndef __TAO_BNCG_H */
54: PETSC_INTERN PetscErrorCode TaoBNCGEstimateActiveSet(Tao, PetscInt);
55: PETSC_INTERN PetscErrorCode TaoBNCGBoundStep(Tao, PetscInt, Vec);
56: PETSC_INTERN PetscErrorCode TaoBNCGComputeScalarScaling(PetscReal, PetscReal, PetscReal, PetscReal*, PetscReal);
57: PETSC_INTERN PetscErrorCode TaoBNCGConductIteration(Tao, PetscReal);
58: PETSC_INTERN PetscErrorCode TaoBNCGStepDirectionUpdate(Tao, PetscReal, PetscReal, PetscReal, PetscReal, PetscReal, PetscBool);
59: PETSC_INTERN PetscErrorCode TaoBNCGComputeDiagScaling(Tao, PetscReal, PetscReal);
60: PETSC_INTERN PetscErrorCode TaoBNCGResetUpdate(Tao, PetscReal);
61: PETSC_INTERN PetscErrorCode TaoBNCGCheckDynamicRestart(Tao, PetscReal, PetscReal, PetscReal, PetscBool*, PetscReal);