Подтвердить что ты не робот

Почему pow (int, int) так медленно?

Я работаю над несколькими упражнениями Эйлера для улучшения моих знаний о С++.

Я написал следующую функцию:

int a = 0,b = 0,c = 0;

for (a = 1; a <= SUMTOTAL; a++)
{
    for (b = a+1; b <= SUMTOTAL-a; b++)
    {
        c = SUMTOTAL-(a+b);

        if (c == sqrt(pow(a,2)+pow(b,2)) && b < c)
        {
            std::cout << "a: " << a << " b: " << b << " c: "<< c << std::endl;
            std::cout << a * b * c << std::endl;
        }
    }
}

Это вычисляется за 17 миллисекунд.

Однако, если я изменяю строку

if (c == sqrt(pow(a,2)+pow(b,2)) && b < c)

to

if (c == sqrt((a*a)+(b*b)) && b < c)

вычисление происходит через 2 миллисекунды. Есть ли некоторые очевидные детали реализации pow(int, int), которые я пропускаю, что делает первое выражение вычисляемым так медленнее?

4b9b3361

Ответ 1

pow() работает с реальными числами с плавающей запятой и использует под капотом формулу

pow(x,y) = e^(y log(x))

для вычисления x^y. int преобразуются в double перед вызовом pow. (log - натуральный логарифм, основанный на e)

x^2, используя pow(), поэтому медленнее, чем x*x.

Изменить на основе соответствующих комментариев

  • Использование pow даже с целыми показателями может привести к неправильным результатам (PaulMcKenzie)
  • В дополнение к использованию математической функции с двойным типом, pow является вызовом функции (в то время как x*x не является) (jtbandes)
  • Многие современные компиляторы фактически оптимизируют pow с постоянными целыми аргументами, но на это не следует полагаться.

Ответ 2

Вы выбрали один из самых медленных способов проверки

c*c == a*a + b*b   // assuming c is non-negative

Это компилируется с тремя целыми умножениями (один из которых может быть выведен из цикла). Даже без pow() вы все еще конвертируете в double и получаете квадратный корень, что ужасно для пропускной способности. (А также латентность, но предсказание ветвей + спекулятивное выполнение на современных CPU означает, что латентность здесь не является фактором).

Инструкция Intel Haswell SQRTSD имеет пропускную способность одного на 8-14 циклов (источник: таблицы инструкций Agner Fog), поэтому даже если ваш sqrt() поддерживает насыщенность исполняемого блока FP sqrt, он все еще примерно в 4 раза медленнее, чем то, что я получил gcc для испускания (ниже).


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

void foo_optimized()
{ 
  for (int a = 1; a <= SUMTOTAL; a++) {
    for (int b = a+1; b < SUMTOTAL-a-b; b++) {
        // int c = SUMTOTAL-(a+b);   // gcc won't always transform signed-integer math, so this prevents hoisting (SUMTOTAL-a) :(
        int c = (SUMTOTAL-a) - b;
        // if (b >= c) break;  // just changed the loop condition instead

        // the compiler can hoist a*a out of the loop for us
        if (/* b < c && */ c*c == a*a + b*b) {
            // Just print a newline.  std::endl also flushes, which bloats the asm
            std::cout << "a: " << a << " b: " << b << " c: "<< c << '\n';
            std::cout << a * b * c << '\n';
        }
    }
  }
}

Скомпилирует (с gcc6.2 -O3 -mtune=haswell) код с этим внутренним циклом. Полный код проводник компилятора Godbolt.

# a*a is hoisted out of the loop.  It in r15d
.L6:
    add     ebp, 1    # b++
    sub     ebx, 1    # c--
    add     r12d, r14d        # ivtmp.36, ivtmp.43  # not sure what this is or why it in the loop, would have to look again at the asm outside
    cmp     ebp, ebx  # b, _39
    jg      .L13    ## This is the loop-exit branch, not-taken until the end
                    ## .L13 is the rest of the outer loop.
                    ##  It sets up for the next entry to this inner loop.
.L8:
    mov     eax, ebp        # multiply a copy of the counters
    mov     edx, ebx
    imul    eax, ebp        # b*b
    imul    edx, ebx        # c*c
    add     eax, r15d       # a*a + b*b
    cmp     edx, eax  # tmp137, tmp139
    jne     .L6
 ## Fall-through into the cout print code when we find a match
 ## extremely rare, so should predict near-perfectly

В Intel Haswell все эти инструкции по 1 мкп каждый. (И макро-предохранитель cmp/jcc попадает в сопоставления и ветвления.) Таким образом, 10 fused-domain uops, которые могут выдаваться на одной итерации за 2,5 цикла.

Haswell запускает imul r32, r32 с пропускной способностью одной итерации за такт, поэтому два умножения внутри внутреннего цикла не насыщают порт 1 с двумя умножениями на 2.5c. Это оставляет место, чтобы впитать неизбежные конфликты ресурсов из ADD и SUB, крадущего порт 1.

Мы даже не близки к каким-либо другим узким местам исполнения, поэтому узкое место переднего плана является единственной проблемой, и это должно выполняться на одной итерации за 2,5 цикла на Intel Haswell, а затем.

Loop-unrolling может помочь здесь уменьшить количество попыток на проверку. например используйте lea ecx, [rbx+1] для вычисления b + 1 для следующей итерации, поэтому мы можем imul ebx, ebx не использовать MOV, чтобы сделать его неразрушающим.


Возможно также снижение прочности. Учитывая b*b, мы могли бы попытаться вычислить (b-1) * (b-1) без IMUL. (b-1) * (b-1) = b*b - 2*b + 1, поэтому, возможно, мы можем сделать lea ecx, [rbx*2 - 1], а затем вычесть из b*b. (Нет меток адресации, которые вычитают вместо добавления. Хм, возможно, мы могли бы сохранить -b в регистре и подсчитать до нуля, поэтому мы могли бы использовать lea ecx, [rcx + rbx*2 - 1] для обновления b*b в ECX, учитывая -b в EBX).

Если вы на самом деле не столкнулись с пропускной способностью IMUL, это может привести к большему количеству ошибок, а не победе. Было бы забавно видеть, насколько хорошо компилятор будет делать это с уменьшением силы в источнике С++.


Возможно, вы могли бы также векторизовать это с помощью SSE или AVX, одновременно проверяя 4 или 8 последовательных значений b. Поскольку хиты действительно редки, вы просто проверяете, был ли у кого-то из 8 ударов, а затем разобрал, какой из них был в редком случае, когда был матч.

См. также tag wiki для большей информации по оптимизации.