Результаты целочисленного теста

Ну вот, готовы результаты целочисленного теста после исправления досадной, грубой ошибки, допущенной мною по невнимательности. Зато добавился новый участник теста: уже немного старенький, но десктопный Intel i5-6500.

Напомню, что тест основан на вычислении обратного элемента в кольце вычетов по модулю N. И фактически тестирует регистровые целочисленные 64-битные операции в двух вариантах: с использованием операции деления или на основе операций сложения и сдвига.

Запустил сначала однопоточный тест на своем AMD FX-4350, получил результаты и подумал, что, наверное, делать многопоточный тест смысле нет: все же чистая синтетика, 100% регистровых операций, ни одного обращения к памяти... Но потом вспомнил про HyperThreading, CMT, SMT и вот это вот всё. Поэтому все же прогнал многопоточный тест ради интереса.

В качестве базы я, как обычно взял свой AMD FX-4350, поэтому его производительность везде принята за единицу. Производительность сложения и сдвига я считаю в 1.5 раза более ценной, чем операции деления (так как встречаются в коде они гораздо чаще), как и производительность однопоточного выполнения. Сначала интегральная производительность.

Как видно из графика, что-то производительность регистровых целочисленных операций у процессоров Intel особо не задалась. Лишь относительно новый i5-8300H в интегральной производительности слегка обгоняет мой старенький FX-4350.

Недавно читал, что за последние лет пять производительность процессоров выросла примерно вдвое. Ну, почти так и есть. Ryzen 5 3600 в 1.82 раза быстрее моего процессора. Но надо заметить, что в основном ускорение пошло экстенсивным путем ‒ за счет увеличения количества ядер, в то время как однопоточная производительность увеличилась всего лишь на жалкие 16%.

Несколько удивляет крайне низкий результат i5-6500, который лишь чуть-чуть быстрее весьма древнего мобильного AMD A10-4600M. Впрочем, тест проводил не я лично, возможно, тестирование было выполнено не совсем аккуратно.

Как можно видеть, в целочисленных регистровых операциях даже относительно новые, хоть и мобильные процессоры Intel с трудом тягаются со стареньким десктопным AMD. Ну а относительно новый Ryzen 5 просто закрепляет победу.
Так что, похоже, и в самом деле архитектура Ryzen очень хороша. Хоть она, похоже, чуть-чуть хуже в операциях с плавающей точкой, но опережает в целочисленных, что для большинства приложений важнее. За исключением, может быть, трехмерных игр и приложений ИИ.

Посмотрим более детально.

Вариант с использованием операции деления.

Если в функции встречается всего лишь одна операция деления, то тут даже старенький процессор от AMD рвет как Тузик грелку все процессоры от Intel, участвовавшие в тестировании. А Ryzen еще немножко добавляет.

Бинарный вариант на сложениях и сдвигах.

Тут Intel в последних поколениях мобильных процессоров в однопотоке немного обогнала PileDriver, но Ryzen 3000 это компенсировал и на обычных целочисленных регистровых операциях AMD и Intel сейчас на равных, но за счет большего количества ядер в общем зачете все же AMD впереди. Из любопытного: все процессоры, за исключением двух, показывают большее ускорение при параллельном выполнении функции на основе деления, что в общем-то понятно при использовании технологии HyperThreading: пока одно виртуальное ядро выполняет длинную операцию деления, другое ядро может вовсю пользоваться остальной частью физического АЛУ. Аномальное поведение замечено лишь у Ryzen 5 3600 и i5-6500. Не могу даже предположить, с чем это может быть связано.

Ну и напоследок табличка с абсолютной производительностью процессоров (кол-во вычисленных обратных элементов в секунду): 

Строки с префиксом P ‒ многопоточные данные.

Как я тестирую несколько потоков

В ответ Ивану.

Да вроде везде сделано одинаково. Надо бы, конечно, "причесать" и поприличней сделать, а то сплошной Copy-Paste. Но, тем не менее, мог где-то и ошибиться. Так что взгляд со стороны будет не лишним.

В управляющем потоке есть такое объявление:
  PT : array of TSingleCoreTester; // массив потоков для параллельного теста

Далее процедура для запуска этих потоков:
procedure TTestControl.TestInt64Parallel;
var
  j : integer;
begin
  PrepareMessage('Testing thread has started test PARALLEL INT64');
  if not(Terminated) then
  begin
    id := 0;
    Prefix := 'PINT64DIV';
    SetLength(PT, ProcessorCount);
    for j := 0 to High(PT) do
      PT[j] := TSingleCore_Euclid64.Create(Self, MaxInt, MaxInt + LC[id]);
    LC[id] := LC[id] * ProcessorCount;
    for j := 0 to High(PT) do
      PT[j].Start;
    Sleep(500);
    WaitForThreadsFinished;
    Synchronize(ShowIntResults);
    Synchronize(WriteIntResult);
    PrepareMessage('Testing thread has finished testing INT64DIV.');
  end;
// здесь аналогичный код для бинарного варианта
  PT := nil;
end;

WaitForThreadsFinished тупо ожидает, пока все вычислительные процессы закончатся:
procedure TTestControl.WaitForThreadsFinished;
var
  Fin : boolean;
  w : integer;
begin
  Fin := false;
  while not(Fin) and not(Terminated) do
  begin
    Fin := true;
    w := 0;
    repeat
      Fin := Fin and PT[w].Finished;
      inc(w);
    until not(Fin) or (w = ProcessorCount);
    if not(Fin) then
      Sleep(100);
  end;
  for w := 0 to High(PT) do
    FreeAndNil(PT[w]);
end;


Ну и собственно код самого потока:
procedure TSingleCore_Euclid64.Execute;
var
  T1, T2 : Int64;
begin
  QueryPerformanceCounter(T1);
  gr_InvElASM(c_pr_number64_1, St, Fin);
  QueryPerformanceCounter(T2);
  Owner.CS_TS.Acquire;
  if Owner.TS1 > T1 then
    Owner.TS1 := T1;
  if Owner.TS2 < T2 then
    Owner.TS2 := T2;
  Owner.CS_TS.Release;
  ReturnValue := 0;
end;

Здесь процедура gr_InvElASM вызывает расчет обратных элементов по модулю c_pr_number64_1 для чисел от St до Fin-1. Идея была сократить количество сохранений регистров в стек. Она на ассемблере написана, но не сильно ускорила расчет по сравнению с циклом for на Pascal.
Поэтому для бинарного варианта так делать не стал, просто в цикле Pascal вызываю функцию расчета обратных элементов. Все остальное аналогично.

Тестирую целочисленные операции

Ну уж коли я тут как-то намедни написал функции поиска обратного элемента, грех этим было не воспользоваться, что бы оценить целочисленную производительность современных процессоров. Что я собственно и сделал.
Запустил сначала однопоточный тест на своем AMD FX-4350, получил результаты и подумал, что, наверное, делать многопоточный тест смысле нет: все же чистая синтетика, 100% регистровых операций, ни одного обращения к памяти...
Но потом вспомнил про HyperThreading, CMT, SMT и вот это вот всё. Поэтому все же прогнал многопоточный тест и сильно удивился!

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

Мой родной AMD FX-4350 показал в многопотоке на чисто регистровых операциях ускорение всего-то в 1.7 раза на классическом Евклиде и в 1.4 (!) на бинарном. Честно говоря, я был в шоке! Попялился на всякий случай на код и ничего там криминального не нашел.
Надо было не тупо пялиться, а ошибку искать. 😕
Что же меня смутило? Вообще-то процессор четырехядерный (на самом деле не совсем), все операции регистровые, поэтому ускорение теоретически должно было бы равняться числу ядер, то есть четырем.
В реальности мой процессор состоит из двух модулей (то есть физически у него два ядра), но в каждом модуле по два полноценных АЛУ, которые представлены как два виртуальных ядра. Это так называемая технология CMT (Clustered Multithreading) от AMD, некий аналог HyperThreading от Intel.
Но даже если с CMT все так плохо, то все же двукратное ускорение должно же было бы получиться. Но ведь и его не выходит! Причем бинарный вариант, использующий самые простые арифметические операции: сложение и сдвиги, показывает наихудший коэффициент ускорения в многопоточном варианте!
Честно говоря, я не знаю, в чем тут дело. Складывается впечатление, что одно АЛУ вообще общее на все ядра. Поэтому одновременное исполнение максимум четырех команд возможно лишь когда активно одно виртуальное ядро. А другие ядра в это время простаивают или ждут завершение длинной операции, такой как обращение к памяти, умножение или деление.
Похоже, именно поэтому вариант поиска обратного элемента на классическом алгоритме Евклида показывает более высокое ускорение: пока один поток ждет завершение операции деления, остальные успевают сделать несколько простых вспомогательных арифметических операций.

Впрочем, такое поведение характерно не только для PileDriver. Наименьшее ускорение для бинарного варианта (1.1) показали мобильные AMD A10-4600M и Intel i3-3227U.
Впрочем, более мощные процессоры не намного лучше: Intel i5-8300H и Intel i7-6700HQ показывают ускорение 2.1 и 2.0 соответственно (по 4 физических ядра), а Ryzen 5 3600 ‒ 2.8 при 6 физических ядрах.
И все они показывают лучшее ускорение при использовании алгоритма с делением.

В связи с чем я с ностальгией вспоминаю свой старенький шестиядерный Phenom II. Похоже, это был последний процессор с относительно честной многоядерностью.
К сожалению, протестировать его нет возможности, так как несколько лет назад сгорела мать, при ее замене я недооценил степень развития собственной дальнозоркости и со слепу нечаянно спалил и Phenom II при замене материнки.
При замене на AMD FX-4350 я сразу заметил, что система стала менее отзывчивой. Я связывал это с тем, что поменялся южный мост, но, возможно дело не только в нем. Конечно, пиковая производительность у FX-4350 выше, чем у Phenom II, но вот насчет многозадачности я сейчас уже не уверен.

Ну и наконец-то собственно сама производительность. В качестве базы я, как обычно взял свой AMD FX-4350, поэтому его производительность везде принята за единицу.
Производительность сложения и сдвига я считаю в 1.5 раза более ценной, чем операции деления (так как встречаются в коде они гораздо чаще), как и производительность однопоточного выполнения.
Графики тоже удалил, так как они ошибочны.

Как можем видеть, в целочисленных регистровых операциях даже относительно новые, хоть и мобильные процессоры Intel с трудом тягаются со стареньким десктопным AMD. Ну а относительно новый Ryzen 5, теперь уже предыдущего поколения, вообще для них недостижим.
Так что и в самом деле архитектура Ryzen очень хороша. Хоть она, похоже, чуть-чуть хуже в операциях с плавающей точкой, но опережает в целочисленных, что для большинства приложений важнее. За исключением, может быть, трехмерных игр и приложений ИИ.

Посмотрим более детально. Вариант с использованием операции деления.

Если в функции встречается операция деления, то тут даже старенький процессор от AMD предпочтительней новых изделий от Intel. А Ryzen так вообще рвет все интеловские процессоры, как Тузик грелку.

Бинарный вариант на сложениях и сдвигах.

Тут Intel в последних поколениях в однопотоке немного обогнала PileDriver, но Ryzen 3000 это компенсировал и на обычных целочисленных регистровых операциях AMD и Intel сейчас на равных, но за счет большего количества ядер в общем зачете все же AMD впереди.

Продолжение про Евклида

Собственно, не рекурсивный алгоритм:

function ExtEuclid(X, Y : Int64; var A, B : Int64) : Int64;
// возвращает наибольший общий делитель X и Y, X >= Y
var
  B1, B1p, XS, YS : Int64;
  Q, T : Int64;
begin
  XS := X;
  YS := Y;
  B1p := 0;
  B1 := 1;
  while Y <> 0 do
  begin
    Q := X div Y;
    T := Y;
    Y := X - Q * Y;
    X := T;
    if Y <> 0 then
    begin
      T := B1;
      B1 := B1p - Q * B1;
      B1p := T;
    end;
  end;
  Result := X;
  B := B1;
  A := (X - B1*YS) div XS;
end;

Ну и его ассемблерный вариант:

function ExtEuclidASM(X, Y : UInt64; var A, B : Int64) : UInt64;
// возвращает наибольший общий делитель X и Y, X >= Y
asm
  push RDI; // B1Lo
  push RSI; // B1Hi
  push RBX; // Y;
  push R12; //XS;
  push R13; //YS
  push R14; // Временный
  push R15; // Временный
  mov R12, RCX;
  mov R13, RDX;
  mov RBX, RDX;
  xor R10, R10;
  xor R11, R11;
  xor RSI, RSI;
  xor RDI, RDI;
  inc DI;
@Loop:
  xor RDX, RDX;
  mov RAX, RCX;
  div RBX;
  mov RCX, RBX;
  mov RBX, RDX;
  test RDX, RDX;
  jz @FinLoop;
  mov R14, RAX
  mul RDI;
  xchg R14, RAX;
  mov R15, RDX;
  imul RSI;
  add R15, RAX;
  sub R10, R14;
  sbb R11, R15;
  xchg RDI, R10;
  xchg RSI, R11;
  jmp @Loop;
@FinLoop:
  mov R14, B;
  mov [R14], RDI;
  mov RAX, RDI;
  mul R13;
  mov R14, RAX;
  mov R15, RDX;
  mov RAX, RSI;
  imul R13;
  add R15, RAX;
  mov RAX, RCX;
  xor RDX, RDX;
  sub RAX, R14;
  sbb RDX, R15;
  idiv R12;
  mov R14, A;
  mov [R14], RAX;
  mov RAX, RCX;
  pop R15;
  pop R14;
  pop R13;
  pop R12;
  pop RBX;
  pop RSI;
  pop RDI;
end;

Бинарный вообще не стал переделывать по следующим соображениям: для моих задач нужна версия, которая работает на полном беззнаковым 64-битном целом, что требует 128-битных вычислений. Соответственно, отказ от второго if во внешнем цикле с обменом уравнений (с целью упростить условие в цикле repeat) потребует пять регистровых обменов, а в результате условие упроститься всего на 4 команды, так что особого выигрыша не будет.
Ресурсы по использованию conditional move исчерпаны: я использую все 16 регистров в теле функции, использовать же память для хранения временных значений, которые, возможно, понадобятся,  очень плохая идея с точки зрения производительности. А если вспомнить, что я использую 128-битную арифметику ‒ то просто катастрофичная.

Теперь результаты ассемблерных вариантов расширенного алгоритма Евклида (паскалевские чуть медленнее, да и сравнивать их не корректно: они используют более короткую 64-битную арифметику):

• рекурсивный ‒ 7.05 млн./сек;
• обычный ‒ 9.7 млн./сек;
• бинарный ‒ 1.47 млн./сек.
  

Видимо, что бы бинарный работал быстрее обычного, нужно что бы одна итерация бинарного выполнялась примерно в 4 раза быстрее.
Так как полный вариант деления и умножения и так уже почти полностью 128-битный у процессоров x64, то тут дополнительных издержек у обычного алгоритма нет, чего не скажешь про бинарный.
Но если потребуется выполнять те же действия с более длинными данными, например, 256-битными, то, видимо, бинарный алгоритм снова будет впереди.

Расширенный алгоритм Евклида

Ситуация сложилась так, что временно три компьютера мне не доступны для теста производительности на решении СЛАУ методом Гаусса-Жордана. Поэтому пока я жду их доступности, что бы получить результаты, и тогда проведу детальный анализ.

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

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

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

function ExtEuclidRec(X, Y : Int64; var A, B : Int64) : Int64;
var
  OldA : Int64;
begin
  OldA := X mod Y;
  if OldA = 0 then
  begin
    Result := Y;
    A := 0;
    B := 1;
  end
  else
  begin
    Result := ExtEuclidRec(Y, OldA, A, B);
    OldA := A;
    A := B;
    B := OldA - B*(X div Y);
  end;
end;

Уточню обозначения. Расширенный алгоритм Евклида фактически решает линейное диофантово уравнение Ax + By = D.
Встречаются разные трактовки переменных в этом выражении. У меня x и у ‒ заданные значения, а A и B ‒ те константы, которые находит алгоритм Евклида, D же, соответственно ‒ наибольший общий делитель x и у.
Часто встречает обратный подход, когда A и B задано, а x и y необходимо найти.
Не тестил, но, возможно при значения x, y, близких к 2^63 будет работать неверно, из-за переполнения при операциях вычитания. При меньших значения работает всегда верно.

Есть еще бинарный вариант расширенного алгоритма Евклида. Вроде как в обычном варианте, по утверждению Кнута, бинарный вариант быстрее на 60% классического. Но с того момента много воды утекло, и операции деления, хоть и остаются медленными, но стали все же намного быстрее, чем это было 20 лет назад, поэтому это утверждение требует проверки на современных компьютерах. Тем не менее, вполне допускаю это для обычного алгоритма Евклида, так как он вполне себе компактен и красив.
Но вот та версия расширенного бинарного алгоритма, что мне удалось найти, не так как красива, как обычного. Приведу ее рабочий и протестированный вариант:

function ExtEuclidBin(X, Y : Int64; var A, B : Int64) : Int64;
var
  g : int64;
  u, v : Int64;
  A1, B1, C, D : Int64;
begin
  g := 1;
  while (X and 1 = 0) and (Y and 1 = 0) do
  begin
    X := X shr 1;
    Y := Y shr 1;
    g := g shl 1;
  end;
  u := X;
  v := Y;
  A1 := 1;
  B1 := 0;
  C := 0;
  D := 1;
  repeat
    while (u and 1 = 0) do
    begin
      u := u shr 1;
      if (A1 and 1 = 0) and (B1 and 1 = 0) then
      begin
        A1 := A1 div 2;
        B1 := B1 div 2;
      end
      else
      begin
        A1 := (A1+Y) div 2;
        B1 := (B1-X) div 2;
      end;
    end;
    while (v and 1 = 0) do
    begin
      v := v shr 1;
      if (C and 1 = 0) and (D and 1 = 0) then
      begin
        C := C div 2;
        D := D div 2;
      end
      else
      begin
        C := (C+Y) div 2;
        D := (D-X) div 2;
      end;
    end;
    if u > v then
    begin
      u := u-v;
      A1 := A1-C;
      B1 := B1-D;
    end
    else
    begin
      v := v-u;
      C := C-A1;
      D := D-B1;
    end;
  until (u = 0) or (v=0);
  if u = 0 then
  begin
    Result := v*g;
    A := C;
    B := D;
  end
  else
  begin
    Result := u*g;
    A := A1;
    B := B1;
  end;
end;

Здесь такое же обозначение переменных, как и в предыдущем примере.
Это просто жесть, конечно, по сравнению с компактной рекурсивной записью расширенного алгоритма Евклида в классическом варианте. Как думаете, будет ли такой алгоритм работать быстрее классического?
Кроме сложности алгоритма, у него есть и другие особенности. В частности, линейное диофантово уравнение допускает множество решений, и классический расширенный алгоритм Евклида дает минимальный вариант коэффициентов А и B, в то время как бинарный дает произвольное правильное решение, то есть найденные А и B могут быть не минимальными по абсолютному значению.
В связи с этим у приведенного бинарного алгоритма есть более серьезное ограничение, чем у классического варианта: он гарантированно корректно работает только при x и у < 2^32. Вызвано это возможным переполнением при вычислении A1, B1, C, D. Может быть, и есть вариант, у которого данное ограничение отсутствует, но мне про него не известно.

Кстати, любопытно, но последние версии Delphi не умеют оптимизировать деление знаковых целых на два, хотя, например, Delphi 7 делала это влет.

Lazarus и AVX

Lazarus - это IDE а-ля Delphi, только под Free Pascal (FPC). Установил, на первый взгляд все было очень ничего. Внешний вид очень похож на Delphi 7.0 - самый удобный интерфейс для разработки. Тот интерфейс, что появился в более поздних версиях Delphi, лучше использовал площадь широкоформатного экрана, но был гораздо менее удобен.
Единственно, что удивило, в настройках интерфейса по умолчанию вся пунктуация (это называется "языковые символы") - красного цвета. Жесть.

К сожалению, на второй и третий взгляд все оказалось не айс. В конфах пишут, что умеет в AVX прямо с паскаля конвертировать, но слегка с глюками, то есть все операции с массивами должны быть кратны 4, иначе не компилит. Не знаю, не проверял.
Меня интересовало как спортируется код Delphi и как ассемблер встроенный там пашет. Портировалось нормально, с минимальными переделками все даже заработало. Но на паскале.

С ассемблером оказалось сложнее. Споткнулся компилятор вот на таком фрагменте: 

procedure SolveTestSingleLESP_SSE(var P : TLESInfoP; LESNum : UInt64); // P - RCX, LESNum - RDX;
asm
// R8 - P, R9w - n, R10w - index
  push RBX;
  mov R8, P.Systems;
  movzx R9, word ptr P.LESSize; // в R9w - n;
  mov R8, [R8+LESNum*8]; //
  ...

Параметры передаются в регистре RCX и RDX соответственно (это, оказывается, в Windows так принято, поэтому работает одинаково что в Delphi, что в FPC). И если с первыми строчками все нормально, то вот с последней компилятор забавно загнал. Он правильно распознал, что LESNum находится в RDX, но почему-то транслировал после замены переменной на регистр в mov R8, [R8*8+RDX].
Если переменную заменить явно на имя регистра, то транслирует нормально.
Второй затык случился на команде такого типа:

movss XMM0, dword ptr c_single1[0]; 

Почему-то компилятор выдал предупреждение, что он сгенерировал абсолютную 32-битную ссылку на константу, поэтому программа обязательно должна быть загружена в младшие четыре гигабайта ОЗУ!
Едрит-Мадрид! Это как и почему? Вроде же используется виртуальная память, которой все равно, где физически располагается переменная. И зачем специально генерировать абсолютную ссылку?
Естественно, в этом месте возникает ошибка. Наверное, это можно как-то порешать, используя какие-то специальные конструкции для ассемблера, имеющиеся в FPC. Но мне стало скучно и лень. Тем более, нормальной документации нету, что бы можно было легко докопаться до таких тонкостей.

Отдельно про debbuger. Я выяснил, что в Delphi он глючит, когда FPU (x87) занято данными более чем наполовину (просто данные верхней, самой новой половины стека не показывает). В Lazarus такого глюка нет. Просто потому, что там нет ресурсов FPU в отладчике как такового.
Ну и вообще там очень хиленькие и неудобные возможности по отладке в ассемблере по сравнению с Delphi. Регистров AVX там, кстати, тоже нет.
Синтаксис AT&T ассемблера отдельно достает. Нет, если впитал его с молоком матери, то нормально. Но у меня обратная ситуация. И если при разработке можно переключить ассемблер на синтаксис Intel, то отладчик всегда показывает только AT&T синтаксис.
В целом Lazarus произвел впечатление полудоделанного и слегка подзаброшенного продукта.

В общем, плюнул я на это дело и решил попробовать запасной вариант: написать нужные мне процедуры на чистом ассемблере.
Идеальный для меня вариант - TASM. Но беда в том, что Embarcadero его давно забросил. И поэтому он поддерживает ровно тот же набор команд, что и компилятор Delphi: то есть AVX там нет.
Почитал в инете, вроде прямым наследником TASM, который более-менее развивается, является FASM. Это как бы и так, но не совсем. FASM настолько далеко ушел в своем развитии, что я не смог сгенерировать объектный файл с самыми простыми процедурами. Отказывается последняя версия FASM делать такой объектный файл. Хотя вроде долго смотрел в примеры, писал все по образцу, но не получается. Правда, не догадался попробовать сами примеры скачать и скомилировать... Надо будет как-нибудь ради прикола сделать.
Короче, заменил процедуры метками, благо, у меня все передается в регистрах, а локальные переменные мне не нужны. Видимо, если бы были нужны, надо было стек вручную создавать? Ну да бог с ним.

Быстренько скопировал код SSE, адаптировал его под AVX и на типе single (float), удивительное дело, все правильно сработало с первого раза, без отладки. Но, ради интереса, решил посмотреть под отладчиком. И не смог: при попытке зайти в процедуру с кодом под AVX отладчик Delphi просто "умирает".
Так что отладки под AVX у меня нету. С кодом под double переход был не таким гладким, в паре мест накосячил, пришлось повозится, выискиваю баги без отладчика. Но тем не менее все получилось.

Теперь осталось собрать тестер и оттестировать доступные мне процессоры. Но предварительно уже могу сказать: на PileDriver код на AVX немного медленнее, чем на SSE. Впрочем, это ожидаемо: давно известно, что PileDriver очень медленно выполняется обращение к памяти по 256 бит. А большинство команд AVX именно такие.
И даже команды FMA PileDriver не помогли. Работает ровно с той же скоростью, что и просто AVX. Но точность чуть выше, чем отдельно умножение + сложение. Впрочем, для метода Гаусса-Жордана повышение точности не очень существенное.
Наверное, если использовать команды AVX на половинных регистрах (XMM вместо YMM), то может быть AVX слегка и обгонит SSE на моем FX-4350, особенно если FMA использовать. Но, мне кажется, это будет не совсем честно по отношению к другим процессорам 😉.

Маленький размер

Заболел. В перерывах постельного режима решил ускорить работу СЛАУ как раз самого маленького размера: в 3 уравнения. Все же похожие системы очень часто встречаются в реальности. Практически все трехмерные расчеты так или иначе связаны примерно с этим размером. В то время как решение очень больших СЛАУ - это скорее очень редкие исключения. По крайней мере мне так кажется.

 

Начал, естественно, с FPU (x87). У него есть 8 вещественных регистров. А коэффициентов в СЛАУ 12. Впрочем, точно хранить нужно только элементы над верхней диагональю и свободные члены, которых всего 7. Но еще нужны свободные регистры для временных операций...

Короче, программирование матричных операций на регистрах, завязанных в стек ‒ та еще радость! И самое главное, не очень ясно, велик ли выигрыш от хранения коэффициентов в регистрах.

Попутно полдня мучил свой воспаленный мозг, пытаясь понять, что же не так то!? Через несколько часов выяснил, что "не так" - это глюк в дебаггере Delphi, который по крайней мере в 64-битном режиме неверно отображает заполненный более чем на половину регистровый стек FPU!

Да, видимо придется с Delphi все же уходить, как не жаль. Тем более, что Lazarus (точнее, Free Pascal) умеет вроде даже в AVX.

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

SSE

Ну вот, наконец-то добрался до SSE. И предварительные результаты несколько обескураживающие. Я, конечно не ожидал, что использование векторных команд увеличит производительность пропорционально размеру вектора, все же даже решение СЛАУ, которое прямо просятся решаться векторными операциями, на 100% не векторизуется.

Но все же я надеялся хотя бы на двукратный рост на типе single (float), которому на SSE соответствует вектор из 4 элементов. На практике же максимальное ускорение относительно решения на Scalar SSE составило всего 1.66 раза, что на AMD FX-4350 соответствует решению 158 уравнений в секунду.
Очевидно, что на маленьких СЛАУ вышло даже небольшое замедление. На СЛАУ из 3 уравнений ускорение получилось 0.93.

Но самая печаль на типе double. В SSE ему соотвествует вектор всего из двух элементов. Максимальное ускорение достигнуто на СЛАУ из 192 переменных и составляет 1.18 раза. На маленьких размерностях (3 и 6 уравнений) ускорение меньше единицы, далее растет, после максимума на 192 уравнениях плавно снижается, приближаясь к единице 1536 уравнениях.

Одна из возможных причин состоит в том, что SSE не умеет выполнять арифметические операции с невыровненными на границу 16 байт данными. А при решении СЛАУ они гарантированно не выровнены.
AMD рекомендует в таких случая выполнять сохранение в две команды, записывая за раз по половине векторного регистра. Попробовал, ощутимой разницы не увидел: на каких-то размерностях стало быстрее, на других же наоборот.
Вроде бы в AVX эту "недоработку" пофиксили и арифметические команды выполняют и над невыровненными данными.
Тут кстати обнаружил, что моя версия Delphi не умеет в AVX. Надо подумать, как красивее обойти это органичение.

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

procedure MakeDoubleNextRowsP_SSE;
asm
// R8 - P, R9w - n, R10w - index
  lea EBX, [R10d+1]; // BX - i
  cmp BX, R9w;
  jae @exit;
@Loop1:
  mov RCX, [R8 + RBX*8]; // RCX - S
  lea EDX, [R10d+1];
  lea EAX, [R9d-7];
  movsd XMM0, [RCX+R10*8];
  unpcklpd XMM0, XMM0;
  cmp EDX, EAX;
  jg @Loop2_4;
@Loop2:
  movapd XMM1, XMM0;
  movupd XMM2, [R11+RDX*8];
  mulpd XMM1, XMM2;
  movupd XMM2, [RCX+RDX*8];
  subpd XMM2, XMM1;
  movupd [RCX+RDX*8], XMM2;
  movapd XMM1, XMM0;
  movupd XMM2, [R11+RDX*8+16];
  mulpd XMM1, XMM2;
  movupd XMM2, [RCX+RDX*8+16];
  subpd XMM2, XMM1;
  movupd [RCX+RDX*8+16], XMM2;
  movapd XMM1, XMM0;
  movupd XMM2, [R11+RDX*8+32];
  mulpd XMM1, XMM2;
  movupd XMM2, [RCX+RDX*8+32];
  subpd XMM2, XMM1;
  movupd [RCX+RDX*8+32], XMM2;
  movapd XMM1, XMM0;
  movupd XMM2, [R11+RDX*8+48];
  mulpd XMM1, XMM2;
  movupd XMM2, [RCX+RDX*8+48];
  subpd XMM2, XMM1;
  movupd [RCX+RDX*8+48], XMM2;
  add DX, 8;
  cmp EDX, EAX;
  jle @Loop2;
@Loop2_4:
  add EAX, 4;
  cmp EDX, EAX;
  jg @Loop2_2;
  movapd XMM1, XMM0;
  movupd XMM2, [R11+RDX*8];
  mulpd XMM1, XMM2;
  movupd XMM2, [RCX+RDX*8];
  subpd XMM2, XMM1;
  movupd [RCX+RDX*8], XMM2;
  movapd XMM1, XMM0;
  movupd XMM2, [R11+RDX*8+16];
  mulpd XMM1, XMM2;
  movupd XMM2, [RCX+RDX*8+16];
  subpd XMM2, XMM1;
  movupd [RCX+RDX*8+16], XMM2;
  add DX, 4;
@Loop2_2:
  add EAX, 2;
  cmp EDX, EAX;
  jg @Loop2_1;
  movapd XMM1, XMM0;
  movupd XMM2, [R11+RDX*8];
  mulpd XMM1, XMM2;
  movupd XMM2, [RCX+RDX*8];
  subpd XMM2, XMM1;
  movupd [RCX+RDX*8], XMM2;
  add DX, 2;
@loop2_1:
  cmp DX, R9w;
  ja @Loop2_Fin;
  movapd XMM1, XMM0;
  mulsd XMM1, [R11+RDX*8];
  movsd XMM2, [RCX+RDX*8];
  subsd XMM2, XMM1;
  movsd [RCX+RDX*8], XMM2;
@Loop2_Fin:
  inc BX;
  cmp BX, R9w;
  jb @Loop1;
@exit:
end;
 

Продолжаю оптимизировать

Точнее, пытаюсь закончить оптимизировать Scalar SSE в тех местах, которые до конца не проверил, или проверил не правильно. 

1. Замена вычитания сложением. В чем суть этой попытки оптимизации? У нас есть вот такое выражение:
  S[j] := S[j]-R[j]*S[index];
На ассемблере я его транслировал в такой набор команд
  movapd XMM1, XMM0;
  mulsd XMM1, [R11+RDX*8];
  movsd XMM2, [RCX+RDX*8];
  subsd XMM2, XMM1;
  movsd [RCX+RDX*8], XMM2;
в двух операторной форме команд, принятой в x86, я вынужден сначала загрузить S[j] в регистр, а уже потом выполнять вычитание.
Если же у S[index] перед внутренним циклом поменять знак, то операцию вычитания можно заменить сложением без предварительной загрузки из памяти в регистр XMM2. В результате получится вот такой код:
  movapd XMM1, XMM0;
  mulsd XMM1, [R11+RDX*8];
  addsd XMM1, [RCX+RDX*8];
  movsd [RCX+RDX*8], XMM1;
Как видим, он на одну команду короче, и, возможно, быстрее. Не стоит забывать, что нам пришлось во внешней цикле добавить команду изменения знака S[index].
Но будет ли этот код реально быстрее? Зависит от того, какая команда выполняется быстрее: умножение XMM1 на значение в память или загрузка в XMM2 значения из памяти. Дело в том, что команда addsd не может начаться, пока не закончится умножение XMM1.
Поэтому, если умножение медленное, то XMM2 успеет загрузиться из памяти (так как эти команды независимы) и затем можно будет вычесть из XMM2 XMM1 без потери производительности. То есть в этом случае замена вычитания сложением не нужна, так как ничего не даст.
В том случае, если умножение выполняется быстрее, чем загрузка XMM2, то ситуация будет обратная: выгоднее будет второй вариант со сложением.
На мой взгляд, разницы тут не должно быть много, но давайте посмотрим на результаты.
На типе single (float) ускорение от замены вычитания на сложение составило от сотых долей процента на системах из 3 уравнений, до 4-5% на больших системах до 48 переменных включительно. А вот на больших размерностях результат оказался парадоксальный: незначительное снижение производительности на системах из 96 уравнений и потери до 3.5% и выше на больших размерностях.
На типе double ускорение не превышает 3%, достигая максимума на 24 уравнениях, затем постепенно снижается и заменяется замедлением на системах в 1536 уравнений.
В общем, я не увидел особого смысла в такой оптимизации.

2. Развертывание цикла. Я прошлый раз ограничился четырьмя операциями, но Иван предложил попробовать увеличить их до восьми.
Результат моих размышлений на эту тему.
Ускорение от развертывания можно получить за счет снижения доли служебных операций на организацию цикла. Но это работает только в том случае, когда тело цикла (рабочая нагрузка) очень проста и по времени сравнима с работой по организации самого цикла.
В моем случае сам цикл организуется всего тремя командами: двумя короткими регистровыми операциями и условным переходом. А полезной нагрузкой является 5 команд, в том числе три с памятью. В принципе, доля служебных команд достаточна велика и такой цикл имеет смысл разложить на две итерации. Дальнейшее разложение приведет лишь к незначительному росту производительности, тем более, что хорошо работает блок предсказаний переходов, команды организации цикла и полезной нагрузки исполняются на разных блоках, так что частично могут выполняться параллельно.
Второй источник ускорения связан с отсутствием сбросов конвейера. Чем длиннее линейная часть программы, без условных переходов, тем полнее заполнен конвейер и выше производительность процессора. Наличие условных переходов сбрасывает конвейер,что приводит к сильному снижению производительности.
Вышеописанная страшилка более актуальная для очень старых процессоров. Современные версии процессоров имеют очень хороший блок предсказания переходов, который работает совместно с процессором и остановка конвейера происходит только в том случае, если блок предсказаний ошибся. В частности, для for-циклов блок предсказаний ошибается очень редко, в идеальном случае лишь один раз ‒ когда цикл заканчивается. Думаю, что именно этот источник ускорения при решении моей задачи имеет очень низкое значение.
Третий источник ускорения ‒ суперскалярная архитектура, причем именно она, на мой взгляд, дает наибольший эффект в моем случае. Посмотрим на полезную нагрузку цикла. В случае одной операции на итерацию у нас требуется выполнить одно умножение и одно сложение, причем сложение зависит от умножение. Фактически здесь нет операций, которые могли бы выполнятся параллельно.
Если выполнять две операции за итерацию, то нужно выполнить 2 умножения и 2 сложения, тут уже есть что распараллелить. 2 умножения можно выполнять параллельно, но со сложениям из-за зависимости от умножений ситуация сложнее, скорее всего параллельно умножениям удастся выполнить лишь одну команду сложения.
При увеличении количества операций на итерацию ситуация с распараллеливанием становится все лучше и лучше. Теоретически могут выполняться 4 инструкции за такт, но практически, есть подозрение, их количество лежит в пределах от 2 до 3, причем ближе к двум, по всей видимости. Это связано с тем, что умножение все же выполняется медленнее сложения, поэтому для сумматоров данные всегда поставляются с задержкой.
Тем не менее, попробуем оценить, насколько увеличивается производительность при переходе от 4 операций на итерацию к 8.
На типе single (float) ускорение составило от 0% (3 переменных) до 11% (192 переменных), затем уменьшилось до 4% на 384 переменных.
На типе double ускорение максимум составило около 7% на 192 переменных, на всех остальных размерах ускорение составляет около 1%, что ниже погрешности измерений.
В результате решил оставить развертывание на восемь операций для типа single (float), а для типа double остановиться на четырех.

3. Решил побаловаться с командами PREFETCH. Эффект для моей задачи равен нулю, так как, видимо, механизм автоматической предзагрузки кэша в этом случае справляется на ура.

Оптимизирую Scalar SSE

Первый вариант моей оптимизации путем развертывания цикла оказался не очень удачным, как совершенно справедливо указал Иван Колесников в комментариях к посту Scalar SSE. Так что исправляю.

1. Да, организация цикла были весьма корява, исправление существенно увеличило производительность.

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

3. Ручное перемешивание кода, опять же, не давало нормально работать механизму внеочередного исполнения команд. После того, как поставил команды в естественном порядке, код стал работать быстрее.

4. А вот предложение заменить вычитание умножением на -1 и сложением наоборот, замедляло производительность.

5. Не стал реализовать концовку аналогом оператора case (switch). Вот мои соображения:
во-первых, такой оператор не очень быстр: тут либо несколько условных переходов, либо, более быстро косвенный безусловный переход по адресу в памяти (8-ми байтному); не уверен, что адрес будет в кэше, предыдущий длинный цикл по строке скорее всего его вытеснит; и переход такой будет лишь один раз на все итерации внутреннего цикла;
во-вторых, такой вариант не очень хорошо ложится на последующий переход к векторному варианту решения.
Спорно, конечно, но решил что этот вариант не сильно ускорит решение СЛАУ.

Насколько удалось ускорить, покажут будущие тесты.

Вскоре возьмусь за векторный вариант, для начала на SSE.

Ну и напоследок код главного цикла (цикл для верхней строки тоже аналогично оптимизировал; ну и то же самое сделал для FPU):

procedure MakeDoubleNextRowsPASM;
asm
// R8 - P, R9w - n, R10w - index
  lea EBX, [R10d+1]; // BX - i
  cmp BX, R9w;
  jae @exit;
@Loop1:
  mov RCX, [R8 + RBX*8]; // RCX - S
  lea EDX, [R10d+1];
  mov EAX, R9d;  // --
  movsd XMM0, [RCX+R10*8];
  sub EAX, 3;
  cmp EDX, EAX; // --
  jg @Loop2_2;
@Loop2:
  movapd XMM1, XMM0;
  mulsd XMM1, [R11+RDX*8];
  movsd XMM2, [RCX+RDX*8];
  subsd XMM2, XMM1;
  movsd [RCX+RDX*8], XMM2;
  movapd XMM1, XMM0;
  mulsd XMM1, [R11+RDX*8+8];
  movsd XMM2, [RCX+RDX*8+8];
  subsd XMM2, XMM1;
  movsd [RCX+RDX*8+8], XMM2;
  movapd XMM1, XMM0;
  mulsd XMM1, [R11+RDX*8+16];
  movsd XMM2, [RCX+RDX*8+16];
  subsd XMM2, XMM1;
  movsd [RCX+RDX*8+16], XMM2;
  movapd XMM1, XMM0;
  mulsd XMM1, [R11+RDX*8+24];
  movsd XMM2, [RCX+RDX*8+24];
  subsd XMM2, XMM1;
  movsd [RCX+RDX*8+24], XMM2;
  add DX, 4;
  cmp EDX, EAX; // --
  jle @Loop2;
@Loop2_2:
  add EAX, 2;
  cmp EDX, EAX;
  jg @Loop2_1;
  movapd XMM1, XMM0;
  mulsd XMM1, [R11+RDX*8];
  movsd XMM2, [RCX+RDX*8];
  subsd XMM2, XMM1;
  movsd [RCX+RDX*8], XMM2;
  movapd XMM1, XMM0;
  mulsd XMM1, [R11+RDX*8+8];
  movsd XMM2, [RCX+RDX*8+8];
  subsd XMM2, XMM1;
  movsd [RCX+RDX*8+8], XMM2;
  add DX, 2;
@loop2_1:
  cmp DX, R9w; //--
  ja @Loop2_Fin;
  movapd XMM1, XMM0;
  mulsd XMM1, [R11+RDX*8];
  movsd XMM2, [RCX+RDX*8];
  subsd XMM2, XMM1;
  movsd [RCX+RDX*8], XMM2;
@Loop2_Fin:
  inc BX;
  cmp BX, R9w;
  jb @Loop1;
@exit:
end;

Scalar SSE. Предварительные результаты.

Оттестировал разные процессоры на пока еще не до конца оптимизированном коде. Думаю, после оптимизации он станет существенно быстрее. Но так как оптимизация будет инвариантна для всех процессоров, то относительный результат не сильно изменится.

Во-первых, развертывание циклов даже не в самом оптимальном варианте показало все же более серьезное ускорение, особенно на больших размерностях.
У AMD FX-4350 это порядка 30% для single (float) и 36% для double. Мобильный AMD A10-4600M показывает чуть худший результат: 20% и 33% соответственно. Рекордсмен среди процессоров от AMD, конечно Ryzen 5, у него 41% и 67%. Только надо учитывать, что это относительно производительности не развернутого цикла FPU (x87), который в Ryzen, как выяснил Иван (см. комментарии к прошлому посту), зарезан в 2 раза по отношению к Scalar SSE, тут все четко получилось.
Intel i3-3227U показывает 29% и 27%. Старшие Intel тут менее результативны: Intel i7-6700HQ 20% для обоих типов и Intel i5-8300H 18% и 17%.

Теперь сравнение производительности Scalar SSE для разных процессоров относительно AMD FX-4350.

 

Как видим, появился тест, в котором AMD Ryzen 5 3600 смог догнать интеловские процессоры. Правда, радость омрачает тот факт, что догнать он смог только мобильные процессоры.
Люблю я процессоры AMD, в первую очередь, конечно, за цены. Еще они первыми добавили поддержку PCIe 4.0, что тоже плюс, правда, не очень пока ясно с чем его кушать, этот плюс.
Но справедливости ради должен сказать, что по однопоточной производительности они все еще далеки до отчетливого лидерства. Но надо проверить и другие дисциплины.😉

Еще бросается в глаза тот факт, что на младших мобильных, что Intel, что AMD, в многозадачном режиме SSE приходится жить похуже, а вот старшие мобильные Intel гораздо лучше с этим справляются. Даже лучше, чем с FPU. На Ryzen многопоточный режим одинаково эффективен что на FPU, что на Scalar SSE, ну и банально: этот процессор однозначный лидер в многопоточной производительности тупо из-за числа ядер.

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

Scalar SSE

Начал реализовывать. В принципе, на что рассчитывал? FPU (x87) осуществляет операции, расширяя предварительно операнды в памяти до 80 бит, в то время как Scalar SSE этого вроде бы не делает, выполняя операции точно в размер типа.
Соответственно, на типе single (float) в 32 бита теоретически можем получить ускорение в 2.5 раза. Заменяем команды FPU на соответствующие команды Scalar SSE и получаем... ту же самую производительность, по крайней мере на AMD FX-4350. То есть производительность FPU и Scalar SSE отличается на плюс/минус погрешность измерения.
Хм... Крайне странно. Попробую что ли развернуть циклы?

Лирическое отступление про развертывание. Оптимизация здесь достигается за счет того, что не сбрасывается конвейер и полностью используются суперскалярные возможности процессора. В полной мере воспользоваться преимуществом развертывания легко на коротких циклах.
Тем не менее, на современных процессорах влияние развертывания на производительность заметно меньше, в первую очередь за счет очень качественного механизма предсказания переходов, а также использования теневых регистров.
И еще меньше оно становится, если параллельно используются два совершенно независимых устройства, в нашем случае блоки целочисленных операции и операций с плавающей точкой.
Ведь пока выполняется относительно длинная операция с плавающей точкой, мы можем параллельно выполнить пару коротких целочисленных операций, практически без потери производительности.
Еще один нюанс связан с количеством итераций цикла. Если количество итераций заранее известно, то можно цикл развернуть максимально оптимальным способом. Если же такой информации нет и пишется универсальный цикл, то при развертывании придется за счет дополнительных условий пожертвовать производительностью коротких циклов.

Теперь вернемся к решению СЛАУ на Scalar SSE. Документация AMD говорит о том, что при использовании SSE можно рассчитывать на увеличение производительности до 4 раз.
От чего это зависит? По моим соображениям, от количества АЛУ в блоке SSE. То есть 4-х кратного ускорения можно достигнуть, если у нас есть 4 сумматора, например. Но только на операциях сложения.
В случае наличия по два мультипликатора и сумматора можно достичь двукратного ускорения. Или 4-х кратного, если у нас есть данные, которые позволяют параллельно выполнить два умножения и два сложения.
При таком построение блока SSE возможно получить существенный выигрыш и от раскрытия цикла Scalar SSE в случае суперскалярной архитектуры.

В принципе, возможно и другое аппаратное решение: одно АЛУ, которое выполняет действие над всем регистром SSE, в результате которого получается сразу четыре результата (для типа single). В этом случае раскрытие цикла на Scalar SSE не даст слишком большого выигрыша.

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

В общем, такое ощущение, что Scalar SSE ничего не знает про суперскалярность! 😉 Или я просто не умею его готовить?
Ниже пример кода, реализующего внутренний цикл, на этот раз, правда, для типа double. В принципе, ничем от single не отличается, кроме размера операндов. Результаты тоже аналогичны.

procedure MakeDoubleNextRowsPASM;
asm
// R8 - P, R9w - n, R10w - index
  lea EBX, [R10d+1]; // BX - i
  cmp BX, R9w;
  jae @exit;
@Loop1:
  mov RCX, [R8 + RBX*8]; // RCX - S
  lea EDX, [R10d+1];
  movsd XMM0, [RCX+R10*8];
  mov ax, R9w;  // --
  sub ax, R10w; // --
@Loop2:
  cmp ax, 4;
  jb @Loop2_2;
  movsd XMM1, XMM0;
  movsd XMM2, [RCX+RDX*8];
  movsd XMM3, XMM0;
  mulsd XMM1, [R11+RDX*8];
  movsd XMM5, XMM0;
  movsd XMM4, [RCX+RDX*8+8];
  subsd XMM2, XMM1;
  mulsd XMM3, [R11+RDX*8+8];
  movsd XMM7, XMM0;
  movsd XMM6, [RCX+RDX*8+16];
  movsd [RCX+RDX*8], XMM2;
  subsd XMM4, XMM3;
  mulsd XMM5, [R11+RDX*8+16];
  movsd [RCX+RDX*8+8], XMM4;
  subsd XMM6, XMM5;
  mulsd XMM7, [R11+RDX*8+24];
  movsd XMM8, [RCX+RDX*8+24];
  movsd [RCX+RDX*8+16], XMM6;
  subsd XMM8, XMM7;
  sub ax, 4
  movsd [RCX+RDX*8+24], XMM8;
  add dx, 4;
  jmp @Loop2_Fin;
@Loop2_2:
  cmp ax, 2;
  jb @Loop2_1;
  movsd XMM1, XMM0;
  movsd XMM2, [RCX+RDX*8];
  movsd XMM3, XMM0;
  mulsd XMM1, [R11+RDX*8];
  movsd XMM4, [RCX+RDX*8+8];
  subsd XMM2, XMM1;
  mulsd XMM3, [R11+RDX*8+8];
  movsd [RCX+RDX*8], XMM2;
  subsd XMM4, XMM3;
  sub ax, 2;
  movsd [RCX+RDX*8+8], XMM4;
  add dx, 2;
  jmp @Loop2_Fin;
@loop2_1:
  movapd XMM1, XMM0;
  mulsd XMM1, [R11+RDX*8];
  movsd XMM2, [RCX+RDX*8];
  subsd XMM2, XMM1;
  movsd [RCX+RDX*8], XMM2;
  inc DX; //--
@Loop2_Fin:
  cmp DX, R9W;
  jbe @Loop2;
  inc BX;
  cmp BX, R9w;
  jb @Loop1;
@exit:
end;

Многопоточное ускорение

Небольшой анализ, насколько хорошо ускоряется решение множества СЛАУ в параллельном режиме. Сначала график.

Мелкая деталь, которая обращает на себя внимание: AMD FX-4350 и A10-4600M обеспечивают почти одинаковое ускорение с небольшим преимуществом десктопного процессора. Что не удивительно ‒ в основе лежит одинаковая архитектура PileDriver. Видимо, небольшой проигрыш определяется отсутствием кэша третьего уровня у мобильной версии.
Тем не менее, AMD здесь постаралась на славу, если вспомнить, что на 4 ядра в этой архитектуре имеется всего два полноценных АЛУ с плавающей точкой. Несмотря на это, в параллельном режиме обеспечивается трехкратное ускорение!

Но еще более впечатляющих результатов AMD добилась в AMD Ryzen 5 3600. Несмотря на 6 физических ядер, на небольших размерностях ускорение превышает 6 раз. Старшие мобильные Intel вроде тоже так умеют, но на более меньших размерностях, которые в этот раз в тест не попали.

Есть и общая закономерность, которую демонстрируют все процессоры: замедление ускорения от параллельного решения СЛАУ при росте размерности задачи. Особенно это заметно у десктопного Ryzen 5 3600, а также мобильных Intel i5-8300H и Intel i7-6700HQ. Причем последний начинает сдуваться даже раньше, на 384 переменных.
В итоге ускорение параллельной обработки у мобильных процессоров Intel становится меньше единицы, т. е. один поток выполнит работу быстрее, чем восемь!
Десктопные процессоры не демонстрируют такую особенность, видимо, это произойдет на одной из следующих размерностях, которые я не тестировал.

Объясняется эта закономерность несоответствием скорости процессора и пропускной способностью памяти. При решении СЛАУ каждая операция слишком проста, и когда строки памяти становятся слишком длинные, то каждое ядро слишком часто обращается к памяти, увеличивается промахи мимо кэша, а пропускной способности шины памяти становится недостаточно, что бы успеть обслужить запросы всех ядер.
Видимо, при решение подобных задач, обрабатывающих простым алгоритмом огромный объем памяти, использование технологии HyperThreading выглядит не очень оправданным. Т. е. ее видимо, либо необходимо отключить, либо просто использовать в два раза меньше вычислительных потоков, чем ядер.
Впрочем, это не гарантирует, что на очень больших размерностях этого будет достаточно. Все зависит от скорости шины памяти. Если она достаточна медленна, а обработка одного элемента слишком быстра, то количество потоков придется еще больше ограничить.

Кстати, самые слабенькие процессоры из теста, Intel i3-3227U и AMD A10-4600M демонстрируют крайне слабое падение ускорения в зависимости от размерности задачи. Они просто достаточно медленные, и ядер не так много, поэтому подсистема памяти успевает обслуживать их запросы.

Новые результаты

Оттестировал новой 64-битной версией программки FPU (x87) на нескольких компьютерах.
Попутно понял, что не совсем правильно построил логику работы теста. Сейчас у меня, если памяти не хватает, то тест на этой размерности пропускается и берется, следующая, большая.
Логика тут простая: при увеличении размерности в 2 раза требуемый объем памяти увеличивается в 4 раза, а производительность падает (теоретически) в 8 раз, то есть при увеличении размерности нужно меньше СЛАУ, что бы провести тест.
Если же использовать только доступную память, то время решения будет слишком маленьким, а погрешность теста ‒ высокая.

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

Параллельно этим размышлениям появилась следующая мысль: видимо, при покупке нового компьютера, если есть цель обеспечить максимальную производительность современного процессора, нужно много памяти. Иначе процессор будет простаивать, ожидая загрузки страниц виртуальной памяти с внешнего ЗУ. Конечно, я тут не имею в виду процессорные обрубки для офисных пиш. машинок.
Минимум сейчас требуется, кмк, 16 ГиБ, рекомендуется 32 ГиБ, ну а если брать топовые процессоры с количеством ядер 8 и больше, то тут нужно все 64 ГиБ.
Для рабочих станции все цифры, мне кажется, надо увеличить минимум вдвое.

А теперь итоги теста.

Как видно, старенький мобильный i3-3227U с 2 ядрами выступает вровень с AMDшным A10-4600M с четырьмя ядрами. Более новый мобильный i5-8300H обгоняет i7-6700HQ какого-то более старого поколения.
А вот свежий бюджетный AMD Ryzen 5 3600 откровенно разочаровывает, так как новый десктопный процессор проигрывает обоим более пожилым мобильным конкурентам от Intel. Всегда почему-то думал, что десктопный процессор должен быть быстрее аналогичного по позиционированию мобильного.

После массы новостей, заполонивших тематические сайты, у меня создалось впечатление, что все же Ryzen должен выступать как минимум вровень с процессорами от Intel. Но нет. По крайней мере по FPU он явно проигрывает.
Правда, есть робкая надежда, что AMDшные процессоры хороши в каких-то иных дисциплинах. Проверю это чуть позже. Ну и понятно, что они точно дешевле Intel, а кроме того, у них есть PCIe 4.0.

Теперь посмотрим на многоядерную производительность в параллельном режиме.

Вот здесь процессор AMD наконец-то продемонстрировал на что он способен. Конечно, он одержал уверенную победу за счет большего количества ядер (6 физических против 4 у мобильных i5 и i7). Надо сказать, что AMD в Ryzen существенно улучшила многозадачность. Более подробно я напишу об этом в отдельном посте.

Интересно посмотреть на аутсайдеров. У i3 2 физических ядра и 4 виртуальных, которые дает Hyper Threading. У A10 вроде как 4 физических ядра и никакого HT. Но сделан A10 по технологии Piledriver, в которой лишь одно FPU на каждые два ядра (кстати, как и у FX-4350). Поэтому эти два процессора, i3 и A10, фактически и выступают на равных в многопоточном тесте, AMD даже оказался чуть-чуть быстрее. А extended (long double), в котором у Intel традиционно преимущество (видимо, память быстрее работе при доступе к невыровненным данным) я в этот раз не тестирую.

 

Ну и наконец диаграмма интегральной производительности.