1% Copyright (C) 2009 International Business Machines 2% All Rights Reserved. 3% This code is published under the Eclipse Public License. 4% 5% $Id: hs071_c.c 699 2006-04-05 21:05:18Z andreasw $ 6% 7% Author: Andreas Waechter IBM 2009-04-02 8% 9% This file is part of the Ipopt tutorial. It is a version with 10% mistakes for the matlab implemention of the coding exercise problem 11% (in AMPL formulation): 12% 13% param n := 4; 14% 15% var x {1..n} <= 0, >= -1.5, := -0.5; 16% 17% minimize obj: 18% sum{i in 1..n} (x[i]-1)^2; 19% 20% subject to constr {i in 2..n-1}: 21% (x[i]^2+1.5*x[i]-i/n)*cos(x[i+1]) - x[i-1] = 0; 22% 23% The constant term "i/n" in the constraint is supposed to be input data 24 25function [x, info] = TutorialMatlab 26 27 % Size of the problem 28 n = 5; 29 30 % Problem data 31 a = ((1:n)/n)'; 32 33 % Starting point 34 x0 = -0.5*ones(n,1); 35 36 % Lower and upper bounds for the variables 37 options.lb = -1.5*ones(n,1); 38 options.ub = zeros(n,1); 39 40 % Constraint bounds 41 options.cl = zeros(n-2,1); 42 options.cu = zeros(n-2,1); 43 44 % Set the Ipopt options 45 % options.ipopt.mu_strategy = 'adaptive'; 46 options.ipopt.derivative_test = 'second-order'; 47 options.ipopt.max_iter = 10; 48 49 % Set the callback functions 50 funcs.objective = @eval_f; 51 funcs.constraints = @eval_g; 52 funcs.gradient = @eval_grad_f; 53 funcs.jacobian = @eval_jac_g; 54 funcs.jacobianstructure = @eval_jac_g_struct; 55 funcs.hessian = @eval_hess; 56 funcs.hessianstructure = @eval_hess_struct; 57 58 [x info] = ipopt(x0, funcs, options); 59 60 61 % End of main function 62 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 63 64 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 65 % Evaluate value of objective function 66 function f = eval_f(x) 67 68 tmp = (x - 1).^2; 69 f = sum(tmp); 70 71 end 72 73 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 74 % Evaluate gradient of objective function 75 function df = eval_grad_f(x) 76 77 df = (x - 1); 78 79 end 80 81 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 82 % Evaluate value of constraint bodies 83 function g = eval_g(x) 84 85 g = (x(2:n-1).^2 + 1.5*x(1:n-2) - a(2:n-1)).*cos(x(3:n)) - x(1:n-2); 86 87 end 88 89 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 90 % Return constraint Jacobian strcture 91 function A = eval_jac_g_struct 92 93 % tri-diagonal structure 94 Diags = ones(n, 3); 95 A = spdiags(Diags, 0:2, n-2, n); 96 97 end 98 99 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 100 % Evaluate constraint Jacobian 101 function A = eval_jac_g(x) 102 103 Diags = -1*ones(n-2,3); 104 Diags(:,2) = (2*x(2:n-1)+1.5).*cos(x(3:n)); 105 Diags(:,3) = -(x(2:n-1).^2 + 1.5*x(2:n-1) - a(2:n-1)).*sin(x(3:n)); 106 A = spdiags(Diags, 0:2, n-2, n); 107 108 end 109 110 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 111 % Return Hessian of Lagrangian function structure 112 function H = eval_hess_struct 113 114 Diags = ones(n, 2); 115 H = spdiags(Diags, -1:0, n, n); 116 117 end 118 119 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 120 % Evaluate Hessian of Lagrangian function 121 function H = eval_hess(x, sigma, lambda) 122 123 Diags = zeros(n,2); 124 125 % part from the objective function 126 Diags(:,2) = sigma * 2; 127 128 % (x_i , x_i) part 129 Diags(2:n-1,2) = Diags(2:n-1,2) + lambda.*cos(x(3:n)); 130 131 % (x_{i+1}, x_{i+1}) part 132 Diags(3:n,2) = Diags(3:n,2) ... 133 - lambda.*(x(2:n-1).^2 + 1.5*x(1:n-2) - a(2:n-1)).*cos(x(3:n)); 134 135 % (x_i, x_{i+1}) part 136 Diags(2:n-1,1) = -lambda.*sin(x(3:n)).*(2*x(2:n-1)+1.5); 137 138 H = spdiags(Diags, -1:0, n, n); 139 140 end 141 142end 143