Asked by alberto martin
on 8 Sep 2019 at 23:04

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

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

alberto martin
on 9 Sep 2019 at 0:35

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';

The movie is like this https://www.youtube.com/watch?v=9L_Fdqtsa_U

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

Steven Lord
on 9 Sep 2019 at 17:27

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.

alberto martin
on 14 Sep 2019 at 0:49

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.

Opportunities for recent engineering grads.

Apply Today
## 0 Comments

Sign in to comment.