Как улучшить скорость работы этого кода Matlab?

Question

Я пытаюсь смоделировать движение некоторых частиц. Программа кажется очень медленной. Я только начал бережливое программирование, поэтому я не знаю, в чем именно проблема, которая делает его медленным, я думаю, что построение графика занимает много времени. может кто-нибудь подсказать мне, как его улучшить?

function folks(N , T)

 if nargin < 1
     N = 50; 
     T = 100; 
 end

A=20;R=50;a=100;r=2;

x0=10*randn(3,N);
v0=0*randn(3,N);

clear c
% Initilazing plot
color = hsv(N);
xh=zeros(1,N);
f=2*max(max(x0));ff=f/1000000;
figure(2),clf
set(gcf,'doublebuffer','on')
hold on, grid on, axis([-1 1 -1 1 -.5 .5]*f)

 for j = 1:N
      xh(j) = line(x0(1,j),x0(2,j),x0(3,j),'color',color(j,:), ...
         'marker','.','markersize',15);
 end
title('t = 0','fontsize',18)
rotate3d; 
view([1.8,-1.8,1])

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Solving the ode system
tol = 3e-14;
opts = odeset('reltol',tol,'abstol',tol);

X0=[x0(:);v0(:)];
dt=T/1000;

for tnext = dt:dt:10*T
    tspan = [tnext-dt tnext];
    [T1,X1] = ode45(@odefolk,tspan,X0,opts,r,R,a,A);

    X0 = X1(end,:);
    if max(abs(X0(1:3*N)))>f
        f=1.1*f;
        axis([-1 1 -1 1 -1 1]*f)
    end

    for j=0:N-1
        set(xh(j+1),'xdata',X0(3*j+1),'ydata',X0(3*j+2),'zdata',X0(3*j+3))
    end

    title(sprintf('t = %3.0f      %5.0e',tnext),'fontsize',18),
    drawnow
end

Я пытался сохранить данные в матрице, а затем отобразить все в новом цикле, но это сделало его намного медленнее. А также не относится к решателю ode45 или системе odefolk

function dx = fkt(~,x,r,R,a,A)
n = length(x);  % n = 90 = 6 *N 
    M = n/6;  % M= 15
    dx = zeros(n,1);

    for i = 1 : 3 : n/2

        vx = 0 ; vy = 0 ; vz = 0;

        for j = 1 : 3 : n/2

            if (i~=j)

                vx = x(i)   - x(j)  ;
                vy = x(i+1) - x(j+1);
                vz = x(i+2) - x(j+2);
%                 ex = expo (vx , vy , vz ,r , a);
%                 
%                 val = (ex(1) * R - ex(2) * A);
                leng = sqrt(vx^ 2 +vy^ 2 + vz^ 2);
                expo1  = exp(-1 * leng / r) / r;
                expo2  = exp(-1 * leng / a) / a;

                val = (expo1 * R - expo2 * A);

                vx = vx * val;
                vy = vy * val;
                vz = vz * val;

            end
        end

        vx = vx/M ;
        vy = vy/M;
        vz = vz/M;

        dx(i)     = x(i + 3 * M  );
        dx(i + 1) = x(i + 3 * M + 1);
        dx(i + 2) = x(i + 3 * M + 2);

        dx(i + 3 * M  ) = vx;
        dx(i + 3 * M + 1 ) = vy;
        dx(i + 3 * M + 2 ) = vz;
    end

end

Answer 1

Несколько советов:

• используйте tic и toc для измерения времени выполнения вашего кода (строк, функций и т. Д.)
• предварительно назначайте массивы и не позволяйте им расти в петлях
• не используйте plot в цикле. Попробуйте использовать команды низкого уровня и измените значение существующих линий, ссылаясь на их дескрипторы (например, handle = line([0 1][2 3][4 5]), а затем set(handle,'XDATA',[4 7],'YDATA',[2 6],'ZDATA',[2 4]))
• MATLAB - интерпретируемый язык, поэтому ваш код будет работать медленнее в большинстве случаев по сравнению со скомпилированными программами

Answer 2

Чтобы выполнить код, рекомендуется проанализировать его с помощью Profiler:

https://es.mathworks.com/help/matlab/ref/profile.html

Во всяком случае, я вижу некоторые вещи, которые вы могли бы изменить. Во-первых, избегайте рисовать в каждой итерации. Сделайте это в конце или нарисуйте каждое M взаимодействий.

Во-вторых, избегайте использования цикла для определения графики. Построить его с помощью входных данных векторного типа.

В-третьих, ode45 может выполнить итерации за вас и вернуть вам вектор со всеми взаимодействиями. Используйте это, чтобы построить.