Same code but different results

4 views (last 30 days)
Using the same code, I get different results than those obtained on a different pc and version of Matlab (the versions involved are R2020b and R2012b). Does anyone have any idea why this is happening?
  2 Comments
Alessandro Gottardi
Alessandro Gottardi on 24 Dec 2020
Edited: Alessandro Gottardi on 24 Dec 2020
clear all
close all
n = 20000;
maxiter = 100000;
epsilon = 1.0e-10;
A = sparse(-2*eye(n)+diag(ones(n-1,1),1)+diag(ones(n-1,1),-1));
y = linspace(1,5,n)';
b = A*y;
tcpu = cputime;
[x, z, k] = coniugato(-A,-b,epsilon,maxiter);
tcpu = cputime - tcpu;
if k==maxiter
disp('Attenzione: raggiunto il numero massimo di iterazioni');
end
if isempty(x); return; end
r = A*x-b;
fprintf('Dimensione %d, Iter %d, tcpu = %f, err = %e \n', n,k,tcpu,z(k))
and the function coniugato
function [x, z, k] = coniugato(A,b,epsilon,maxiter,x)
n = length(b);
[k,l] = size(A);
if k~=l || n~=k
x = [];
fprintf('Dimensioni di A e/o b non corrette \n');
return;
end
z = zeros(maxiter,1);
if nargin<5; x=zeros(n,1); end;
r = A*x-b;
p = r;
rr = r'*r;
for k=1:maxiter;
Ap = A*p;
pAp = p'*Ap;
alfa = rr/pAp;
x = x - alfa * p;
r = r - alfa * Ap;
rrn = r'*r;
z(k) = sqrt(rrn);
if z(k)<epsilon; break; end
beta = rrn/rr;
p = r + beta * p;
rr = rrn;
end
fprintf('Gradiente coniugato: ');
end

Sign in to comment.

Accepted Answer

Walter Roberson
Walter Roberson on 24 Dec 2020
It looks to me as if you are using mrdivide ( / ) a couple of times in your code. The / operator has been improved to be able to select more efficient algorithms in some cases, and the underlying high performance math libraries have been upgraded at least twice since 2012 (my memory is saying three times for Windows, but that would have to be confirmed.)
A few years ago, the upgrades to the underlying math libraries started causing error messages in code that had run before. What had happened was that more accurate algorithms were used, but the more accurate algorithms were also more likely to detect that the problem was nearly singular — that really the old libraries should have warned about singularity too but had been a little sloppy and had permitted calculations that were not numerically meaningful.

More Answers (0)

Categories

Find more on Signal Integrity Kits for Industry Standards in Help Center and File Exchange

Tags

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!