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