About preallocating for speed

I'm trying to display the TSE and CR Linear and Non Linear compression graph but it's been 4 hours and the graph and calculation still hasn't been shown
clc
clear all
close all
[f,Fs] = audioread('C:\Users\fadli\Downloads\voice_1090883.wav');
N = length(f);
n = 0:N-1;
DC = mean(f);
M = N/2;
for K = 1:M
m = 1 : K;
G = cos(2*pi*(n')*m/N);
H = sin(2*pi*(n')*m/N);
A = (2/N)*(G'*f);
B = (2/N)*(H'*f);
fhat = DC + G*A + H*B;
e = f-fhat;
TSE_Lin(K) = norm(e)^2;
CR(K) = (2*K)*100/N;
end
plot(CR,TSE_Lin,'r--o')
m = 1 : M;
G = cos(2*pi*(n')*m/N);
H = sin(2*pi*(n')*m/N);
A = (2/N)*(G'*f);
B = (2/N)*(H'*f);
[~,i] = sort(abs(A),'descend');
[~,j] = sort(abs(B),'descend');
for K = 1:M
Gc = G(:,i(1:K));
Ac = A(i(1:K));
Hc = H(:,j(1:K));
Bc = B(j(1:K));
fhat = DC + Gc*Ac + Hc*Bc;
e = f-fhat;
TSE_Non_Lin(K) = norm(e)^2;
end
hold on
plot(CR,TSE_Non_Lin,'b--*')
legend('Linear Compression', 'Nonlinear Compression')
xlabel('Compression Ratio (%)')
ylabel('TSE')
xticks(4:8:100)
I've looked up threads like https://www.mathworks.com/matlabcentral/answers/258627-about-preallocating-for-speed and I still don't quite understand what and where i should put the zero in my code.
I apologize if this sounds like a silly question but I sincerely thank you anyone that'd kindly help me

2 Comments

Have you tried running it with a smaller file first? That way you can run the profiler. I doubt preallocation will be the solution.
I’ll try using a smaller file, but the audio file is only like two-three seconds long i don’t see how i can shortened it anymore haha. but yes, i’ll try using a shorter/smaller file

Sign in to comment.

Answers (1)

Jan
Jan on 30 Oct 2021
Edited: Jan on 30 Oct 2021
for K = 1:M
...
TSE_Lin(K) = norm(e)^2;
CR(K) = (2 * K) * 100 / N;
end
The variables TSE_lin and CR grow in each iteration. Preallocate the final size bfeore the loop:
TSE_lin = zeros(1, M);
CR = zeros(1, M);
for K = 1:M
...
TSE_Lin(K) = norm(e)^2;
CR(K) = (2 * K) * 100 / N;
end
Now the arrays are allocated once only and do not grow anymore.
But you will see, that this does not reduce the processing time remarkably. Use some small test data to run the profiler as Rik has suggested:
f = rand(1e3, 1) * 2 - 1; % Instead of the wav file
Most of the time is spent in the lines:
G = cos(2 * pi * (n') * m / N);
H = sin(2 * pi * (n') * m / N);
For N=1000 these are 250'500'000 sine and cosine operations, which are expensive. Creating the argument consumes some time also and at least this can be saved, because it is equal:
TSE_lin = zeros(1, M);
CR = zeros(1, M);
c1 = 2 * pi * n.' / N
for K = 1:M
m = 1 : K;
Arg = c1 * m;
G = cos(Arg);
H = sin(Arg);
...
TSE_Lin(K) = e.' * e; % Faster than: norm(e)^2;
CR(K) = (2 * K) * 100 / N;
end
For N=1000 this reduces the runtime by 10%. Not impressing.
The main problem remains the computations of the trigonometric functions. But your code computes the Argument of the 1st iteration in all following iterations again. The sine and cosine determined for m = 1 are calculated again for m = 1:2, m = 1:3 and so on. So just append the new values to the formerly calculated ones:
TSE_Lin = zeros(1, M);
c1 = 2 * pi * n.' / N;
fn = 2 * f / N;
ff = f - DC;
G = [];
H = [];
for K = 1:M
Arg = c1 * K;
G = [G, cos(Arg)];
H = [H, sin(Arg)];
A = G.' * fn;
B = H.' * fn;
e = ff - G*A - H*B;
TSE_Lin(K) = e.' * e;
end
CR = 200 / N * (1:M);
CR can be evaluated outside the loop easily.
Now G and H are growing iteratively again. But this is cheaper than cropping copies for the repeated calculations. The idea is to avoid repeated calculations either by moving it out of the loop (see c1, fc, ff), or re-using it in the following loops.
For N=1e3 the processing time is reduced from 2.79 to 0.14 seconds on my i7, Win10, Matlab R2018b system.
Unfortunately the computation time grows fast with the input. For the double size N=2000 it takes roughly 10 times longer.

Categories

Find more on Mathematics in Help Center and File Exchange

Products

Release

R2019a

Asked:

on 30 Oct 2021

Edited:

on 30 Oct 2021

Community Treasure Hunt

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

Start Hunting!