Het artikel legt uit hoe trage integerdelingen ($x/y$) en modulo-operaties ($x\%y$) kunnen worden versneld door ze te verplaatsen naar floating-point hardware, die over het algemeen een lagere latentie en hogere doorvoersnelheid heeft.
De kernpunten zijn:
- Implementatie: Voor integers binnen de precisiegrenzen (53 bits voor double, 24 bits voor single) kan de integerdeling worden berekend met
d = trunc(x/y) en de rest via de Fused Multiply-Add operatie: m = -fma(d, y, -x).
- Wiskundige onderbouwing: De auteur bewijst dat bij de standaard afrondingsmodus (round-to-nearest) er geen fouten optreden. Omdat 'ties' onmogelijk zijn en het fractionele deel nooit groot genoeg is om een onterechte afronding naar boven te veroorzaken, blijft het integer-deel identiek.
- Praktijk: Er wordt gewezen op de invloed van conversiekosten (int naar float), de voordelen van SIMD-implementaties en specifieke hardware-ondersteuning zoals AVX-512.
Integerdeling omzetten naar floating-point is triviaal
Er zijn daarom gevallen waarin het interessant kan zijn om bepaalde integer div/mod-operaties te verplaatsen naar floating-point. Hoewel dit complex lijkt, is het in werkelijkheid vrij eenvoudig. De wiskunde is rechttoe rechtaan.
De claim
Voor twee integers $x$ en $y$ (signed of unsigned) die passen binnen 53 bits voor double precision of 24 bits voor single precision, waarbij beide zijn gepromoveerd naar floating-point, geldt het volgende in de standaard afrondingsmodus (round-to-nearest, ties to even):
// Floating-point berekening:
// d is dezelfde integer als de integerdeling x / y
// m is dezelfde integer als de integerrest x % y
d = trunc(x/y); // floor werkt voor unsigned
m = -fma(d, y, -x); // fma (Fused Multiply-Add) is vereist
Opmerking: Als alleen d nodig is en deze vervolgens wordt geconverteerd naar een integer, dan is de truncate- of floor-operatie gratis inbegrepen bij de conversie van float naar integer.
In het vervolg wordt er alleen gekeken naar unsigned integers, omdat dit het moeilijkste scenario is (signed integers hebben een kleinere maximale magnitude), maar houd er rekening mee dat floating-point effectief een signed magnitude quantity opslaat.
Praktische overwegingen
- Afrondingsmodi: Dit werkt door de afrondingsmodus in te stellen op
TOWARD_ZERO en deze daarna terug te zetten, maar het wijzigen van control words is vaak erg duur.
- Hardware-ondersteuning: Sommige hardware heeft divisie-opcodes waarmee de afrondingsmodus per operatie kan worden gekozen. AVX-512 heeft bijvoorbeeld de intrinsic
mmdivroundsd met een geliste latentie van 14 cycli.
- Conversiekosten: Het converteren van int naar float en terug kost tijd. (Let op: de meeste x64 CPU's hebben geen specifieke op voor unsigned naar float en terug, wat een implementatievraagstuk is.)
- Precisie: De floating-point precisie $p$ bepaalt de bovengrens voor de breedte van de integer. Als het werkt voor $p$, dan werken alle kleinere breedtes triviaal. Voorbeelden zijn 32-bit integers in doubles ($p=53$) en 16-bit integers in singles ($p=24$). In deze gevallen is het mogelijk dat er "trucs" zijn om de int-naar-float conversie te omzeilen.
- SIMD: De meest veelbelovende toepassing is in SIMD, om de overhead te amortiseren.
- Constante delers: Bij constante delers zou de compiler dit moeten optimaliseren. Voor delers die tijdens runtime bekend zijn en een klein aantal keren worden hergebruikt, bestaan software-oplossingen zoals libdivide.
Er bestaan overigens reeds bewezen methoden voor dit proces (zie: Formally verified 32- and 64-bit integer division using double-precision floating-point arithmetic).
Wiskundige onderbouwing
De standaard afrondingsmodus is round-to-nearest (ties to even). De eerste belangrijke observatie betreft het "tie"-gedeelte. Een tie treedt op wanneer het exacte resultaat van een operatie exact in het midden ligt tussen twee floating-point getallen. Floating-point deling (in ons geval met basis-2 en dezelfde werkprecisie) heeft echter geen midpunten. Daarom treedt er nooit een tie op en is alleen het round-to-nearest gedeelte van toepassing.
Laten we een voorbeeld nemen met een floating-point formaat van 4-bit precisie, waarbij we kijken naar alle configuraties waarin het exacte resultaat kleiner is dan één, maar zo dicht mogelijk bij het afronden naar boven ligt:
| GRS | Resultaat | Actie |
|---|
| `.1111\ | 0xx` | .1111 | Geen afronding |
| `.1111\ | 100` | - | Tie case is onmogelijk |
| `.1111\ | 101` | 1.000 | Afronden naar boven |
| `.1111\ | 11x` | 1.000 | Afronden naar boven |
Vóór de afronding berekent de hardware drie extra bits: guard bit (G), round bit (R) en sticky bit (S). Omdat de tie case onmogelijk is, hebben we alleen de G-bit nodig. In elk formaat waar het exacte resultaat $r$ bits heeft voor het fractionele deel, is voor afronding naar de volgende integer vereist dat het fractionele deel ten minste $r+1$ leidende enen heeft.
De claim is dat dit onmogelijk is met legale inputs.
Analyse van het fractionele deel
Om te bepalen wat er gebeurt met het fractionele deel, splitsen we het exacte resultaat van $x/y$ op in een integer-deel $n$ en een fractioneel deel:
$$\frac{x}{y} = n + \frac{a}{b}$$
Waarbij $b=y$, dus duidelijk is dat $\frac{a}{b} \in [0, 1)$.
Bij een binair floating-point formaat met precisie $p$ produceert de deling een integer van $d$-bits, waardoor er $r$-bits overblijven voor de rest:
- $d = \lfloor \log_2(n) + 1 \rfloor$
- $r = p - d$
Voor een fractioneel deel van $r$-bits kunnen we de onder- en bovengrens definiëren voor de waarde van $b$:
- $\text{lower}(r) = 2^r$
- $\text{upper}(r) = 2^{r+1} - 1$
We willen weten wat het maximale fractionele deel is dat mogelijk is voor $r$-bits. De waarde die met de ondergrens het dichtst bij 1 komt is:
$$\frac{2^r-1}{2^r} = \sum_{i=1}^r \frac{1}{2^i}$$
Deze som laat zien dat dit het kleinste getal is dat $r$ leidende enen produceert in de binaire expansie. Als we nu kijken naar de waarde van de bovengrens die het dichtst bij 1 ligt, en deze vergelijken met het kleinste getal met $r+1$ leidende enen (de ondergrens van $r+1$), zien we:
$$\frac{2^{r+1}-2}{2^{r+1}-1} < \frac{2^{r+1}-1}{2^{r+1}-0}$$
Deze ongelijkheid betekent dat de bovengrens beperkt is tot $r$ leidende enen. We kunnen specifiek zijn over de binaire expansie van de grenzen: de ondergrens bestaat uit $r$ enen en stopt daar; de exacte waarde van de bovengrens is puur cyclisch met een periode van $r+1$, bestaande uit $r$ enen gevolgd door één nul.
Overzicht voor de eerste vier waarden van $r$:
| $r$ | Ondergrens | Binaire expansie | Bovengrens | Binaire expansie |
|---|
| 1 | $\frac{1}{2}$ | .1 | $\frac{2}{3}$ | .(10) |
| 2 | $\frac{3}{4}$ | .11 | $\frac{6}{7}$ | .(110) |
| 3 | $\frac{7}{8}$ | .111 | $\frac{14}{15}$ | .(1110) |
| 4 | $\frac{15}{16}$ | .1111 | $\frac{30}{31}$ | .(11110) |
Het belangrijkste punt is: voor legale inputs kan $\frac{x}{y}$ geen fractioneel deel produceren dat groot genoeg is om afronding naar de volgende integer te veroorzaken. Hierdoor is het integer-deel van het floating-point resultaat identiek aan het integer-resultaat.
Verificatie met delers
Om dit samen te vatten, kijken we naar enkele specifieke delers:
- Delen door 1: Is correct.
- Delen door 2: Modificeert alleen de exponent; er vindt geen afronding plaats.
- Delen door 3: Hier is $r=1$. De grootste waarde voor $x \% 3$ is 2, dus de maximale fractie is $\frac{2}{3}$. Dit kan niet worden afgerond naar de volgende integer.
- Delen door $\{4,5,6,7\}$: Hier is $r=2$. De exacte deler maakt niet uit, omdat ze allemaal beginnen met
.110x.
Dit geldt voor alle mogelijke delers.
Integerdeling omzetten naar floating-point is triviaal
Er zijn daarom gevallen waarin het interessant kan zijn om bepaalde integer div/mod-operaties te verplaatsen naar floating-point. Hoewel dit complex lijkt, is het in werkelijkheid vrij eenvoudig. De wiskunde is rechttoe rechtaan.
De claim
Voor twee integers $x$ en $y$ (signed of unsigned) die passen binnen 53 bits voor double precision of 24 bits voor single precision, waarbij beide zijn gepromoveerd naar floating-point, geldt het volgende in de standaard afrondingsmodus (round-to-nearest, ties to even):
// Floating-point berekening:
// d is dezelfde integer als de integerdeling x / y
// m is dezelfde integer als de integerrest x % y
d = trunc(x/y); // floor werkt voor unsigned
m = -fma(d, y, -x); // fma (Fused Multiply-Add) is vereist
Opmerking: Als alleen d nodig is en deze vervolgens wordt geconverteerd naar een integer, dan is de truncate- of floor-operatie gratis inbegrepen bij de conversie van float naar integer.
In het vervolg wordt er alleen gekeken naar unsigned integers, omdat dit het moeilijkste scenario is (signed integers hebben een kleinere maximale magnitude), maar houd er rekening mee dat floating-point effectief een signed magnitude quantity opslaat.
Praktische overwegingen
- Afrondingsmodi: Dit werkt door de afrondingsmodus in te stellen op
TOWARD_ZERO en deze daarna terug te zetten, maar het wijzigen van control words is vaak erg duur.
- Hardware-ondersteuning: Sommige hardware heeft divisie-opcodes waarmee de afrondingsmodus per operatie kan worden gekozen. AVX-512 heeft bijvoorbeeld de intrinsic
mmdivroundsd met een geliste latentie van 14 cycli.
- Conversiekosten: Het converteren van int naar float en terug kost tijd. (Let op: de meeste x64 CPU's hebben geen specifieke op voor unsigned naar float en terug, wat een implementatievraagstuk is.)
- Precisie: De floating-point precisie $p$ bepaalt de bovengrens voor de breedte van de integer. Als het werkt voor $p$, dan werken alle kleinere breedtes triviaal. Voorbeelden zijn 32-bit integers in doubles ($p=53$) en 16-bit integers in singles ($p=24$). In deze gevallen is het mogelijk dat er "trucs" zijn om de int-naar-float conversie te omzeilen.
- SIMD: De meest veelbelovende toepassing is in SIMD, om de overhead te amortiseren.
- Constante delers: Bij constante delers zou de compiler dit moeten optimaliseren. Voor delers die tijdens runtime bekend zijn en een klein aantal keren worden hergebruikt, bestaan software-oplossingen zoals libdivide.
Er bestaan overigens reeds bewezen methoden voor dit proces (zie: Formally verified 32- and 64-bit integer division using double-precision floating-point arithmetic).
Wiskundige onderbouwing
De standaard afrondingsmodus is round-to-nearest (ties to even). De eerste belangrijke observatie betreft het "tie"-gedeelte. Een tie treedt op wanneer het exacte resultaat van een operatie exact in het midden ligt tussen twee floating-point getallen. Floating-point deling (in ons geval met basis-2 en dezelfde werkprecisie) heeft echter geen midpunten. Daarom treedt er nooit een tie op en is alleen het round-to-nearest gedeelte van toepassing.
Laten we een voorbeeld nemen met een floating-point formaat van 4-bit precisie, waarbij we kijken naar alle configuraties waarin het exacte resultaat kleiner is dan één, maar zo dicht mogelijk bij het afronden naar boven ligt:
| GRS | Resultaat | Actie |
|---|
| `.1111\ | 0xx` | .1111 | Geen afronding |
| `.1111\ | 100` | - | Tie case is onmogelijk |
| `.1111\ | 101` | 1.000 | Afronden naar boven |
| `.1111\ | 11x` | 1.000 | Afronden naar boven |
Vóór de afronding berekent de hardware drie extra bits: guard bit (G), round bit (R) en sticky bit (S). Omdat de tie case onmogelijk is, hebben we alleen de G-bit nodig. In elk formaat waar het exacte resultaat $r$ bits heeft voor het fractionele deel, is voor afronding naar de volgende integer vereist dat het fractionele deel ten minste $r+1$ leidende enen heeft.
De claim is dat dit onmogelijk is met legale inputs.
Analyse van het fractionele deel
Om te bepalen wat er gebeurt met het fractionele deel, splitsen we het exacte resultaat van $x/y$ op in een integer-deel $n$ en een fractioneel deel:
$$\frac{x}{y} = n + \frac{a}{b}$$
Waarbij $b=y$, dus duidelijk is dat $\frac{a}{b} \in [0, 1)$.
Bij een binair floating-point formaat met precisie $p$ produceert de deling een integer van $d$-bits, waardoor er $r$-bits overblijven voor de rest:
- $d = \lfloor \log_2(n) + 1 \rfloor$
- $r = p - d$
Voor een fractioneel deel van $r$-bits kunnen we de onder- en bovengrens definiëren voor de waarde van $b$:
- $\text{lower}(r) = 2^r$
- $\text{upper}(r) = 2^{r+1} - 1$
We willen weten wat het maximale fractionele deel is dat mogelijk is voor $r$-bits. De waarde die met de ondergrens het dichtst bij 1 komt is:
$$\frac{2^r-1}{2^r} = \sum_{i=1}^r \frac{1}{2^i}$$
Deze som laat zien dat dit het kleinste getal is dat $r$ leidende enen produceert in de binaire expansie. Als we nu kijken naar de waarde van de bovengrens die het dichtst bij 1 ligt, en deze vergelijken met het kleinste getal met $r+1$ leidende enen (de ondergrens van $r+1$), zien we:
$$\frac{2^{r+1}-2}{2^{r+1}-1} < \frac{2^{r+1}-1}{2^{r+1}-0}$$
Deze ongelijkheid betekent dat de bovengrens beperkt is tot $r$ leidende enen. We kunnen specifiek zijn over de binaire expansie van de grenzen: de ondergrens bestaat uit $r$ enen en stopt daar; de exacte waarde van de bovengrens is puur cyclisch met een periode van $r+1$, bestaande uit $r$ enen gevolgd door één nul.
Overzicht voor de eerste vier waarden van $r$:
| $r$ | Ondergrens | Binaire expansie | Bovengrens | Binaire expansie |
|---|
| 1 | $\frac{1}{2}$ | .1 | $\frac{2}{3}$ | .(10) |
| 2 | $\frac{3}{4}$ | .11 | $\frac{6}{7}$ | .(110) |
| 3 | $\frac{7}{8}$ | .111 | $\frac{14}{15}$ | .(1110) |
| 4 | $\frac{15}{16}$ | .1111 | $\frac{30}{31}$ | .(11110) |
Het belangrijkste punt is: voor legale inputs kan $\frac{x}{y}$ geen fractioneel deel produceren dat groot genoeg is om afronding naar de volgende integer te veroorzaken. Hierdoor is het integer-deel van het floating-point resultaat identiek aan het integer-resultaat.
Verificatie met delers
Om dit samen te vatten, kijken we naar enkele specifieke delers:
- Delen door 1: Is correct.
- Delen door 2: Modificeert alleen de exponent; er vindt geen afronding plaats.
- Delen door 3: Hier is $r=1$. De grootste waarde voor $x \% 3$ is 2, dus de maximale fractie is $\frac{2}{3}$. Dit kan niet worden afgerond naar de volgende integer.
- Delen door $\{4,5,6,7\}$: Hier is $r=2$. De exacte deler maakt niet uit, omdat ze allemaal beginnen met
.110x.
Dit geldt voor alle mogelijke delers.
Groetjes,