%% Data:
fs = 35000;
T = 100;
Nt = fs*T;
M = 60;
a1 = randn(M,Nt);
a2 = randn(M,Nt);
lambda = randn(1,Nt);
B = randn(M,Nt);
result = zeros(2,Nt);
est = zeros(M,Nt);
A = [a1(:,1) a2(:,1)];
%% Vectorised method:
tic
num = (repmat(A'*A, 1, 1, Nt) + repmat(reshape([zeros(size(lambda)); lambda], 1, 2, []), 2, 1)); % \ (A' * B)
dem = reshape((A' * B), 2, 1, []);
resultVectorised = reshape(pagemldivide(num,dem), 2, []);
estVectorised = a1 .* result2(1,:) + a2 .* result2(2,:);
timeVectorised = toc;
%% Loop method:
tic
for i = 1 : Nt
result(:,i) = (A'*A + lambda(1,i) * [0 1 ; 0 1]) \ (A' * B(:,i));
est(:,i) = a1(:,i) * result(1,i) + a2(:,i) * result(2, i);
end
timeLoop = toc;
%% Output:
fprintf('Vectorisation: %.3f s\n', timeVectorised)
fprintf('Loop: %.3f s\n', timeLoop)
fprintf('Improvement: %.3f\n', timeLoop / timeVectorised)
fprintf('Absolute error: %e\n', max(max(abs(est - estVectorised))))
Vectorisation: 1.080 s
Loop: 20.117 s
Improvement: 18.632
Absolute error: 1.776357e-15
I think the vectorised method can be improved since, personally, I found hard to implement it without precisely know what I am doing.
Moreover, I think the small numerical error occurs due to pagemldivide which might use a different method from \, but I am not sure.
