Gradient projection method.
45 views (last 30 days)
Show older comments
Hi everyone,
Does somebody implemented the gradient projection method?
I have difficulties to define a constrained set in matlab (where I have to project to).
The constraint that I wanted to implement is D/A <= 1.
Thank u very much in advance.
1 Comment
AMEHRI WALID
on 28 Jun 2019
Edited: AMEHRI WALID
on 4 Sep 2019
Hi,
I have developped a code for the problem shown in this video https://www.youtube.com/watch?v=1mCyC1FyMhk.
The problem is:
minimize f(x) NL
subject to A*X <= B
Here is the code:
% Example 1 :
% minimize: f(X) = x1^2 + x2^2 + x3^2 + x4^2 - 2x^1 - 3x4
% g(X) = AX <= B
% Where A = [Matrix] B = [7; 6; 0; 0; 0; 0]
% and all X >= 0, size(X) = (4, 1)
clear;
close;
clc;
%% Initial State : Choose an initial point that satisfy all the constraints
X0 = zeros(4, 1);
x10 = X0(1);
x20 = X0(2);
x30 = X0(3);
x40 = X0(4);
% Evaluate the Objective function, the Gradient of the Objective function and the constraint function at X = X0
ObjFun_X0 = x10^2 + x20^2 + x30^2 + x40^2 - 2*x10 - 3*x40;
Minus_Grad_ObjFun_X0 = [2 - 2*x10; -2*x20; -2*x30; 3-2*x40];
d_X0 = Minus_Grad_ObjFun_X0;
B = [7; 6; 0; 0; 0; 0];
A = [2 1 1 4; % A is the gradient of g(X)
1 1 2 1;
-1 0 0 0;
0 -1 0 0;
0 0 -1 0;
0 0 0 -1];
g_X0 = A*X0 - B; % the constraints equations are described by g(X) = AX <= B
% Find LC and TC
Test_C = A * X0;
index_AT = 0;
index_AL = 0;
LC = 0;
TC = 0;
for i = 1:6
% Tight Constraints
if ( Test_C(i) >= B(i) )
index_AT = index_AT + 1;
TC(index_AT) = i;
% Loose Constraints
else
index_AL = index_AL + 1;
LC(index_AL) = i;
end
end
TC = TC';
LC = LC';
% Get AL and AT
AL = zeros(size(LC, 1), 4);
AT = zeros(size(TC, 1), 4);
for j = 1:size(LC, 1)
AL (j, :) = A(LC(j), :);
end
for k = 1:size(TC, 1)
AT (k, :) = A(TC(k), :);
end
% find the Direction vector
Test_G = AT * d_X0;
dir_X0 = d_X0;
for e = 1:size(Test_G, 1)
if ( Test_G (e) > 0 )
Projection = eye(4) - AT' * (inv(AT * AT')) * AT;
dir_X0 = Projection * d_X0;
break;
end
end
% Compute the Langrange Multipliers for the tight Constraints : They must all postives > 0 to satify the optimality of the solution
LM = inv(AT * AT') * AT * d_X0;
count_LM = 0;
for i = 1:size(LM, 1)
if ( LM(i) <= 0 )
LM_Condition = 0;
count_LM = 0;
break;
end
if ( LM(i) > 0)
count_LM = count_LM + 1;
if ( count_LM >= size(LM, 1))
LM_Condition = 1;
end
end
end
count = 1;
Tolerance = 1e-8;
All_X(:, count) = X0;
%% ***** Begin the Main Loop of the Gradient Prjection method ***** %%
while ( norm(dir_X0) > Tolerance || LM_Condition~= 1 )
count = count + 1;
% Compute Lamda using the Table
All_Lamda = zeros(size(LC, 1), 1);
for i = 1 : size(LC, 1)
index = LC(i);
All_Lamda(i) = -g_X0(index) / ( AL(i, :) * dir_X0 );
end
All_Lamda_corrected = All_Lamda(All_Lamda >= 0);
Lamda = min(All_Lamda_corrected);
% Compute the New X
X = X0 + Lamda * dir_X0;
X0 = X;
% Store the solution at each step
All_X(:, count) = X0;
% Evaluate the ObjFun, Grad of the ObjFun and the Constraint Func
x10 = X0(1);
x20 = X0(2);
x30 = X0(3);
x40 = X0(4);
ObjFun_X0 = x10^2 + x20^2 + x30^2 + x40^2 - 2*x10 - 3*x40;
Minus_Grad_ObjFun_X0 = [2 - 2*x10; -2*x20; -2*x30; 3-2*x40];
d_X0 = Minus_Grad_ObjFun_X0;
g_X0 = A*X0 - B;
% Find LC and TC, we always compare to initial original constraints
Test_C = A * X0;
index_AT = 0;
index_AL = 0;
LC = 0;
TC = 0;
for i = 1:6
% Tight Constraints
if ( round(Test_C(i), 4) >= B(i) )
index_AT = index_AT + 1;
TC(index_AT) = i;
% Loose Constraints
else
index_AL = index_AL + 1;
LC(index_AL) = i;
end
end
TC = TC';
LC = LC';
% Get AL and AT
AL = zeros(size(LC, 1), 4);
AT = zeros(size(TC, 1), 4);
for j = 1:size(LC, 1)
AL (j, :) = A(LC(j), :);
end
for k = 1:size(TC, 1)
AT (k, :) = A(TC(k), :);
end
% Test if we need to project the gradient or not
Test_G = AT * d_X0;
dir_X0 = d_X0;
for e = 1:size(Test_G, 1)
if ( Test_G (e) > 0 )
Projection = eye(4) - AT' * (inv(AT * AT')) * AT;
dir_X0 = Projection * d_X0;
break;
end
end
% Compute the Langrange Multipliers for the tight Constraints
LM = inv(AT * AT') * AT * d_X0;
count_LM = 0;
for i = 1:size(LM, 1)
if ( LM(i) <= 0 )
LM_Condition = 0;
count_LM = 0;
break;
end
if ( LM(i) > 0)
count_LM = count_LM + 1;
if ( count_LM >= size(LM, 1))
LM_Condition = 1;
end
end
end
end % End While
Answers (1)
See Also
Categories
Find more on Nonlinear Optimization in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!