Почему эта программа рисования эллипса очень медленная? - PullRequest
0 голосов
/ 26 апреля 2019

У меня есть программа для рисования сетки эллипсов с равномерным распределением фаз. Тем не менее, это очень медленно.

Я бы хотел, чтобы мой код был быстрее, чтобы я мог использовать, например, N = 150, M = 150. Как я могу ускорить этот код?

N = 10;
M = 10;
y = 1;
x = 1;
a = 1;
b = 2;
for k = 1:N
  for m = 1:N
    w = rand(1,1);
    for l = 1:N
      for s = 1:N
        if(((l-x)*cos(w*pi)+(s-y)*sin(w*pi)).^2/a^2 + (-(l-x)*sin(w*pi) + (s-y)*cos(w*pi)).^2/b.^2 <= 1)
          f(l,s) = 1*(cos(0.001)+i*sin(0.001));
        end
      end
    end
    y = y+4;
  end
  y = 1;
  x = x+5;
end
image(arg(f),'CDataMapping','scaled');

Вот что выдает код:

Picture

Обновлен:

N = 10;
M = 10;
y = 1;
x = 1;
a = 1;
b = 2;
for x = 1:5:N
  for y = 1:4:N
    w = rand(1);
    for l = 1:N
      for s = 1:N
        if(((l-x).*cos(w.*pi)+(s-y).*sin(w.*pi)).^2/a.^2 + (-(l-x).*sin(w.*pi) + (s-y).*cos(w.*pi)).^2/b.^2 <= 1)
          f(l,s) = cos(0.001)+i.*sin(0.001);
        end
      end
    end
  end
  y = 1;
end
image(arg(f),'CDataMapping','scaled');

Ответы [ 2 ]

1 голос
/ 27 апреля 2019

Есть много вещей, которые вы можете сделать, чтобы ускорить вычисления.Одним из важных является удаление циклов и замена их векторизованным кодом.Octave работает намного быстрее, когда выполняет много вычислений по одному, а не по одному за раз.

Например, вместо

for l = 1:N
  for s = 1:N
    if(((l-x).*cos(w.*pi)+(s-y).*sin(w.*pi)).^2/a.^2 + (-(l-x).*sin(w.*pi) + (s-y).*cos(w.*pi)).^2/b.^2 <= 1)
      f(l,s) = cos(0.001)+i.*sin(0.001);
    end
  end
end

можно написать

l = 1:N;
s = (1:N).';
index = ((l-x).*cos(w.*pi)+(s-y).*sin(w.*pi)).^2/a.^2 + (-(l-x).*sin(w.*pi) + (s-y).*cos(w.*pi)).^2/b.^2 <= 1;
f(index) = cos(0.001)+i.*sin(0.001);

Однако здесь мы все еще выполняем слишком много работы, потому что мы вычисляем index в тех местах, которые, как мы знаем, находятся за пределами эллипса.В идеале мы бы нашли меньшую область вокруг каждой точки (x,y), в которую, как мы знаем, эллипс вписывается.

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

Также выполняется много избыточных вычислений.Например, w.*pi вычисляется несколько раз, и cos и sin тоже.Вы также назначаете cos(0.001)+i.*sin(0.001) для вывода пикселей снова и снова, это может быть константа, вычисляемая один раз.

Следующий код запускается в MATLAB за доли секунды (хотя Octave будет намного медленнее),Я также правильно разделил N и M (поэтому выходные данные не всегда квадратные) и сделал размеры шагов переменными для лучшего понимания кода.Я присваиваю 1 элементам эллипса, вы можете заменить их своей константой, умножив f = f * (cos(0.001)+i*sin(0.001)).

N = 150;
M = 200;
a = 5;
b = 10;
x_step = 25;
y_step = 25;
f = zeros(N,M);
for x = x_step/2:x_step:M
  for y = y_step/2:y_step:N
    phi = rand(1)*pi;
    cosphi = cos(phi);
    sinphi = sin(phi);
    l = (1:M)-x;
    s = (1:N).'-y;
    index = (l*cosphi+s*sinphi).^2/a.^2 + (-l*sinphi + s*cosphi).^2/b.^2 <= 1;
    f(index) = 1;
  end
end

output of algorithm

0 голосов
/ 27 апреля 2019

Я не на 100% уверен, что ты пытаешься сделать. Посмотрите код ниже и дайте мне знать. Мне потребовалось около 30 секунд, чтобы запустить дело 150 х 150. Не уверен, что это достаточно быстро

[h,k] = meshgrid(0:150, 0:150);

a = 0.25;
b = 0.5;

phi = reshape( linspace( 0 , 2*pi , numel(h) ), size(h));

theta = linspace(0,2*pi,50)';

x = a*cos(theta);
y = b*sin(theta);

h = h(:);
k = k(:);
phi = phi(:);

ellipseX = arrayfun(@(H,K,F) x*cos(F)-y*sin(F) + H , h,k,phi, 'uni', 0);
ellipseY = arrayfun(@(H,K,F) x*sin(F)+y*cos(F) + K , h,k,phi, 'uni', 0);

ellipse = [ellipseX, ellipseY, repmat({'r'}, size(ellipseX))]';

fill(ellipse{:})
...