Performance analysis of series summation in Matlab

问题:

I'm writing a Matlab program to compute pi through the summation of series

A = Sum of a_i from i=1 to N

where

pi/4 = 1 - 1/3 + 1/5 - 1/7 + 1/9 - 1/11 + 1/13 ...

To compute pi through the series summation, the suggested approach is to set

a_i = (-1)^(i+1)/(2i-1)

In order to do this, I wrote the program below

n=100;
f=[];
for jj=1:n
    ii=1:jj;
    f=[f 4*sum( ((-1).^(ii+1))./(2.*ii-1)  )];
end;
hold on
plot(f)
title('Computing of \pi using a finite sum')
xlabel('Number of summations') 
ylabel('Estimated value of \pi') 
plot(1:size(f,2),ones(size(f))*pi)

This program shows that the series approximation is somewhat accurate near N=80.

I am now attempting to adjust my program so that the y-axis displays total calculation time T_N and the x-axis displays N (the number of summations). The total calculation time T_N should increase as N increases. Ideally, I am looking to have the graph display something close to a linear relationship between T(N) and N

To do this, I have adjusted my original program as follows

n=100;
f=[];
tic
for jj=1:n
    ii=1:jj;
    f=[f 4*sum( ((-1).^(ii+1))./(2.*ii-1)  )];
end;
hold on
plot(f)
title('Time it takes to sum \pi using N summations')
xlabel('Number of summations (N)') 
ylabel('Total Time (T_N)') 
plot(1:size(f,2),toc)
slope = polyfit(1:size(f,2),toc,1);

This looks wrong. I must have incorrectly applied the built-in timing functions in Matlab (tic and toc). So, I am going to analyze my code and ask two questions -

1. How could I adjust my code above so that the y-axis correctly displays the total calculation time per the summation N? It looks like I did something wrong in plot(1:size(f,2),toc).

2. After I get the y-axis to display the correct total calculation time (T_N), I should be able to use the polyfit command to find the slope of T(N)/N. This will give me a linear relationship between T(N) and N. I could then use the value of slope = polyfit(1:size(f,2),toc,1) to compute

   t_N = a + b*N

where t_N is computed for every value of N and b is the slope calculated through the polyfit command.

I think that I should be able to find the values of a and b after I correctly display the y-axis and correctly reference the polyfit command.

回答1:

If I understand your problem correctly, I think there are two different issues here. First, you plot your result function then the elapsed time which is several orders of magnitude smaller than pi:

 hold on
 plot(f)  % <---- Comment me out!
 ...stuff...
 plot(1:size(f,2),toc)

Secondly, you need to store the execution time of each pass of the loop:

n=100;
f=[];
telapsed = zeros(1,n);
tic
for jj=1:n
    ii=1:jj;
    f=[f 4*sum( ((-1).^(ii+1))./(2.*ii-1)  )];
    telapsed(jj) = toc;
end
hold on
% plot(f)
title('Time it takes to sum \pi using N summations')
xlabel('Number of summations (N)') 
ylabel('Total Time (T_N)') 
plot(1:n,telapsed)
slope = polyfit(1:n,telapsed,1);

Note the new polyfit expression for slope of the execution time. Does that help?

回答2:

There are several things that can be improved in your code:

• f should be preallocated, so as not to waste time in repeatedly assigning memory.
• tic should be called within the loop, in order to restart the stopwatch timer.
• When you call toc you get the current time from the last tic. The time spent should be stored in a vector (also preallocated).
• Since the computations you want to time are very fast, measuring the time they take is very unrealiable. The computations should be repeated many times, so the measured time is larger and you get better accuracy. Even better would be to use timeit (see below).
• You cannot plot the time and the results in the same figure, because the scales are too different.

The code incorporating these changes is:

n = 100;
f = NaN(1,n); % preallocate
times = NaN(1,n); % preallocate
repeat_factor = 1e4; % repeat computations for better time accuracy
for jj=1:n
    tic % initiallize time count
    for repeat = 1:repeat_factor % repeat many times for better time accuracy
        ii=1:jj;
        f(jj) = 4*sum( ((-1).^(ii+1))./(2.*ii-1) ); % store value
    end
    times(jj) = toc; % store time
end
times = times / repeat_factor; % divide by repeat factor
plot(f)
title('Time it takes to sum \pi using N summations')
xlabel('Number of summations (N)') 
ylabel('Total Time (T_N)')
figure % new figure for time
plot(1:size(f,2), times)
p = polyfit(1:size(f,2),times,1);
slope = p(1);

---

Using timeit for measuring the time will probably give improved accuracy (but not very good because, as mentioned above, the computations you want to time are very fast). To use timeit you need to define a function with the code to be timed. The simplest way is to use an anonymous function without input arguments. See code below.

n = 100;
f = NaN(1,n); % preallocate
times = NaN(1,n); % preallocate
for jj=1:n
    ii=1:jj;
    fun = @() 4*sum( ((-1).^(ii+1))./(2.*ii-1) );
    f(jj) = fun(); % store value
    times(jj) = timeit(fun); % measure and store time
end
plot(f)
title('Time it takes to sum \pi using N summations')
xlabel('Number of summations (N)') 
ylabel('Total Time (T_N)')
figure % new figure for time
plot(1:size(f,2), times)
p = polyfit(1:size(f,2),times,1);
slope = p(1);