MATLAB, умножение матриц против цикла для каждого столбца

Question

При умножении двух матриц, я попробовал следующие два варианта:

1)

res = X*A; 

2)

for i = 1:size(A,2) 
    res(:,i) = X*A(:,i); 
end 

Я выделенная память для res в обоих. И удивительно, что я нашел вариант 2 быстрее.

Может кто-нибудь объяснить, как это происходит?

редактировать: Я попытался

K=10000; 
clear t1 t2 
t1=zeros(K,1); 
t2=zeros(K,1); 

for k=1:K 
    clear res 
    x = rand(100,100); 
    a = rand(100,100); 
    tic 
    res = x*a; 
    t1(k) = toc; 
end 

for k=1:K 
    clear res2 
    res2 = zeros(100,100); 
    x = rand(100,100); 
    a = rand(100,100); 
    tic 
    for i = 1:100 
     res2(:,i) = x*a(:,i); 
    end 
    t2(k) = toc; 
end 

Answer 1

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

Перед тем, как версия Matlab 2008a (или версия рядом с этим выпуском), для циклов сильно пострадали в любом коде Matlab, потому что интерпретатор (слой между очень читаемым скриптом и реализация кода более низкого уровня) должен был бы повторно интерпретировать код каждый раз через цикл for.

С этого выпуска интерпретатор стал прогрессивно лучше, поэтому при запуске современной версии Matlab интерпретатор может посмотреть ваш код и сказать «Ах, я знаю, что он делает, позвольте мне оптимизировать его просто бит "и избегать попадания, которое в противном случае это потребовало бы путем переинтерпретации кода.

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

Один широкий урок, который мы должны извлечь из этого, не все версии равны. Я работаю над несколькими случаями кровотечения, используя два дополнения Matlab, SimBiology и Parallel Computing Toolboxes, оба из которых (особенно если вы хотите, чтобы они работали вместе) зависят от скорости исполнения и время от времени другие проблемы стабильности. Таким образом, я сохраняю три последних выпуска Matlab, проверю, что получаю одинаковые ответы из каждой версии, и я иногда буду возвращаться к более ранней версии, если найду проблемы с некоторыми функциями. Это, вероятно, слишком много для большинства людей, но дает вам представление о различиях в версии.

Надеюсь, это поможет.

редактирует:

Чтобы уточнить, код векторизации по-прежнему важно. Но данный сценарий, как:

x_slow = zeros(1,1e5); 
x_fast = zeros(1,1e5); 


tic; 
for i=1:1e5 
    x_slow(i) = log(i); 
end 
time_slow = toc; % evaluates for me in .0132 seconds 

tic; 
x_fast = log(1:1e5); 
time_fast = toc; % evaluates for me in .0055 seconds 

Несоответствие между time_slow и time_fast уменьшила в последние несколько версий, основанных на совершенствовании переводчика. Пример, который я видел, я считаю, был в 2000a против 2008b, но это зависит от моего воспоминания.

Есть что-то еще, что может произойти, на что обратились Оли и Юк. Существует часто разница между time_1 и time_2 в:

tic; x = log(1:1e5); time_1 = toc 
tic; x = log(1:1e5); time_2 = toc 

Таким образом, тест одного миллиона оценок против одной оценки является ценным, в зависимости от того, где в памяти х является (в кэше или нет).

Надеюсь, что это поможет снова.

Answer 2

Мне кажется, что вы не умножение матрицы должным образом, вы должны суммировать все продукты из й строки X матрицы и -го столбце матрицы, которые могли бы быть причиной. Посмотрите here, чтобы посмотреть, как это делается.

Answer 3

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

Answer 4

Я запускаю оба кода в цикле 1000 раз. В среднем (но не всегда) первый векторизованный код был в 3-4 раза быстрее. Я очистил переменные результата и предварительно выделил до запуска таймера.

x = rand(100,100); 
a = rand(100,100); 

K=1000; 
clear t1 t2 
t1=zeros(K,1); 
t2=zeros(K,1); 

for k=1:K 
    clear res 
    tic 
    res = x*a; 
    t1(k) = toc; 
end 

for k=1:K 
    clear res2 
    res2 = zeros(100,100); 
    tic 
    for i = 1:100 
     res2(:,i) = x*a(:,i); 
    end 
    t2(k) = toc; 
end 

Итак, никогда не делайте временного заключения на основе одного прогона.