Need little help in making matlab code related with iteration methods

I have made following matlab code for computing the Moore-Penrose inverse of a given matrix A.

A = randn(100)            % given matrix

beta = 1/norm(A,2)^2

x0 = beta.*A'             % initial approximation, A' is the transpose of matrix A

k = 0;

iter = 0

f = 1 ;

I = eye(100);

while (f > 1.0e-007)

x1 = 2*x0 - x0*A*x0;   % x1 is approximation of Moore-Penrose inverse of matrix A

iter = iter+1

a = norm(A*x1*A-A, 2);  %  error norm

b = norm(x1*A*x1 - x1);   % error norm

c = norm ((A*x1)' - A*x1, 2 );   % error norm

d = norm ((x1*A)' - x1 *A, 2);   % error norm

B = [a, b, c, d];

f = max(B);

x0 = x1;

end


Since matrix A is changing after every iterations hence values of a,b ,c and d are also changing. I am trying to compute average of values a, b, c and d after hundred repetitions but I am unable to do so. Could anybody help me? I would be very much thankful to you.

-
Change your code to "B(iter,:) = [a,b,c,d];" and "f = max(B(iter,:));" This will store the value of a,b,c,d for each iteration as successive rows of the matrix B. After the loop ends, you can find the averages using "mean_a = mean(B(:,1));" and similarly for the others. – in_wolframAlpha_we_trust Sep 25 '12 at 10:37
@in_wolfram_we_trust i did this although value of a is in the power s of $10^{-13}$ but mean a is giving me results like 1.33...while it should be similar to a..since it is simply mean of values of a after 100 repetitions... – srijan Sep 25 '12 at 10:50
I've just run the code as well at found mean_a = 0.8126. The difference in values is due to a different initial seed, but both these answers are acceptable. If you use plot(B(:,1)) it shows you how the value of a has changed - while the value of a at the last iteration is very small (in my case ~10^-11) it starts off quite large (in my case ~10). Averaging over the whole set of values gives the answer you think is too big. – in_wolframAlpha_we_trust Sep 25 '12 at 11:04
@srijan wait you mean after 100 iterations or after rerunning the code with 100 different matrices A ? – k1next Apr 11 '13 at 11:32
So lets say for a single matrix your method converges, you get $a=10^{-8}$ and for the next matrix your method converges with $a=10^{-9}$ you would want to get $mean(a) = 5 \cdot 10^{-8}$ ? – k1next Apr 11 '13 at 11:45

I think the reason, why the mentioned approaches didn't work was because everyone was assuming ( me too at first ) that you want to have the mean value for a during the iterations with one single matrix. However, you seem to want the mean value of a at the end of the algorithm, averaged for several runs with different matrices.

If I understand this correct, the below code should work. (Note $M=10$ instead of $M=100$, line 2)

clc              % cleans the output screen, I don't like old information...
M = 10;          % number of matrices for which perform the calculation
E = zeros(M,4);  % we will store for all M iterations all the 4 values

fprintf('\nPerforming %i iterations in total.\n',M);  % this is just output..

%now the important part, we will perform your calculation M times
%therefore, we loop over your calculation with a for loop
for i = 1:M

fprintf('Iteration %i is running...\n' ,i ) % output..

A = rand(100);            % given matrix
beta = 1/norm(A,2)^2;
x0 = beta.*A';             % initial approximation, A' is the transpose of matrix A
k = 0;
iter = 0;
f = 1 ;
I = eye(100);

while (f > 1.0e-7)

x1 = 2*x0 - x0*A*x0;   % x1 is approximation of Moore-Penrose inverse of matrix A
iter = iter+1;
a = norm(A*x1*A-A, 2);  %  error norm
b = norm(x1*A*x1 - x1);   % error norm
c = norm ((A*x1)' - A*x1, 2 );   % error norm
d = norm ((x1*A)' - x1 *A, 2);   % error norm

B = [a, b, c, d];

f = max(B);

x0 = x1;

end

% at this point, your calculation has been performed for a matrix A.
% We will now store  the values for a,b,c,d of the last iteration in the matrix E
E(i,:) = B;
end

% here we print the mean value of a, which is mean(E(:,1))
fprintf('\nMean value of a =%s\n', mean(E(:,1)));
`
-
I have no words to thank you. I just want your little explanation about the code. I am not very much good in matlab coding. Thanks again. :) – srijan Apr 11 '13 at 12:05
No hurry. I will add some comments to the code soon. – k1next Apr 11 '13 at 12:14
Thanks again sir. – srijan Apr 11 '13 at 12:19
@srijan Comments added. I think that everybody else also could have helped, but the problem was that we first missinterpreted your question ;) – k1next Apr 11 '13 at 12:20
Dear sir everything is clear now. No queries. Heartily thanks. :) – srijan Apr 11 '13 at 12:28