Самый эффективный способ вычесть один массив из другого

Question

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

var 
    Array1, Array2 : array of integer; 

..... 
// Code that fills the arrays 
..... 

for ix := 0 to length(array1)-1 
    Array1[ix] := Array1[ix] - Array2[ix]; 

end; 

Есть ли у кого-нибудь предложение?

Answer 1

Мне очень понравилась оптимизация скорости в этом простом случае. Итак, я сделал 6 простых процедур и измерил тик процессора и время при размере массива 100000;

1. процедура Паскаль с диапазоном опций компилятора и переполнением Проверки на
2. Pascal процедуры с Range Опции компилятора и перелив Проверки от
3. Классической процедуры x86 ассемблера.
4. Процедура ассемблера с инструкциями SSE и несимметричный 16-байтовый ход.
5. Процедура ассемблера с инструкциями SSE и выровнено по 16 байт.
6. Ассемблер 8 раз развернутый контур с инструкциями SSE и выровненный 16-байтовый ход.

Проверьте результаты на картинке и код для получения дополнительной информации.

Чтобы получить выравнивание памяти 16 байт первого Delite точки в файле директиве 'FastMM4Options.inc' {$ .define Align16Bytes} !

program SubTest; 

{$APPTYPE CONSOLE} 

uses 
//In file 'FastMM4Options.inc' delite the dot in directive {$.define Align16Bytes} 
//to get 16 byte memory alignment! 
    FastMM4, 
    windows, 
    SysUtils; 

var 
    Ar1     :array of integer; 
    Ar2     :array of integer; 
    ArLength    :integer; 
    StartTicks   :int64; 
    EndTicks    :int64; 
    TicksPerMicroSecond :int64; 

function GetCpuTicks: int64; 
asm 
    rdtsc 
end; 
{$R+} 
{$Q+} 
procedure SubArPasRangeOvfChkOn(length: integer); 
var 
    n: integer; 
begin 
    for n := 0 to length -1 do 
    Ar1[n] := Ar1[n] - Ar2[n]; 
end; 
{$R-} 
{$Q-} 
procedure SubArPas(length: integer); 
var 
    n: integer; 
begin 
    for n := 0 to length -1 do 
    Ar1[n] := Ar1[n] - Ar2[n]; 
end; 

procedure SubArAsm(var ar1, ar2; length: integer); 
asm 
//Length must be > than 0! 
    push ebx 
    lea ar1, ar1 - 4 
    lea ar2, ar2 - 4 
@Loop: 
    mov ebx, [ar2 + length * 4] 
    sub [ar1 + length * 4], ebx 
    dec length 
    jnz @Loop 
@exit: 
    pop ebx 
end; 

procedure SubArAsmSimdU(var ar1, ar2; length: integer); 
asm 
//Prepare length 
    push  length 
    shr  length, 2 
    jz  @Finish 
@Loop: 
    movdqu xmm1, [ar1] 
    movdqu xmm2, [ar2] 
    psubd xmm1, xmm2 
    movdqu [ar1], xmm1 
    add  ar1, 16 
    add  ar2, 16 
    dec  length 
    jnz  @Loop 
@Finish: 
    pop  length 
    push  ebx 
    and  length, 3 
    jz  @Exit 
//Do rest, up to 3 subtractions... 
    mov  ebx, [ar2] 
    sub  [ar1], ebx 
    dec  length 
    jz  @Exit 
    mov  ebx, [ar2 + 4] 
    sub  [ar1 + 4], ebx 
    dec  length 
    jz  @Exit 
    mov  ebx, [ar2 + 8] 
    sub  [ar1 + 8], ebx 
@Exit: 
    pop  ebx 
end; 

procedure SubArAsmSimdA(var ar1, ar2; length: integer); 
asm 
    push  ebx 
//Unfortunately delphi use first 8 bytes for dinamic array length and reference 
//counter, from that reason the dinamic array address should start with $xxxxxxx8 
//instead &xxxxxxx0. So we must first align ar1, ar2 pointers! 
    mov  ebx, [ar2] 
    sub  [ar1], ebx 
    dec  length 
    jz  @exit 
    mov  ebx, [ar2 + 4] 
    sub  [ar1 + 4], ebx 
    dec  length 
    jz  @exit 
    add  ar1, 8 
    add  ar2, 8 
//Prepare length for 16 byte data transfer 
    push  length 
    shr  length, 2 
    jz  @Finish 
@Loop: 
    movdqa xmm1, [ar1] 
    movdqa xmm2, [ar2] 
    psubd xmm1, xmm2 
    movdqa [ar1], xmm1 
    add  ar1, 16 
    add  ar2, 16 
    dec  length 
    jnz  @Loop 
@Finish: 
    pop  length 
    and  length, 3 
    jz  @Exit 
//Do rest, up to 3 subtractions... 
    mov  ebx, [ar2] 
    sub  [ar1], ebx 
    dec  length 
    jz  @Exit 
    mov  ebx, [ar2 + 4] 
    sub  [ar1 + 4], ebx 
    dec  length 
    jz  @Exit 
    mov  ebx, [ar2 + 8] 
    sub  [ar1 + 8], ebx 
@Exit: 
    pop  ebx 
end; 

procedure SubArAsmSimdAUnrolled8(var ar1, ar2; length: integer); 
asm 
    push  ebx 
//Unfortunately delphi use first 8 bytes for dinamic array length and reference 
//counter, from that reason the dinamic array address should start with $xxxxxxx8 
//instead &xxxxxxx0. So we must first align ar1, ar2 pointers! 
    mov  ebx, [ar2] 
    sub  [ar1], ebx 
    dec  length 
    jz  @exit 
    mov  ebx, [ar2 + 4] 
    sub  [ar1 + 4], ebx 
    dec  length 
    jz  @exit 
    add  ar1, 8      //Align pointer to 16 byte 
    add  ar2, 8      //Align pointer to 16 byte 
//Prepare length for 16 byte data transfer 
    push  length 
    shr  length, 5     //8 * 4 subtructions per loop 
    jz  @Finish      //To small for LoopUnrolled 
@LoopUnrolled: 
//Unrolle 1, 2, 3, 4 
    movdqa xmm4, [ar2] 
    movdqa xmm5, [16 + ar2] 
    movdqa xmm6, [32 + ar2] 
    movdqa xmm7, [48 + ar2] 
// 
    movdqa xmm0, [ar1] 
    movdqa xmm1, [16 + ar1] 
    movdqa xmm2, [32 + ar1] 
    movdqa xmm3, [48 + ar1] 
// 
    psubd xmm0, xmm4 
    psubd xmm1, xmm5 
    psubd xmm2, xmm6 
    psubd xmm3, xmm7 
// 
    movdqa [48 + ar1], xmm3 
    movdqa [32 + ar1], xmm2 
    movdqa [16 + ar1], xmm1 
    movdqa [ar1], xmm0 
//Unrolle 5, 6, 7, 8 
    movdqa xmm4, [64 + ar2] 
    movdqa xmm5, [80 + ar2] 
    movdqa xmm6, [96 + ar2] 
    movdqa xmm7, [112 + ar2] 
// 
    movdqa xmm0, [64 + ar1] 
    movdqa xmm1, [80 + ar1] 
    movdqa xmm2, [96 + ar1] 
    movdqa xmm3, [112 + ar1] 
// 
    psubd xmm0, xmm4 
    psubd xmm1, xmm5 
    psubd xmm2, xmm6 
    psubd xmm3, xmm7 
// 
    movdqa [112 + ar1], xmm3 
    movdqa [96 + ar1], xmm2 
    movdqa [80 + ar1], xmm1 
    movdqa [64 + ar1], xmm0 
// 
    add  ar1, 128 
    add  ar2, 128 
    dec  length 
    jnz  @LoopUnrolled 
@FinishUnrolled: 
    pop  length 
    and  length, $1F 
//Do rest, up to 31 subtractions... 
@Finish: 
    mov  ebx, [ar2] 
    sub  [ar1], ebx 
    add  ar1, 4 
    add  ar2, 4 
    dec  length 
    jnz  @Finish 
@Exit: 
    pop  ebx 
end; 

procedure WriteOut(EndTicks: Int64; Str: string); 
begin 
    WriteLn(Str + IntToStr(EndTicks - StartTicks) 
    + ' Time: ' + IntToStr((EndTicks - StartTicks) div TicksPerMicroSecond) + 'us'); 
    Sleep(5); 
    SwitchToThread; 
    StartTicks := GetCpuTicks; 
end; 

begin 
    ArLength := 100000; 
//Set TicksPerMicroSecond 
    QueryPerformanceFrequency(TicksPerMicroSecond); 
    TicksPerMicroSecond := TicksPerMicroSecond div 1000000; 
// 
    SetLength(Ar1, ArLength); 
    SetLength(Ar2, ArLength); 
//Fill arrays 
//... 
//Tick time info 
    WriteLn('CPU ticks per mikro second: ' + IntToStr(TicksPerMicroSecond)); 
    Sleep(5); 
    SwitchToThread; 
    StartTicks := GetCpuTicks; 
//Test 1 
    SubArPasRangeOvfChkOn(ArLength); 
    WriteOut(GetCpuTicks, 'SubAr Pas Range and Overflow Checking On, Ticks: '); 
//Test 2 
    SubArPas(ArLength); 
    WriteOut(GetCpuTicks, 'SubAr Pas, Ticks: '); 
//Test 3 
    SubArAsm(Ar1[0], Ar2[0], ArLength); 
    WriteOut(GetCpuTicks, 'SubAr Asm, Ticks: '); 
//Test 4 
    SubArAsmSimdU(Ar1[0], Ar2[0], ArLength); 
    WriteOut(GetCpuTicks, 'SubAr Asm SIMD mem unaligned, Ticks: '); 
//Test 5 
    SubArAsmSimdA(Ar1[0], Ar2[0], ArLength); 
    WriteOut(GetCpuTicks, 'SubAr Asm with SIMD mem aligned, Ticks: '); 
//Test 6 
    SubArAsmSimdAUnrolled8(Ar1[0], Ar2[0], ArLength); 
    WriteOut(GetCpuTicks, 'SubAr Asm with SIMD mem aligned 8*unrolled, Ticks: '); 
// 
    ReadLn; 
    Ar1 := nil; 
    Ar2 := nil; 
end. 

...

Самая быстрая процедура ASM с 8 раз инструкции SIMD развернутой принимает только 68us и составляет около 4 раз быстрее, чем процедуры Паскаля.

Как мы видим, процедура цикла Pascal, вероятно, не является критичной, требуется только около 277us (проверка переполнения и диапазона) на частоте 2,4 ГГц при 100000 вычитаниях.

Таким образом, этот код не может быть узким местом?

Answer 2

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

Answer 3

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

Вы отключили настройки компилятора, проверку переполнения и диапазона?

Вы можете попробовать использовать ассемблерный rutine, это очень просто ...

Что-то вроде:

procedure SubArray(var ar1, ar2; length: integer); 
asm 
//length must be > than 0! 
    push ebx 
    lea ar1, ar1 -4 
    lea ar2, ar2 -4 
@Loop: 
    mov ebx, [ar2 + length *4] 
    sub [ar1 + length *4], ebx 
    dec length 
//Here you can put more folloving parts of rutine to more unrole it to speed up. 
    jz @exit 
    mov ebx, [ar2 + length *4] 
    sub [ar1 + length *4], ebx 
    dec length 
// 
    jnz @Loop 
@exit: 
    pop ebx 
end; 


begin 
    SubArray(Array1[0], Array2[0], length(Array1)); 

Это может быть намного быстрее ...

EDIT: Добавлена с инструкциями SIMD. Эта процедура запрашивает поддержку процессора SSE. Он может принимать 4 целых числа в регистре XMM и вычитать сразу. Существует также возможность использовать movdqa вместо movdqu быстрее, но сначала вы должны обеспечить 16-байтовый алимент. Вы также можете развернуть XMM-пар, как в моем первом asm-случае. (Я интересно о скорости Измерение. :))

var 
    array1, array2: array of integer; 

procedure SubQIntArray(var ar1, ar2; length: integer); 
asm 
//prepare length if not rounded to 4 
    push  ecx 
    shr  length, 2 
    jz  @LengthToSmall 
@Loop: 
    movdqu xmm1, [ar1]   //or movdqa but ensure 16b aligment first 
    movdqu xmm2, [ar2]   //or movdqa but ensure 16b aligment first 
    psubd xmm1, xmm2 
    movdqu [ar1], xmm1   //or movdqa but ensure 16b aligment first 
    add  ar1, 16 
    add  ar2, 16 
    dec  length 
    jnz  @Loop 
@LengthToSmall: 
    pop  ecx 
    push  ebx 
    and  ecx, 3 
    jz  @Exit 
    mov  ebx, [ar2] 
    sub  [ar1], ebx 
    dec  ecx 
    jz  @Exit 
    mov  ebx, [ar2 + 4] 
    sub  [ar1 + 4], ebx 
    dec  ecx 
    jz  @Exit 
    mov  ebx, [ar2 + 8] 
    sub  [ar1 + 8], ebx 
@Exit: 
    pop  ebx 
end; 

begin 
//Fill arrays first! 
    SubQIntArray(Array1[0], Array2[0], length(Array1)); 

Answer 4

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

как
Thread1: SubArray (ar1 [0], ar2 [0], 50);
Thread2: SubArray (ar1 [50], ar2 [50], 50);

procedure SubArray(var Array1, Array2; const Length: Integer); 
var 
    ap1, ap2 : PInteger; 
    i : Integer; 
begin 
    ap1 := @Array1; 
    ap2 := @Array2; 
    i := Length; 
    while i > 0 do 
    begin 
    ap1^ := ap1^ - ap2^; 
    Inc(ap1); 
    Inc(ap2); 
    Dec(i); 
    end; 
end; 

// similar assembly version 
procedure SubArrayEx(var Array1, Array2; const Length: Integer); 
asm 
    // eax = @Array1 
    // edx = @Array2 
    // ecx = Length 
    // esi = temp register for array2^ 
    push esi 
    cmp ecx, 0 
    jle @Exit 
    @Loop: 
    mov esi, [edx] 
    sub [eax], esi 
    add eax, 4 
    add edx, 4 
    dec ecx 
    jnz @Loop 
    @Exit: 
    pop esi 
end; 


procedure Test(); 
var 
    a1, a2 : array of Integer; 
    i : Integer; 
begin 
    SetLength(a1, 3); 
    a1[0] := 3; 
    a1[1] := 1; 
    a1[2] := 2; 
    SetLength(a2, 3); 
    a2[0] := 4; 
    a2[1] := 21; 
    a2[2] := 2; 
    SubArray(a1[0], a2[0], Length(a1)); 

    for i := 0 to Length(a1) - 1 do 
    Writeln(a1[i]); 

    Readln; 
end; 

Answer 5

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

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