I wrote the code below to invert an upper triangular matrix, avoiding any possible multiplication/subtraction by zero. It just uses $\frac{1}{6}n^3+\ldots$ flops instead of $n^3+\ldots$ flops.
function B = InvTrMat(A,n,type)
% This code inverts an upper triangular matrix of order n
% without doing any multiplication or subtraction by zero.
B = zeros(n);
for j = 1:n
B(j,j) = 1 / A(j,j);
for i = j-1:-1:1
for k = i+1:j
B(i,j) = B(i,j) - A(i,k) * B(k,j);
end
B(i,j) = B(i,j) * B(i,i);
end
end
return
end
The problem is it is slower than the backslash command of Matlab. I used the following code to test and compare both:
n = 128;
t = zeros(1,100);
u = zeros(1,100);
for i = 1:100
A = rand(n,n);
A = A+A';
A = A + n*eye(n);
R = chol(A);
tic;
B = R \ eye(n);
t(i) = toc;
tic;
B = InvTrMat(R,n,'u');
u(i) = toc;
end
sum(t)/100
sum(u)/100
Does anyone know how to make my code faster than Matlab's code? It must be possible since I'm using just a fraction of flops from the other code.
My Matlab version is 7.10.0 (R2010a), installed on a Windows 8.0 64-bits PC.
Thank you, everyone.