1*a4963045SJacob Faibussowitsch #pragma once 2af0996ceSBarry Smith #include <petsc/private/taoimpl.h> 3a7e14dcfSSatish Balay 4a7e14dcfSSatish Balay /* 5a7e14dcfSSatish Balay Context for Interior-Point Method 6a7e14dcfSSatish Balay */ 7a7e14dcfSSatish Balay 8a7e14dcfSSatish Balay typedef struct { 9a7e14dcfSSatish Balay PetscInt mi, me, n, nxb, nib, nb, nslack; 10a7e14dcfSSatish Balay PetscInt nuser_inequalities; 11a7e14dcfSSatish Balay PetscInt nxlb, nxub, niub, nilb; 12a7e14dcfSSatish Balay PetscScalar sig, mu, taumin, dec; 13a7e14dcfSSatish Balay PetscScalar muaff; 14a7e14dcfSSatish Balay TaoLineSearch lag_ls; 15a7e14dcfSSatish Balay Vec work, rhs_x, save_x; 16a8d3b578SPierre Jolivet Vec lambdai, dlambdai, rhs_lambdai, save_lambdai; 17a8d3b578SPierre Jolivet Vec lambdae, dlambdae, rhs_lambdae, save_lambdae; 18a7e14dcfSSatish Balay Vec s, ds, rhs_s, save_s; 19a7e14dcfSSatish Balay Vec ci; 20a7e14dcfSSatish Balay Vec Zero_nb, One_nb, Inf_nb; 21a7e14dcfSSatish Balay PetscScalar kkt_f; /* d'*x + (1/2)*x'*H*x; */ 22a8d3b578SPierre Jolivet Vec rd; /* H*x + d + Ae'*lambdae - Ai'*lambdai */ 23a7e14dcfSSatish Balay Vec rpe; /* residual Ae*x - be */ 24a7e14dcfSSatish Balay Vec rpi; /* Ai*x - yi - bi */ 25a8d3b578SPierre Jolivet Vec complementarity; /* yi.*lambdai */ 26a7e14dcfSSatish Balay PetscScalar phi; 27a8d3b578SPierre Jolivet Mat L; /* diag(lambdai) */ 28a7e14dcfSSatish Balay Mat Y; /* diag(yi) */ 29a7e14dcfSSatish Balay Mat Ai; /* JacI (lb) 30a7e14dcfSSatish Balay -JacI (ub) 31a7e14dcfSSatish Balay I (xlb) 32a7e14dcfSSatish Balay -I (xub) */ 33a7e14dcfSSatish Balay Mat K; /* [ H , 0, Ae',-Ai']; 34a7e14dcfSSatish Balay [Ae , 0, 0 , 0]; 35a7e14dcfSSatish Balay [Ai ,-Imi, 0 , 0]; 36a7e14dcfSSatish Balay [ 0 , L , 0 , Y ]; */ 37a7e14dcfSSatish Balay 38a8d3b578SPierre Jolivet Vec bigrhs; /* rhs [x; lambdae; yi; lambdai] */ 39a8d3b578SPierre Jolivet Vec bigstep; /* [dx; dyi; dlambdae; dlambdai] */ 40a7e14dcfSSatish Balay PetscBool monitorkkt; 41a7e14dcfSSatish Balay PetscScalar alpha1, alpha2; 42a7e14dcfSSatish Balay PetscScalar pushs, pushnu; 43a7e14dcfSSatish Balay IS isxl, isxu, isil, isiu; 44a7e14dcfSSatish Balay VecScatter ci_scat, xl_scat, xu_scat; 45a7e14dcfSSatish Balay VecScatter step1, step2, step3, step4; 46a7e14dcfSSatish Balay VecScatter rhs1, rhs2, rhs3, rhs4; 47a7e14dcfSSatish Balay } TAO_IPM; 48