MATLAB Answers

0

'Long-runtime problem', improve code, partial sums.

Asked by alberto martin on 8 Sep 2019 at 23:04
Latest activity Edited by alberto martin on 14 Sep 2019 at 15:51
Hello, i am using Matlab R2019a and i am trying to make a video of the dynamics of the partial sums of the Riemman Zeta function using matlab plot frames.
I have a few problems:
1) Memory grows until matlab crash before 2 or 5 hours depending on the size.
2)Slow execution , i tried to vectorize the code but due to the size I must use loops to not overflow the memory, the goal is n=500k.
I am looking for some way to improve the code or look for...something more efficient to carry out the task.
Thanks
n=50000;%The sum will be from 1 to n
m=1:1:n;%vector from 1 to n
partes=80;% divide the calculations in 'partes' parts
na=2;%the real part of 's' will be from 0 to 1 in na parts
b=10;%the imaginary part of 's' will be from 0 to b
nb=50;%imaginary part from 0 to b and nb steps
f=zeros(1,length(m));%preallocate
M=zeros(n/partes,1);%preallocate
CS=zeros(n/partes,n);%preallocate
count=0;%savenamecount...
for u=1:1:na%loop real part
for j=1:b*nb% loop imaginary part
tic
f=(exp(-1i*(j/nb)*log(m))./((m).^(u/na))).';%array with each value of the n component of the sum
for pp=0:partes-1%divide the (n,n) matrix into (n/partes,n)
CS=partess(partes,n,pp);%Generate parts of the full matrix
M=CS*f;%matrix*vector
%------plots--------
plot([-10;10],[0;0],'g');
hold on;
plot([0;0],[-10;10],'g');
axis([-10 10 -10 10]);
xlabel 'n';
xlabel 'real';
ylabel 'imag';
title({ ['n=' ,num2str(n)];['a=',num2str(u/na)];['b=' ,num2str(j/nb)]});
grid 'on';
grid 'minor';
plot(real(M),imag(M),'o','MarkerSize',1.5,'MarkerFaceColor','r','MarkerEdgeColor','r');
hold on;
plot(real(M),imag(M),'b');
set(gca,'color','k');
end
plot(real(M(n/partes)),imag(M(n/partes)),'o','MarkerSize',1.5,'MarkerFaceColor','g','MarkerEdgeColor','w');
drawnow;
clf;
cla;
%-----save the frame------
count=count+1;
temp=['~/matlab/imagenes/zeta/videos/frames/','zeta','-',num2str(count),'.jpg'];
saveas(gcf,temp)
toc
end
end
partess:
function C=partess(partes,n,pp)
C=zeros(n/partes,n);
for i=1:n/partes;
for k=1:n
if i>=k-(n/partes)*pp
C(i,k)=1;
end
end
end

  0 Comments

Sign in to comment.

Products


Release

R2019a

1 Answer

Answer by Steven Lord
on 8 Sep 2019 at 23:34

n=50000;%The sum will be from 1 to n
m=1:1:n;%vector from 1 to n
partes=80;% divide the calculations in 'partes' parts
na=2;%the real part of 's' will be from 0 to 1 in na parts
b=10;%the imaginary part of 's' will be from 0 to b
nb=50;%imaginary part from 0 to b and nb steps
f=zeros(1,length(m));%preallocate
M=zeros(n/partes,1);%preallocate
CS=zeros(n/partes,n);%preallocate
These last three lines aren't necessary. You're performing non-subscripted assignment to f, M, and CS in the loops which overwrite the vectors/matrices you've created here.
for u=1:1:na%loop real part
So the outermost loop runs 2 iterations
for j=1:b*nb% loop imaginary part
This one runs 500.
tic
f=(exp(-1i*(j/nb)*log(m))./((m).^(u/na))).';%array with each value of the n component of the sum
for pp=0:partes-1%divide the (n,n) matrix into (n/partes,n)
This one runs 80.
CS=partess(partes,n,pp);%Generate parts of the full matrix
M=CS*f;%matrix*vector
%------plots--------
plot([-10;10],[0;0],'g');
This line of code gets executed 2*500*80 = 80000 times, creating the same line each time! While you won't see any difference, each line is a separate object in memory. Move this outside the loops so it runs only once.
hold on;
plot([0;0],[-10;10],'g');
axis([-10 10 -10 10]);
xlabel 'n';
xlabel 'real';
ylabel 'imag';
These lines should be outside the loop as well, since they don't depend on any of your loop variables.
title({ ['n=' ,num2str(n)];['a=',num2str(u/na)];['b=' ,num2str(j/nb)]});
grid 'on';
grid 'minor';
plot(real(M),imag(M),'o','MarkerSize',1.5,'MarkerFaceColor','r','MarkerEdgeColor','r');
hold on;
plot(real(M),imag(M),'b');
Rather than replotting these lines (adding two new line objects to your axes each and every iteration of your inner loop) you probably want to call plot with an output argument before the loop (plotting NaN vectors of the correct size) and update the XData and YData properties of those lines inside the loop. That way you only have one or two lines that keep getting reused rather than creating tens of thousands of lines stacked on top of one another. This is one of the techniques described in the animation section of the documentation.
set(gca,'color','k');
Outside the loops. Generally, if a section of your code doesn't depend on the loop variables, it probably shouldn't be inside the loops. [One exception is the drawnow statement you have a couple lines down in the code. Leave that to give MATLAB an opportunity to update the graphics with the new data.]
end
plot(real(M(n/partes)),imag(M(n/partes)),'o','MarkerSize',1.5,'MarkerFaceColor','g','MarkerEdgeColor','w');
drawnow;
%-----save the frame------
temp=['~/matlab/imagenes/zeta/videos/frames/','zeta','-',num2str(count),'.jpg'];
saveas(gcf,temp)
This line is only inside two of your three loops, meaning it only runs 1000 times. This is probably still going to be a bottleneck, but if you want each of your images saved it's an unavoidable one.
Except ... this is the only place in your code count is referenced. I assume you cut out some code when you posted it that created and initialized the count variable and updated it inside the middle loop? If not you could simplify this code a whole lot by simply running the body of the middle loop once for the last element of the vectors over which the two outermost loops iterate.
toc
end
end

  3 Comments

Sorry!
First, thanks for your reply.
Second,i forgot some lines when I cleaned the code, now it's edited.
As you can see now, i put these lines inside the loop because i need to delete the draw, if not they overlap.
%------plots--------
plot([-10;10],[0;0],'g');
hold on;
plot([0;0],[-10;10],'g');
axis([-10 10 -10 10]);
xlabel 'n';
xlabel 'real';
ylabel 'imag';
I will read the Animation section.
The count variable is started with the other variables, but delete it when I clean the code to post it, I only use it as a save counter.
Thanks
So I see you're using clf and cla in your updated code. That means you're not actually creating thousands and thousands of line objects, but you are recreating the same lines (visually) over and over. Rather than clearing the figure and axes, if you call plot with an output argument you could delete those output arguments to eliminate just the lines being plotted, leaving these axes lines in place. That would save you the time to clear the figure and axes, replot the lines in the axes, and reconfigure the axes each time.
Actually, now that I see the movie, there are easier ways to achieve your goal. When you create your axes, set its XAxisLocation and YAxisLocation properties each to 'origin' rather than their default values of 'bottom' and 'left' respectively. Just like the plot function accepts property names 'MarkerSize', 'MarkerFaceColor', and 'MarkerEdgeColor' the axes function accepts property names 'XAxisLocation' and 'YAxisLocation'.
Alternately, if you need the lines and the tick labels to be different colors you could use the xline and yline functions to simulate axis lines through the origin.
Hello
i wrote the code following your tips and now the plot part is a lot faster.
I reallized that instead of vectorize the problem,iterated sums are ridiculously faster.
Here is the new code
clf
cla
n=num;%The sum will be from 1 to n
m=1:1:n;%vector from 1 to n
na=1;%the real part of 's' will be from 0 to 1 in na parts
b=2;%the imaginary part of 's' will be from 0 to b
nb=5;%imaginary part from 0 to b and nb steps
count=1;%savenamecount...
nan=NaN(n,1);%initiallize plots
S=zeros(n,1);
%------0plot-----
plot([-10;10],[0;0],'g');
hold on;
plot([0;0],[-10;10],'g');
axis([-10 10 -10 10]);
xlabel 'n';
xlabel 'real';
ylabel 'imag';
set(gca,'color','k');
grid 'on';
grid 'minor';
h3=plot(nan,nan,'o','MarkerSize',1.5,'MarkerFaceColor','g','MarkerEdgeColor','w');
h1=plot(nan,nan,'o','MarkerSize',1.5,'MarkerFaceColor','r','MarkerEdgeColor','r');
h2=plot(nan,nan,'b');
t=0;
tic
for u=1:1:na %loop real part
for j=1:b*nb %loop imaginary part
f=(exp(-1i*(j/nb)*log(m))./((m).^(u/na))).';%array with each value of the n component of the sum
S(1)=f(1);%first element
for pp=2:n
S(pp)=S(pp-1)+f(pp);%ppth partial sum
end
set(h1,'Xdata',real(S),'Ydata',imag(S)); %modify the plot
set(h2,'Xdata',real(S),'Ydata',imag(S));
set(h3,'Xdata',real(S(n)),'Ydata',imag(S(n)));%final point.
title({ ['n=' ,num2str(n),' partes=1',' NOGPU'];['a=',num2str(u/na)];['b=' ,num2str(j/nb)]});
drawnow;
%-----save frame
temp=['~/matlab/imagenes/zeta/videos/frames/','zeta','-',num2str(count),'.jpg'];
export_fig(gcf,temp);
count=count+1;%savenamecount
end
end
toc
t=toc;
1)The function export_fig is really slow ~0.3x sec, maybe there is a faster way to save the frames?
2)There is some way to parallelize the code? i try using parfor with sum() function but it is slower.
I had made plots with time comparision to show the difference, but this is nonsense because with the new code i can calculate each plot....say with n=2millon in 8 sec and with the older one
in 8 sec n~5k.........
Thanks

Sign in to comment.