矢量

乍一看for循环代码告诉我们,由于photAbs是一个二进制数组,每列都根据Data的每个元素进行缩放,因此这个二进制特征可用于矢量化.这在代码中滥用 –

function RNGData = RNGAnsys_vect1(N,Data)

%// Get the 2D Matrix of random ones and zeros

photAbsAll = randint(N,numel(Data));

%// Take care of multData internally by summing along the columns of the

%// binary 2D matrix and then multiply each element of it with each scalar

%// taken from Data by performing elementwise multiplication

sumData = Data.*sum(photAbsAll,1);

%// Divide by n, but account for 0.5 average by n/2

RNGData = (sumData./(N/2))'; %//'

return;

在分析之后,似乎瓶颈是随机二进制数组创建部分.因此,使用this smart solution中建议的更快的随机二进制数组创建器,可以进一步优化上述功能 –

function RNGData = RNGAnsys_vect2(N,Data)

%// Create a random binary array and sum along the columns on the fly to

%// save on any variable space that would be required otherwise.

%// Also perform the elementwise multiplication as discussed before.

sumData = Data.*sum(rand(N,numel(Data))<0.5,1);

%// Divide by n, but account for 0.5 average by n/2

RNGData = (sumData./(N/2))'; %//'

return;

使用智能二进制随机数组创建器,原始代码也可以进行优化,以便稍后在优化的for循环和矢量化代码之间进行公平的基准测试.这里列出了优化的for循环代码 –

function RNGData = RNGAnsys_opt1(N,Data)

multData = zeros(N,numel(Data));

for i = 1:numel(Data)

%// Create N number of random 0's or 1's using a smart approach

%// Then, multiply each element in the molar data by the random numbers

multData(:,i) = Data(i) * rand(N,1)<.5>

end

sumData = sum(multData,1); % sum each individual energy level's data point

RNGData = (sumData/(N/2))'; % divide by n, but account for 0.5 average by n/2

return;

标杆

基准代码

N = 15000; %// Kept at this value as it going out of memory with higher N's.

%// Size of dataset is more important anyway as that decides how

%// well is vectorized code against a for-loop code

DS_arr = [50 100 200 500 800 1500 5000]; %// Dataset sizes

timeall = zeros(2,numel(DS_arr));

for k1 = 1:numel(DS_arr)

DS = DS_arr(k1);

Data = rand(1,DS);

f = @() RNGAnsys_opt1(N,Data);%// Optimized for-loop code

timeall(1,k1) = timeit(f);

clear f

f = @() RNGAnsys_vect2(N,Data);%// Vectorized Code

timeall(2,k1) = timeit(f);

clear f

end

%// Display benchmark results

figure,hold on, grid on

plot(DS_arr,timeall(1,:),'-ro')

plot(DS_arr,timeall(2,:),'-kx')

legend('Optimized for-loop code','Vectorized code')

xlabel('Dataset size ->'),ylabel('Time(sec) ->')

avg_speedup = mean(timeall(1,:)./timeall(2,:))

title(['Average Speedup with vectorized code = ' num2str(avg_speedup) 'x'])

结果

结束语

根据我迄今为止使用MATLAB的经验,循环和矢量化技术都不适合所有情况,但一切都是针对具体情况的.