Войти
ПрограммированиеФорумОбщее

Трюки с float: быстрое вычисление логарифмов. (комментарии) (9 стр)

Страницы: 16 7 8 9 10 11 Следующая »
#120
10:06, 17 янв 2022

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

#121
(Правка: 10:14) 10:13, 17 янв 2022

Но выгода может быть получена уже за счёт того, что не надо вычислять величину (x-1)/(x+1). Если вычисление двух многочленов хорошо распараллелится (не важно, с помощью векторизации или суперскалярного исполнения), то можем получить код, не уступающий по скорости, но превосходящий по точности. Так что в ближайшее время займусь этой идеей.

#122
2:36, 18 янв 2022

FordPerfect
> Если вместо P((x-1)/(x+1)) делать P(x)/Q(x) (в #68 похожее), т. е. деление в конце - может быть лучше latency.
Попытался адаптировать свою процедуру \(L_\infty\) подгона под отношения, чет сходится не туда :/

+ есть у кого мысли?
#123
5:39, 18 янв 2022

}:+()___ [Smile]
У тебя константные x-координаты, а подбираешь ты коэффициенты P, Q?

#124
6:17, 18 янв 2022

FordPerfect
> У тебя константные x-координаты, а подбираешь ты коэффициенты P, Q?
Да, причем P и Q — это линейные комбинации произвольных базисных функций, не обязательно полиномов.

#125
7:21, 18 янв 2022

Мой тупорылый подгонятор
https://rextester.com/TURZP50644
даёт такое
https://rextester.com/ABCAC19569
Можно улучшать.

#126
(Правка: 19 янв 2022, 12:23) 13:31, 18 янв 2022

Выкладываю код функции с коэффициентами, полученными по методу Паде-Чебышёва (подобраны достаточно точно).

_declspec(naked) float _vectorcall ln1(float x)
{
  static const float ct[10] =           // Таблица констант
  {
    0.00830399152f,                     // a3
    0.0903585330f,                      // b3
    1.35235262f,                        // a1
    1.31841099f,                        // b1
    0.377887428f,                       // a2
    0.667115986f,                       // b2
    1.09942472f,                        // a0
    0.762063146f,                       // b0
    1.0f,                               // 1
    0.693147181f                        // ln 2
  };
  _asm
  {
    mov edx,-126                        // В edx - смещение для коррекции порядка нормализованных float
    vmovd eax,xmm0                      // Прочитать в eax двоичное представление числа x
    cmp eax,0x7F800000                  // Сравнить число x со значением Inf
    jge ln_exit                         // Если x=+Inf или +NaN, вернуть аргумент x без изменений
    ror eax,23                          // Циклически сдвинуть число так, чтобы порядок оказался в al
    movzx ecx, al                       // Получить в ecx порядок числа со смещением +127
    jecxz ln_denorm                     // Перейти, если x=0 или x денормализованное (т.е. близко к 0)
      ln_calc:                          // Точка входа для вычисления логарифма после нормализации числа 
    cmp eax,0x570A3E00                  // Сравнить мантиссу с 1.34
    setb al                             // al=1, если мантисса меньше 1.34; иначе al=0
    sub cl,al                           // ecx=k+126, где k - добавка к двоичному логарифму (целая часть)
    or al,126                           // Изменить порядок x так, чтобы привести к диапазону 0.67<=x<1.34
    add ecx,edx                         // ecx=k
    ror eax,9                           // Получить в eax двоичное представление числа и проверить его знак
    jnc ln_cont                         // Перейти, если x>0, т.е. аргумент корректный. Продолжить, если x<0
      ln_nan:                           // Точка входа для формирования результата NaN
    sbb eax,eax                         // eax=-1 (значение NaN) в случае x<0
      ln_inf:                           // Точка входа для формирования результата -Inf в случае x=0
    or eax,0xFF800000                   // Сформировать в eax значение -Inf, если x=0; или NaN, если x<0
    vmovd xmm0,eax                      // Перенести результат из eax в xmm0
      ln_exit:                          // В xmm0 готов результат для нештатного значения x
    ret                                 // Вернуться
      ln_denorm:                        // Точка входа для обработки денормализованного x, в т.ч. x=0
    btr eax,8                           // Прочитать в cf и сбросить знаковый бит числа
    bsr ecx,eax                         // Найти номер старшего единичного бита; установить zf, если x=0
    jz ln_inf                           // Если x=0, перейти для формирования результата -Inf
    jc ln_nan                           // Если x<0 (или -NaN), перейти для формирования результата NaN
    lea edx,[edx+2*ecx-31]              // Обновить смещение порядка для денормализованного числа
    neg ecx                             // 32-ecx - это необходимое количество сдвигов мантиссы влево
    shl eax,cl                          // Сформировать мантиссу нормализованного числа
    jmp ln_calc                         // Число нормализовано. Перейти для вычисления логарифма
      ln_cont:                          // Продолжаем вычисления в штатном случае, когда x>0
    vmovd xmm0, eax                     // Перенести из eax в xmm0 приведённое значение x
    mov eax, offset ct                  // В eax адрес таблицы констант
    vsubss xmm0, xmm0, [eax+32]         // xmm0 = x-1 = t
    vmovups xmm1,[eax]                  // xmm1 = b1 : a1 | b3 : a3
    vshufps xmm0,xmm0,xmm0,0            // xmm0 = t : t | t : t
    vcvtsi2ss xmm4,xmm4,ecx             // xmm4 = k
    vmulps xmm2,xmm0,xmm0               // xmm2 = t^2 : t^2 | t^2 : t^2
    vfmadd213ps xmm1, xmm0, [eax+16]    // xmm1 = b1*t+b0 : a1*t+a0 | b3*t+b2 : a3*t+a2
    vmovhlps xmm3,xmm3,xmm1             // xmm3 = b1*t+b0 : a1*t+a0
    vfmadd213ps xmm1, xmm2, xmm3        // xmm1 = b3*t^3+b2*t^2+b1*t+b0 : a3*t^3+a2*t^2+a1*t+a0 = Q : P
    vmovshdup xmm2,xmm1                 // xmm2 = Q
    vdivss xmm1,xmm1,xmm2               // xmm1 = P/Q
    vfmadd213ss xmm0, xmm1, xmm4        // xmm1 = P/Q*t+k - готовый двоичный логарифм
    vmulss xmm0, xmm0, [eax+36]         // Перевести двоичный логарифм в натуральный
    ret                                 // Вернуться
  }
}
#127
13:36, 18 янв 2022

Подскажите быстрые аппроксимация для pow(x,2.2) и pow(x,1/2.2).
Для 2.2 я нашел такое: https://www.shadertoy.com/view/WlG3zG
Для 1/2.2 лучшее до чего догадался это использовать exp(0.454545 * log(x)), но это, как я понимаю точное решение, и скорее всего стандартный powf как-то также реализован.
Можно ли лучше?

#128
19:25, 18 янв 2022

HolyDel
Тебе гамму?
Моя тема: https://gamedev.ru/flame/forum/?id=226959&page=2&m=5463220#m23 .

#129
19:33, 18 янв 2022

FordPerfect
Ага, спасибо! Интересная тема, читаю.

#130
5:39, 19 янв 2022

Чуть допилил. Точность всё-равно похуже, чем P((x-1)/(x+1)). Скорость - неочевидно.
https://rextester.com/KJAP13381

float log2_v7(float x)
{
    float y=-127.0f;
    uint32_t bits,bias=0x00400000u;
    memcpy(&bits,&x,sizeof(float));
    if((bits>>23)-1u>253u)
    {
        if(!(x>0.0f)) return (x==0.0f?-INFINITY:nanf(""));
        if(x==+INFINITY) return +INFINITY;
        x*=4294967296.0f;
        memcpy(&bits,&x,sizeof(float));
        y-=32.0f;
    }
    bits+=bias;
    y+=float(int32_t(bits>>23));
    bits=((bits&0x007FFFFFu)|0x3F800000u)-bias;
    memcpy(&x,&bits,sizeof(float));
    float z=x-1.0f;
    float P=z*(x+z*(-2.25607605e-02f+z*1.71451341e-01f));
    float Q=6.93147180e-01f
        +z*(1.02408259e+00f
        +z*(3.99834445e-01f
        +z*(3.18533542e-02)));
    return y+P/Q;
}
#131
(Правка: 12:22) 6:21, 19 янв 2022

FordPerfect, производительночть у меня получилась почти идентичная на AMD FX, и немного хуже — на Ryzen 5600X. Точность — тоже немного похуже. Причём основной источник ошибки — это шум из-за погрешностей округления: сначала округляются P и Q, потом их частное. Я коэффициенты подобрал максимально точно, ошибка аппроксимации пренебрежимо мала, тем не менее ошибка вычислений получается больше, чем в сумме многочлена по схеме Горнера.

#132
(Правка: 9:37) 9:37, 19 янв 2022

HolyDel
> Подскажите быстрые аппроксимация для pow(x,2.2) и pow(x,1/2.2).
Если тебе гамму - там сама степень 2.2 является только приближением, в реальной формуле диапазон делится на две части. Если мы хотим более быстрое приближение - то лучше подогнать его к оригинальной формуле, а не к pow(x, 2.2).

#133
6:23, 20 янв 2022

В #130 баг - одна из констант в double.
Поправил (и чуть подвинул):
https://rextester.com/ABXI98604
Точность всё-равно не ахти.

#134
6:39, 20 янв 2022

FordPerfect, ну исправьте прямо в сообщении. Я не совсем понял, по какой формуле вы аппроксимируете. У вас вижу почему-то 6 констант. У меня 8 констант, но одна — свободная, т.е. фактически 7. По моим расчётам, в данном случае это минимальное количество коэффициентов для качественной аппроксимации дробно-рациональной функцией с одинарной точностью.
И, кстати, точность у меня получилось лишь незначительно хуже (точные цифры пока не готов привести), чем через P((x-1)/(x+1)). У вас на сколько хуже точность?

Страницы: 16 7 8 9 10 11 Следующая »
ПрограммированиеФорумОбщее