conjugate_grads_method.m

MATLAB的梯度法源程序

function Conjugate_grads_method
%This function applies only the function like
%this:fx=X1^2+X2^2-X1*X2-10*X1-4*X2+60

syms X1 X2 real;
fx=X1^2+X2^2-X1*X2-10*X1-4*X2+60;
X0=[0;0];
e=0.0001;
k=0;
n=2;
t=1;
N=100;
S0_f=[diff(fx,X1,1);diff(fx,X2,1)];
value_S0_1=inline(vectorize(S0_f(1)),'X1','X2');
value_S0_2=inline(vectorize(S0_f(2)),'X1','X2');
S0=-[value_S0_1(X0(1),X0(2));value_S0_2(X0(1),X0(2))];
options=optimset('Display','off');
fx_value=inline(vectorize(fx),'X1','X2');

while sqrt(sum(S0.^2))>e&k<N
    t=t+1;
    k=k+1;
    if k==(n+1)
        k=0;
        S0_f=[diff(fx,X1,1);diff(fx,X2,1)];
        value_S0_1=inline(vectorize(S0_f(1)),'X1','X2');
        value_S0_2=inline(vectorize(S0_f(2)),'X1','X2');
        S0=-[value_S0_1(X0(1),X0(2));value_S0_2(X0(1),X0(2))];
    end
    syms r real;
    XX=X0+r.*S0;
    r_fx=fx_value(XX(1),XX(2));
    diff_r_fx=inline(vectorize(diff(r_fx,r,1)),'r');
    r=fzero(diff_r_fx,1,options);
    X=X0+r.*S0;
    grad_S1=[value_S0_1(X(1),X(2));value_S0_2(X(1),X(2))];
    bk=sum(grad_S1.^2)/sum(S0.^2);
    S0=-grad_S1+bk*S0;
    X0=X;
end
disp(k);
disp(X0);