Почему этот код превосходит rint () и как мне защитить его от -ffast-math и друзей? - PullRequest
1 голос
/ 07 января 2020

Я ищу способ защитить часть кода от -ffast-math (или msvc / i cc эквивалентов и т. Д. c), который работает с C компиляторами.

My внутренняя l oop ищет данные для чисел, близких к целочисленным значениям (например, в пределах ~ 0,1). Значения данных подписаны, как правило, менее нескольких тысяч без инф / нан. Самая быстрая версия, которую я нашел, использует трюк с большим числом волхвов c:

 remainder = h - ((h+MAGIC)-MAGIC) ;

У кого-нибудь есть идеи, как сохранить порядок приоритетов скобок для ключевой строки выше? Похоже, что этот показатель превышает rint(x) в 3 раза, поэтому мне любопытно, почему он все равно работает. Может быть, это как-то связано с векторизацией?

Большинство компиляторов «упрощают» выражение при использовании -ffast-math или его эквивалента, и оно перестает работать. Я хочу сохранить производительность (3X довольно много), но также сохранить ее неопределенно переносимой (учитывая, что MAGI C зависит от наличия правильного ieee). Если я добавлю летучий, то он замедляется, но, похоже, дает правильные ответы с быстрой математикой, но затем медленнее, чем rint:

 volatile t = h+MAGIC; t-=MAGIC;
 remainder = h - t;

Полный пример приведен ниже. Я попробовал некоторые g cc такие вещи, как __attribute__((optimize("-fno-associative-math"))), но это не похоже на правильный подход, чтобы в конечном итоге работать для icc / gcc / msvc / clang et c. Соответствующая стандартная прагма C99 также, кажется, не является широко доступной.

#include <stdlib.h>
#include <stdio.h>
#include <math.h>
#include <sys/time.h>

/* https://stackoverflow.com/questions/17035464/a-fast-method-to-round-a-double-to-a-32-bit-int-explained */
union i_cast {double d; int i[2];};
#define MAGIC 6755399441055744.0
/* x86_64 for me today */
#define ENDIANLOC 0
#define REMAINDER_LUA  \
    {volatile union i_cast u; u.d = (h) + MAGIC; r = h - u.i[ENDIANLOC]; }

#define REMAINDER_MAGIC  r=(h - ((h+MAGIC)-MAGIC));
#define REMAINDER_RINT   r=(h - rint(h));
#define REMAINDER_TRUNC  r=(h - ( (h>0) ? ((int)(h+0.5)) : ((int)(h-0.5))) );
#define REMAINDER_FLOOR  r=(h - floor(h+0.5));
#define REMAINDER_REMAIN r=(remainder(h, 1.0));
#define REMAINDER_ROUND  r=(h - round(h));
#define REMAINDER_NEARBY r=(h - nearbyint(h));

#define block(MACRO) {                                                  \
    for(i=0 ; i<3 ; i++){                                               \
      gettimeofday(&start, NULL);                                       \
      n = 0;                                                            \
      for (k = 0; k < ng; k++) {                                        \
        h = mul * gv[k];                                                \
        MACRO                                                           \
          if ( (r*r) < tol ) n++;                                       \
      }                                                                 \
      gettimeofday(&end, NULL);                                         \
      double dt = (double)(end.tv_sec - start.tv_sec);                  \
      dt += (1e-6)*(double)(end.tv_usec - start.tv_usec);               \
      if(i==2)                                                          \
      printf("%20s %d indexed in %lf s %f ns/value\n",#MACRO,           \
             n,dt,1e9*dt/ng);                                           \
    }                                                                   \
  } 

int main(){
  struct timeval start, end;
  // Make some test data
  double h, r, tol = 0.02, mul = 123.4;
  int i, n, k, ng = 1024*1024*32;
  srand(42);
  double *gv = (double *) malloc(ng*sizeof(double));
  for(int i=0;i<ng;i++) { gv[i] = ((double)rand())/RAND_MAX*2.-1.; }
  // Measure some timing
  block(REMAINDER_MAGIC);
  block(REMAINDER_LUA);
  block(REMAINDER_RINT);
  block(REMAINDER_FLOOR);
  block(REMAINDER_TRUNC);
  block(REMAINDER_ROUND);
  block(REMAINDER_REMAIN);
  block(REMAINDER_NEARBY);
  free( gv );
  return 0;
}

Для меня сегодня вывод был таким: g cc -O3:

     REMAINDER_MAGIC 9489537 indexed in 0.017953 s 0.535041 ns/value
       REMAINDER_LUA 9489537 indexed in 0.048870 s 1.456439 ns/value
      REMAINDER_RINT 9489537 indexed in 0.050894 s 1.516759 ns/value
     REMAINDER_FLOOR 9489537 indexed in 0.086768 s 2.585888 ns/value
     REMAINDER_TRUNC 9489537 indexed in 0.162564 s 4.844785 ns/value
     REMAINDER_ROUND 9489537 indexed in 0.417856 s 12.453079 ns/value
    REMAINDER_REMAIN 9489537 indexed in 0.517612 s 15.426040 ns/value
    REMAINDER_NEARBY 9489537 indexed in 0.786896 s 23.451328 ns/value

Возможно, какой-нибудь другой язык (rust / go / opencl / что угодно) будет сделать лучше, чем C здесь? Или просто лучше контролировать флаги компилятора и добавить в код тест на время выполнения для корректности?

1 Ответ

2 голосов
/ 07 января 2020

Нет стандартного способа контролировать это нестандартное поведение. Каждый компилятор с параметром -ffast-math-style имеет атрибуты и прагмы для управления этим, но они различаются, как и точные эффекты этого параметра. В этом отношении разные версии некоторых из этих компиляторов имеют замечательно различное поведение fast-math, поэтому речь идет не просто о наборе соответствующих прагм. Стандартный способ получить стандартное поведение - заставить компилятор следовать стандарту языка.

-ffast-math и тому подобное предназначены в первую очередь для программистов, которые не заботятся о деталях математики с плавающей запятой, и просто хотите, чтобы их программы (которые ограниченно и консервативно использовали FP) работали быстрее. В любом случае большинство полезных эффектов -ffast-math можно продублировать с помощью тщательно написанного кода.

...