Actual source code: pdipm.h
2: #ifndef __TAO_PDIPM_H
4: #include <petsc/private/taoimpl.h>
6: /*
7: Context for Primal-Dual Interior-Point Method
8: See the document pdipm.pdf
9: */
11: typedef struct {
12: /* Sizes (n = local, N = global) */
13: PetscInt nx, Nx; /* Decision variables nx = nxfixed + nxub + nxlb + nxbox + nxfree */
14: PetscInt nxfixed, Nxfixed; /* Fixed decision variables */
15: PetscInt nxlb, Nxlb; /* Decision variables with lower bounds only */
16: PetscInt nxub, Nxub; /* Decision variables with upper bounds only */
17: PetscInt nxbox, Nxbox; /* Decision variables with box constraints */
18: PetscInt nxfree, Nxfree; /* Free variables */
19: PetscInt ng, Ng; /* user equality constraints g(x) = 0. */
20: PetscInt nh, Nh; /* user inequality constraints h(x) >= 0. */
21: PetscInt nce, Nce; /* total equality constraints. nce = ng + nxfixed */
22: PetscInt nci, Nci; /* total inequality constraints nci = nh + nxlb + nxub + 2*nxbox */
23: PetscInt n, N; /* Big KKT system size n = nx + nce + 2*nci */
25: /* Vectors */
26: Vec X; /* R^n - Big KKT system vector [x; lambdae; lambdai; z] */
27: Vec x; /* R^nx - work vector, same layout as tao->solution */
28: Vec lambdae; /* R^nce - vector, shares local arrays with X */
29: Vec lambdai; /* R^nci - vector, shares local arrays with X */
30: Vec z; /* R^nci - vector, shares local arrays with X */
32: /* Work vectors */
33: Vec lambdae_xfixed; /* Equality constraints lagrangian multiplier vector for fixed variables */
34: Vec lambdai_xb; /* User inequality constraints lagrangian multiplier vector */
36: /* Lagrangian equality and inequality Vec */
37: Vec ce, ci; /* equality and inequality constraints */
39: /* Offsets for subvectors */
40: PetscInt off_lambdae, off_lambdai, off_z;
42: /* Scalars */
43: PetscReal L; /* Lagrangian = f(x) - lambdae^T*ce(x) - lambdai^T*(ci(x) - z) - mu*sum_{i=1}^{Nci}(log(z_i)) */
44: PetscReal gradL; /* gradient of L w.r.t. x */
46: /* Matrices */
47: Mat Jce_xfixed; /* Jacobian of equality constraints cebound(x) = J(nxfixed) */
48: Mat Jci_xb; /* Jacobian of inequality constraints Jci = [tao->jacobian_inequality ; J(nxub); J(nxlb); J(nxbx)] */
49: Mat K; /* KKT matrix */
51: /* Parameters */
52: PetscReal mu; /* Barrier parameter */
53: PetscReal mu_update_factor; /* Multiplier for mu update */
54: PetscReal deltaw;
55: PetscReal lastdeltaw;
56: PetscReal deltac;
58: /* Tolerances */
60: /* Index sets for types of bounds on variables */
61: IS isxub; /* Finite upper bound only -inf < x < ub */
62: IS isxlb; /* Finite lower bound only lb <= x < inf */
63: IS isxfixed; /* Fixed variables lb = x = ub */
64: IS isxbox; /* Boxed variables lb <= x <= ub */
65: IS isxfree; /* Free variables -inf <= x <= inf */
67: /* Index sets for PC fieldsplit */
68: IS is1, is2;
70: /* Options */
71: PetscBool monitorkkt; /* Monitor KKT */
72: PetscReal push_init_slack; /* Push initial slack variables (z) away from bounds */
73: PetscReal push_init_lambdai; /* Push initial inequality variables (lambdai) away from bounds */
74: PetscBool solve_reduced_kkt; /* Solve Reduced KKT with fieldsplit */
75: PetscBool solve_symmetric_kkt; /* Solve non-reduced symmetric KKT system */
76: PetscBool kkt_pd; /* Add deltaw and deltac shifts to make KKT matrix positive definite */
78: SNES snes; /* Nonlinear solver */
79: Mat jac_equality_trans, jac_inequality_trans; /* working matrices */
81: PetscReal obj; /* Objective function */
83: /* Offsets for parallel assembly */
84: PetscInt *nce_all;
85: } TAO_PDIPM;
87: #endif /* ifndef __TAO_PDIPM_H */