1a4963045SJacob Faibussowitsch #pragma once 2a7e14dcfSSatish Balay 3af0996ceSBarry Smith #include <petsc/private/taoimpl.h> 4aaa7dc30SBarry Smith #include <petscis.h> 5a7e14dcfSSatish Balay #define LCL_FORWARD1 0 6a7e14dcfSSatish Balay #define LCL_ADJOINT1 1 7a7e14dcfSSatish Balay #define LCL_FORWARD2 2 8a7e14dcfSSatish Balay #define LCL_ADJOINT2 3 9a7e14dcfSSatish Balay 10a7e14dcfSSatish Balay typedef struct { 11a7e14dcfSSatish Balay Mat M; /* Quasi-newton hessian matrix */ 12a7e14dcfSSatish Balay Vec dbar; /* Reduced gradient */ 13a7e14dcfSSatish Balay Vec GL; 14a7e14dcfSSatish Balay Vec GAugL; 15a7e14dcfSSatish Balay Vec GL_U; /* Gradient of lagrangian */ 16a7e14dcfSSatish Balay Vec GL_V; /* Gradient of lagrangian */ 17a7e14dcfSSatish Balay Vec GAugL_U; /* Augmented lagrangian gradient */ 18a7e14dcfSSatish Balay Vec GAugL_V; /* Augmented lagrangian gradient */ 19a7e14dcfSSatish Balay Vec GL_U0; /* Gradient of lagrangian */ 20a7e14dcfSSatish Balay Vec GL_V0; /* Gradient of lagrangian */ 21a7e14dcfSSatish Balay Vec GAugL_U0; /* Augmented lagrangian gradient */ 22a7e14dcfSSatish Balay Vec GAugL_V0; /* Augmented lagrangian gradient */ 23a7e14dcfSSatish Balay 24a7e14dcfSSatish Balay IS UIS; /* Index set to state */ 25a7e14dcfSSatish Balay IS UID; /* Index set to design */ 26a7e14dcfSSatish Balay IS UIM; /* Full index set to all constraints */ 27a7e14dcfSSatish Balay VecScatter state_scatter; 28a7e14dcfSSatish Balay VecScatter design_scatter; 29a7e14dcfSSatish Balay 30a7e14dcfSSatish Balay Vec U; /* State variable */ 31a7e14dcfSSatish Balay Vec V; /* Design variable */ 32a7e14dcfSSatish Balay Vec U0; /* State variable */ 33a7e14dcfSSatish Balay Vec V0; /* Design variable */ 34a7e14dcfSSatish Balay Vec V1; /* Design variable */ 35a7e14dcfSSatish Balay 36a7e14dcfSSatish Balay Vec DU; /* State step */ 37a7e14dcfSSatish Balay Vec DV; /* Design step */ 38a7e14dcfSSatish Balay Vec DL; /* Multipliers step */ 39a7e14dcfSSatish Balay 40a7e14dcfSSatish Balay Vec GU; /* Gradient wrt U */ 41a7e14dcfSSatish Balay Vec GV; /* Gradient wrt V */ 42a7e14dcfSSatish Balay Vec GU0; /* Gradient wrt U */ 43a7e14dcfSSatish Balay Vec GV0; /* Gradient wrt V */ 44a7e14dcfSSatish Balay 45a7e14dcfSSatish Balay Vec W; /* work vector */ 46a7e14dcfSSatish Balay Vec X0; 47a7e14dcfSSatish Balay Vec G0; 48a7e14dcfSSatish Balay Vec WU; /* state work vector */ 49a7e14dcfSSatish Balay Vec WV; /* design work vector */ 50a7e14dcfSSatish Balay Vec r; 51a7e14dcfSSatish Balay Vec s; 52a7e14dcfSSatish Balay Vec g1, g2; 53a7e14dcfSSatish Balay Vec con1; 54a7e14dcfSSatish Balay 55a7e14dcfSSatish Balay PetscInt m; /* number of constraints */ 56a7e14dcfSSatish Balay PetscInt n; /* number of variables */ 57a7e14dcfSSatish Balay 58a7e14dcfSSatish Balay Mat jacobian_state0; /* Jacobian wrt U */ 59*7addb90fSBarry Smith Mat jacobian_state0_pre; /* matrix wrt U used to compute the preconditioner */ 60a7e14dcfSSatish Balay Mat jacobian_design0; /* Jacobian wrt V */ 61a7e14dcfSSatish Balay Mat jacobian_state_inv0; /* Inverse of Jacobian wrt U */ 62a7e14dcfSSatish Balay Mat R; 63a7e14dcfSSatish Balay 64a8d3b578SPierre Jolivet Vec lambda; /* Lagrange Multiplier */ 65a8d3b578SPierre Jolivet Vec lambda0; /* Lagrange Multiplier */ 66a8d3b578SPierre Jolivet Vec lambda1; /* Lagrange Multiplier */ 67a7e14dcfSSatish Balay 68a7e14dcfSSatish Balay Vec WL; /* Work vector */ 69a7e14dcfSSatish Balay PetscReal rho; /* Penalty parameter */ 70a7e14dcfSSatish Balay PetscReal rho0; 71a7e14dcfSSatish Balay PetscReal rhomax; 72a7e14dcfSSatish Balay PetscReal eps1, eps2; 73a7e14dcfSSatish Balay PetscReal aug, aug0, lgn, lgn0; 74a7e14dcfSSatish Balay PetscInt subset_type; 75a7e14dcfSSatish Balay PetscInt solve_type; 76a7e14dcfSSatish Balay PetscBool recompute_jacobian_flag; 77a7e14dcfSSatish Balay PetscInt phase2_niter; 78a7e14dcfSSatish Balay PetscBool verbose; 79a7e14dcfSSatish Balay PetscReal tau[4]; 80a7e14dcfSSatish Balay 81a7e14dcfSSatish Balay } TAO_LCL; 82