%% Test different sparse summation implementation
% Run the loop until time between outputs becomes unbearable, then cancel execution and plot results

% starting parameters
N = 10000;
density = .001;

previousNumLoops = 0;
numLoops = 2;

runtimes = struct;
ind = 1;

while 1
    % Setup dummy data matrices
    for i = (previousNumLoops+1):numLoops
        Rcell(i) = {sprand(N,N,density)};
    end

    runtimes(ind).numLoops = numLoops;

    % Evaluate implementation performance matrix multiplication
    runtimes(ind).ConstructSparse = ConstructSparse(N, Rcell, numLoops);
    runtimes(ind).AddSparse = AddSparse(N, Rcell, numLoops);

    disp(runtimes(ind))

    ind = ind + 1;
    previousNumLoops = numLoops;
    numLoops = numLoops + 2;

end

%% Plot runtime results

% Run to create a plot of the results
PlotRuntimes(runtimes, {'ConstructSparse','AddSparse'});

%% Functions to time implementation and plot runtimes according to scaling factor

% measure accumulation of nonzero value properties and construction of the sparse matrix after loop
function runtime = ConstructSparse(N, Rcell, numLoops)
tic

sparseInds = cell(numLoops,3);
for i = 1:numLoops
    [rowInds,colInds,nonzeroVals] = find(Rcell{i});
    sparseInds(i,:) = {rowInds,colInds,nonzeroVals};
end

I = vertcat(sparseInds{:,1});
J = vertcat(sparseInds{:,2});
nonzeroValues = vertcat(sparseInds{:,3});

A = sparse(I,J,nonzeroValues,N,N);

runtime = toc;

end

% Measure summing of sparse matrices within loop
function runtime = AddSparse(N, Rcell, numLoops)
tic

A = sparse(N,N);
for i = 1:numLoops
    A = A + Rcell{i};
end

runtime = toc;

end

function [] = PlotRuntimes(runtimes, testNames)
% Function to plot runtime results
figure
hold on;
for label = testNames
    plot([runtimes(1:(end-1)).numLoops],[runtimes(1:(end-1)).(label{1})])
end
hold off

xlabel('Loop parameter')
ylabel('Runtime(seconds)')

legend(testNames)

end