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