Solving the problem for the free extremum of a function of several variables - numerical determination local minimum.

Hello, I am trying to do an automatization project, which numericaly determine and writes out the local minimum of given method, also it creates a table and a plot as written below. I am also uploading the .m-files which are also listen in one methon. Thank you for every help or solution in advance.

Task: Find the minimum of the objective function 𝒇(𝒙, 𝒚) = x^4 + 5x^2 - 9x^2*y + 2y^2 + 2y^4 + 2x

My given starting point for every method is: [0;0.5]

a) by the Newton and Raphson method

b) Levenberg and Marquardt method. Consider 𝛼 = 8, 𝑐 = 4

c) by the method of conjugate gradients: add function to the konj_tab program (also to graph) and find the minimum for entered starting point and for another selected starting point - at least one different (in each coordinate) from the specified one, I want this: [0.5;1]

d) compare in one image the dependence of the value of the objective function on the iterations (x-axis: number iterations, y-axis: value of the objective function)

For all methods, consider an individually specified starting point. Consider 𝑔1^2 + 𝑔2^2 = 𝑑 ≤ 0.001 as the termination condition of the minimization process, 𝑔1, 𝑔2 are gradient components of the gradient ∇𝑓 = [ g1 ; g2 ].

Solution instructions:

In each method:

- state the basic relationship and briefly describe the procedure, at least 2 steps.

- state the results of the individual calculation steps (in the table).

- draw a contour plot of the function in a range appropriate to the calculated points and plot the points calculated in individual iterations into it.

Evaluate and compare the obtained results. Consider success - whether the method found a local extremum (using the Hessian), accuracy, number of iterations, computational complexity.

[Sam Chak]

%% Main script
type konj_tab.m

disp('');disp('Konjugovane gradienty - Polak a Ribiere');format short e;%clc;
clear;close all
c=input('Vychodiskovy bod x0 = ');d=input('Presnost delta =');x(1,:)=c;o=1;m=0;hold on;

[fl,ql]=e9(c);
g(1,:)=ql;s(1,:)=g(1,:);
   while (g(o,:)*g(o,:)'>d)&(o<44),
     i=1; ml=0;
     jl=-s(o,:)*ql';S=-sign(jl);
     mp=S*d;
     e=c-S*d*s(o,:);
     [fp,q]=e9(e);
     jp=-s(o,:)*q';
       while jl*jp>0,
          jl=jp;fl=fp;ml=mp;
          mp=ml+S*d*2^i;
          e=c-mp*s(o,:);
          [fp,q]=e9(e);
          jp=-s(o,:)*q';
          i=i+1;
       end
       if S<0,
          jo=jp;jp=jl;jl=jo;
          fo=fp;fp=fl;fl=fo;
          mo=mp;mp=ml;ml=mo;
       end 
       j=min([jp,-jl]);
       while j*j>.01*d*d,
          z=jp+jl-3*(fp-fl)/(mp-ml);w=sqrt(z*z-jp*jl);
          m=ml+(mp-ml)*(z+w-jl)/(jp-jl+2*w);
          e=c-m*s(o,:);
          [f,q]=e9(e);j=-s(o,:)*q';
            if j>0,
               fp=f;jp=j;mp=m;
             else
               fl=f;jl=j;ml=m;
            end
        end
        x(o+1,:)=x(o,:)-m*s(o,:);o=o+1;c=x(o,:);
        [fl,ql]=e9(c);
        g(o,:)=ql;
        k(o-1)=g(o,:)*(g(o-1,:)-g(o,:))'/(g(o-1,:)*g(o-1,:)');
        s(o,:)=g(o,:)-k(o-1)*s(o-1,:);
        gm(o-1)=m;
    end;
% kreslenie obrazka 
[X,Y] = meshgrid(-1:.1:2, -1:.1:2);
        Z = X.^4+5*Y.^2-X.*Y-2*Y;
% vykreslenie vrstevnic do dalsieho obrazka
 v=-0.2:.1:20;vc=v.*v;
  contour(X,Y,Z,v)
hold on
title('Metóda konjug. gradientov')

roz=size(x);
a=x(:,1);b=x(:,2);
for i=1:roz(1)-1
P2=plot([a(i),a(i+1)],[b(i),b(i+1)],'g',[a(i),a(i+1)],[b(i),b(i+1)],'ok');
set(P2,'linewidth',2);
pause(0.3)
end
P2=plot(a,b,'g',a,b,'ok');
hold off

u=g(:,1);v=g(:,2);p=s(:,1);n=s(:,2);w=floor(o/22);
x1=a;x2=b;g1=u;g2=v;s1=p;s2=n;
Vysledky=table(x1,x2,g1,g2,s1,s2)
krokov=roz(1)-1

%% Function file
type e9.m

function [p,q]=e9(x)
% x1^2+x2^2-4x1-5x2-5
% p=x(1)^2/4+x(2)^2/9-0.8*x(1)-x(2)-0.3*x(1)*x(2)-3;
% q(1)=0.5*x(1)-.8-0.3*x(2);
% q(2)=2*x(2)/9-1-0.3*x(1);

%x^4+5*y^2-x*y-2*y
p=x(1)^4+5*x(2)^2-x(1)*x(2)-2*x(2);
q(1)=4*x(1)^3-x(2);
q(2)=10*x(2)-x(1)-2;

[Sam Chak]

The minimum point seems to lie at .

X   = linspace(-4, 4, 81);
Y   = linspace( 0, 3, 61);
[x, y]  = meshgrid(X, Y);
z   = x.^4 + 5*x.^2 - 9*(x.^2).*y + 2*y.^2 + 2*y.^4 + 2*x;

figure(1)
surfc(x, y, z),
xlabel('x'), ylabel('y'), zlabel('z')

figure(2)
contour(x, y, z, 50)
xlabel('x'), ylabel('y')

[DS]

Thank you, this is only c) by the method of conjugate gradients I guess, for me I got it completely different (see pictures - graph and table), the minimum was (−0.263, 0.149) , I am attaching the files too I can't guess whats wrong.

Answer 1

[Alan Stevens]

Here's a "starter for ten" with the Newton-Raphson method:

% Newton-Raphson Approach
 
% Function
f = @(x,y)x^4+5*x^2-9*x^2*y+2*y^2+2*y^4+2*x;

% Gradients
gx = @(x,y) 4*x^3+10*x-18*x*y+2;
gy = @(x,y) -9*x^2+4*y+8*y^3;

% Hessian
H = @(x,y) [12*x^2+10-18*y, -18*x;
            -18*x,          4+24*y^2];

% Initial guess
x0 = -3; y0 = 3;

tol = 1E-3; err = 1;
maxsteps = 20; steps = 0;

x = x0; y = y0;
while err>tol && steps<maxsteps
  steps = steps + 1;
  
  X = [x; y] - H(x,y)\[gx(x,y); gy(x,y)];
  err = gx(x,y)^2 + gy(x,y)^2;
  
  x = X(1); y = X(2); 
end

disp(['coordinates are: (', num2str(x),' ', num2str(y),')'])

coordinates are: (-2.3734 1.7606)

disp(['in ', int2str(steps), ' steps'])

in 5 steps

disp(['function value is: ', num2str(f(x,y))])

function value is: -8.6924

Note that the results are very sensitive to the initial guesses as is suggested by the contour plot posted by Sam Chak above.

[DS]

Hello thank you, I got this method also completely different I got values: −0.2638 and 0.1499 with 6 iterations I am attaching the code.

[Torsten]

@DS

Your code for Newton's method is correct. @Alan Stevens chose a different initial point.