CODEWITHSUNDEEP
matlabratenewtons-methodconvergence
Giuliana Paradiso
Giuliana Paradiso

Reputation: 1

Experimental order of convergence for Newton methods - MATLAB

I'm having problems with a piece of MATLAB code that is supposed to calculate the experimental order of convergence for the modified Newton method and the Truncated newton method.

x_min is the real minimum which was given, ek, ekm, ekp are the errors, pseq is the sequence of orders of convergence. Problem is, when testing the MN or TN on the Rosenbrock function it returns orders p very very large or very very small.

Here's the code:

typex_min=[1;1];
if k<=100
    pseq=zeros(k,1);
    for i=2:length(xseq)-1

        ekm=norm(xseq(:,i-1)-x_min);
        ek=norm(xseq(:,i)-x_min);
        ekp=norm(xseq(:,i+1)-x_min);
        pseq(i-1)= log( ekp/ek) ./ log( ek/ekm ) ;
    end
    p=mean(pseq)
else
    xseq_new=xseq(:,length(xseq)-100:length(xseq));
    pseq=zeros(100,1);
    for i=2:length(xseq_new)-1
        ekm=norm(xseq_new(:,i-1)-x_min);
        ek=norm(xseq_new(:,i)-x_min);
        ekp=norm(xseq_new(:,i+1)-x_min);
        pseq(i-1)= log( ekp/ek) ./ log( ek/ekm );
    end
    p=mean(pseq)
end here

Upvotes: 0

Views: 84

Answers (1)

lvand
lvand

Reputation: 365

From what I see, the issue you're encountering with the experimental order being too large or too small could be due to a couple of factors such as minor errors leading to numerical instability or incorrect indexing.

What I would recommend is:

  • Double-check the indexing in your loops is correct. The indexing for ekm, ek and ekp should be consistent with other iterations.
  • Ensure that the errors are not close to machine pression, as when errors are extremely small the logarithm operations can amplify differences leading to large or small values.
  • If your errors are too small consider setting a threshold below the convergence order; so that it isn't computed. You could also add a small epsilon to prevent division from small numbers.

I have altered some of your code to make it work.

% Assuming x_min, k, and xseq are already defined

if k <= 100
    pseq = zeros(k, 1);
    for i = 2:length(xseq)-1
        ekm = norm(xseq(:, i-1) - x_min);
        ek = norm(xseq(:, i) - x_min);
        ekp = norm(xseq(:, i+1) - x_min);
        
        % Check for small values to avoid numerical instability
        if ekm > 1e-12 && ek > 1e-12 && ekp > 1e-12
            pseq(i-1) = log(ekp / ek) / log(ek / ekm);
        else
            pseq(i-1) = NaN; % Ignore this calculation
        end
    end
    % Calculate the mean ignoring NaN values
    p = mean(pseq(~isnan(pseq)));
else
    xseq_new = xseq(:, end-99:end); % Correct the slicing
    pseq = zeros(100, 1);
    for i = 2:length(xseq_new)-1
        ekm = norm(xseq_new(:, i-1) - x_min);
        ek = norm(xseq_new(:, i) - x_min);
        ekp = norm(xseq_new(:, i+1) - x_min);
        
        % Check for small values to avoid numerical instability
        if ekm > 1e-12 && ek > 1e-12 && ekp > 1e-12
            pseq(i-1) = log(ekp / ek) / log(ek / ekm);
        else
            pseq(i-1) = NaN; % Ignore this calculation
        end
    end
    % Calculate the mean ignoring NaN values
    p = mean(pseq(~isnan(pseq)));
end

% Output the experimental order of convergence
disp(['Experimental order of convergence: ', num2str(p)]);

I have just made sure that MATLAB checks to ensure that the ekm, ek and ekp are within thresholds. I have also made sure to ignore NaN values when calculating the mean order of convergence to improve results.

Upvotes: 0

Related Questions