MATLAB 有约束信赖域算法，以四元多项式为算例
functionTRM4% All Rights Reserved lilongduzhi@yahoo.com.cn hehe,:)
% Optimal Function : min = -15;X = ;
% myfun = x1-x2-x3-x1*x3+x1*x4+x2*x3-x2*x4;
% st:x1 + 2x2 <= 8 ;
% 4x1 + x2<= 12 ;
% 3x1 + 4x2 <= 12 ;
% 2x3 - x4<= 8 ;
% x3 + 2x4<= 8 ;
% x3 + x4 <= 5 ;
% x1,x2,x3>= 0 ;
% diff_g(X) = [ 1-x3+x4 ,-1+x3-x4 ,-1-x1+x2 ,x1-x2];
nvars = 4; % 变量个数 ；
A = [ 1,2 ,0 ,0;
4,1 ,0 ,0;
3,4 ,0 ,0;
0,0 ,2 ,-1;
0,0 ,1 ,2;
0,0 ,1 ,1;];
A = [-A;eye(3),zeros(3,1)] ;
b = [ -8 ,-12 ,-12 ,-8 ,-8 ,-5 ,0 ,0 ,0]';
% step 1:
delta = 0.9; %
epsilon = 1e-6;
X = '; % X=;
B = eye(nvars);
k = 0;
%step 2:
grad = diff_g(X);
sqp_A = -A;
sqp_b = A*X-b;
upper = delta*X;
downer = -delta*X;
d = quadprog(B ,grad ,sqp_A ,sqp_b ,[],[] ,-delta*X ,delta*X);
while max(abs(d)) > epsilon
k = k+1 ;
% step 2:
norm_d = max(abs(d))
sqp_b =A*X-b;
lower = -abs(delta*X);
upper =abs(delta*X);
% d = quadprog(B ,grad ,sqp_A,sqp_b,[],[] ,-delta*X ,delta*X);
d = quadprog(B ,grad ,sqp_A,sqp_b,[],[] ,lower ,upper);
% step 2.2:
grad = diff_g(X);
ared = myfun(X)-myfun(X+d);
pred = q(zeros(nvars,1) ,grad ,B)-q(d ,grad ,B);
r = ared/pred;
% step 2.1:
if r > 0
X = X+d;
else
X = X;
end
% step 3:
if r>=0.25 ,
% goto step 4:
else delta = delta/2;
end
if r<0.75 || max(abs(d)) < delta ,
% goto step 5:
else
delta = min(2*delta ,1-1e-10);
end
% step 5:
delta = delta ;
% 5.1 : mod B
B = BFGS(B ,X ,d);
k = k+1
end
disp('求解结果：');
X
min_value = myfun(X)
% =====================> sub function <=====================%
function y = myfun(X)
% myfun = x1-x2-x3-x1*x3+x1*x4+x2*x3-x2*x4;
y = X(1)-X(2)-X(3)-X(1)*X(3)+X(1)*X(4)+ X(2)*X(3)-X(2)*X(4);
function grad = diff_g(X)
% diff_g(X) = [ 1-x3+x4 ,-1+x3-x4 ,-1-x1+x2 ,x1-x2];
grad = [ 1-X(3)+X(4) ,-1+X(3)-X(4) ,-1-X(1)+X(2) ,X(1)-X(2)]';
function Q = q(d ,g ,B)
Q = g'*d+1/2*(d'*B*d);
function B = BFGS(B ,X,d)
y = diff_g(X+d)-diff_g(X);
s = d;
if y'*s > 0
B = B - B*s*s'*B/(s'*B*s)+y*y'/(y'*s);
end
%=============================================================%
PS: 程序仍有一定缺陷，算法依赖迭代初值现象严重。
欢迎大家为我指正，谢谢！
页:
[1]