Branchless clamp to [0..255] [pl]

Parę dni temu na blogu Gynvaela Coldwinda pojawił się wpis o bardzo fajnej tematyce – zrobienia clampa dowolnej liczby całkowitej do zakresu [0..255] bez skoków. Postanowiłem się przyjrzeć temu z bliska, bo problem jest dość ciekawy i nie pozwolę, żeby mi ktoś chleb zabierał ;P.

Napisałem więc dwie funkcje, które także dają poprawny wynik i są branchless:

int my1(int x)
{
  x -= (x - 255) & -(x > 255);
  x += (-x) & -(x < 0);
  return x;
}

Pierwsza instrukcja spowoduje ustawienie x na 255 tylko wtedy, gdy x > 255. Jakim trafem to działa? x > 255 wynosi 1 (w rzeczywistości to wynosi true, ale bool ma niejawną konwersję do int), jeśli, jak łatwo się domyślić, x jest większe od 255 i 0 w przeciwnym wypadku. Jeśli zafundujemy temu negację arytmetyczną, to dostaniemy (x > 255) ? (-1) : (0). Ponieważ (prawie?) wszystkie architektury używają uzupełnienia do dwóch, -1 jest równe 0xFF...FF, czyli ma wszystkie bity równe 1. Pamiętając, że cokolwiek & 0xFF...FF = cokolwiek oraz cokolwiek & 0 = 0 uzyskujemy x -= (x > 255) ? (x - 255) : (0). Pozostaje nam to tajemnicze x - 255. Ile wynosi cokolwiek - (cokolwiek - 255)? Po opuszczeniu nawiasów otrzymujemy cokolwiek - cokolwiek + 255. cokolwiek - cokolwiek wynosi zawsze 0, więc zostaje 255. Ostatecznie daje nam to x = (x > 255) ? (255) : (x).

Druga linia jest podobna – jest to równoważne x += (x < 0) ? (-x) : (0). x + 0 to x, co jest oczywiste, więc pomyślmy nad tym -x. x + (-x) to x - x, czyli ładne 0. Całość sprowadza się więc do x = (x < 0) ? (0) : (x).

Druga funkcja natomiast jest nieco inna (wydaje się dłuższa, ale w rzeczywistości po skompilowaniu jest mniejsza):

int my2(int r)
{
  int min = r ^ ((255 ^ r) & -(r > 255));
  int max = min ^ ((0 ^ min) & -(min < 0));
  return max;
}

W tym przypadku całość jest chyba odrobinę (ale tylko odrobinę, bo ciągle jest mnóstwo operatorów bitowych, więc wskaźnik WTF/min jest wciąż wysoki ;>) czytelniejsza, ze względu na porządne nazwy zmiennych. Pierwsza instrukcja oblicza min(255, r). Druga natomiast max(0, min). Zwracany jest max, czyli efektywnie max(0, min(255, r)). Nie będę teraz tłumaczył czemu to tak wygląda, bo przypadkowo mam w planach pisanie o tym w części drugiej minitutka pisania branchless code /EDIT (16th November 2009): plany oczywiście się nie sprawdziły ;>. Teraz planuję o napisać o tym w trzeciej części, trudno, musisz zaczekać ;>/ ;).

Teraz słowo o wydajności. Przetestowałem cztery funkcje – te dwie pokazane wyżej, funkcję Gynvaela oraz pewną niespodziankę. Co się okazało? Wyniki były dla mnie niespodzianką, bo wyszło na to, że jestem lamką i nie potrafię napisać porządnego branchless code – moje funkcje uzyskiwały zwykle nieco gorsze czasy niż wynalazek Gyna. Najbardziej jednak zdziwiła mnie zwycięska funkcja, która w prawie wszystkich przypadkach była szybsza niż pozostałe… że co? Mam podać jej kod? Ok, proszę bardzo:

int clamp(int r)
{
  if(255 < r)
    return 255;
  else if(r < 0)
    return 0;
  else
    return r;
}

Tak, to jest ta funkcja. Najprostsza możliwa implementacja tego clampa okazała się najszybsza. Nawet już nie wspominam o czytelności. Ech, nie ma to jak niespodzianki sprawiane przez kompilator ;).

…co? Że niby to nie jest branchless? W takim razie sprawdźmy sobie jej kod wynikowy:

__Z5clampi:
; arg = esp + 4
mov eax, [esp + 4] ; eax = arg
xor edx, edx       ; edx = 0
test eax, eax      ;
cmovns edx, eax    ; if(eax >= 0) edx = eax
mov eax, 0FFh      ; eax = 0xFF
cmp edx, 0FFh      ;
cmovle eax, edx    ; if(edx <= 0xFF) eax = edx
ret                ; return eax

Nie ma skoków, mamy tylko conditional moves, więc jest branchless. Amen.

~ - autor: Fanael w dniu Listopad 10, 2009.