Podpora numerických hodnot typu quadruple-precision floating-point v GCC

Obsah

1. Podpora hodnot typu quadruple-precision floating-point v GCC

2. Rozšířený formát extended/temporary

3. Varianta quadruple – hodnoty se čtyřnásobnou přesností

4. Varianta octuple – 256 bitů pro jedno číslo

5. Podpora hodnot typu quadruple v překladači a knihovnách GCC

6. Prototypy všech funkcí definovaných v quadmath.h

7. Tisk hodnot typu __float128

8. Konstanty definované v hlavičkových souborech math.h a quadmath.h

9. Aritmetické operace a speciální hodnoty typu __float128

10. Volání funkce definované v hlavičkovém souboru quadmath.h

11. Rozdílné výsledky operací s hodnotami s plovoucí řádovou čárkou

12. Iterativní výpočet odmocniny realizovaný s typy float, double i quadruple

13. Převod řetězce na hodnotu typu __float128

14. Podpora pro práci s komplexními čísly se 128bitovými reálnými a imaginárními složkami

15. Tisk numerických hodnot typu __complex128

16. Základní aritmetické operace s hodnotami typu __complex128

17. Funkce akceptující parametr typu __complex128 a vracející __float128

18. Funkce akceptující parametr typu __complex128 a vracející __complex128

19. Repositář s demonstračními příklady

20. Odkazy na Internetu

1. Podpora hodnot typu quadruple-precision floating-point v GCC

Naprostá většina mainstreamových programovacích jazyků podporuje práci s numerickými hodnotami s plovoucí řádovou čárkou, které jsou reprezentovány buď ve formátu single/float nebo double. Jedná se o názvy, které do jisté míry vychází z normy IEEE 754, přesněji řečeno z její varianty IEEE 754–1985. Termín single vychází ze sousloví single precision neboli „jednoduchá přesnost“. A double vychází ze sousloví double precision, tj. „dvojitá přesnost“.

Později začaly být tyto formáty označovány i jako binary32 a binary64, a to z toho důvodu, aby došlo k odlišení od decimálních formátů decimal32, decimal64 a později decimal128. A aby to nebylo příliš jednoduché – programovací jazyky mohou používat i jiná označení, například f32, f64, float32, float64 apod.

Formáty single a double definované v normě IEEE 754 patří mezi základní formáty (basic). Do této skupiny formátů dnes řadíme i quadruple-precision floating-point, tj. hodnoty se „čtyřnásobnou přesností“ A právě tímto formátem se budeme zabývat v dnešním článku, protože se jedná o formát, který je podporován překladačem GCC.

2. Rozšířený formát extended/temporary

Kromě obou základních formátů (tj. jednoduché i dvojité přesnosti, ke kterým byl přidán i dnes popisovaný formát se čtyřnásobnou přesností) je v normě IEEE 754 povoleno používat i rozšířené binární formáty (a formáty decimální, ty však v dnešním článku vynecháme). Na platformě x86 je při výpočtech prováděných matematickým koprocesorem používán rozšířený formát nazývaný extended či temporary. Tento formát je zajímavý tím, že pro uložení hodnot s plovoucí řádovou čárkou používá 80 bitů a je do něho možné beze ztráty přesnosti uložit 64bitové hodnoty typu integer (přesně taková je totiž bitová šířka mantisy). Osmdesátibitový vektor je rozdělený do třech skupin (bitových polí) následujícím způsobem:

• 1 bit pro znaménko
• 15 bitů pro exponent (BIAS je roven 16383)
• 64 bitů pro mantisu (maximální hodnota přesahuje 10⁴⁹³²)

U tohoto formátu je poměrně zajímavá funkce bitu s indexem 63. Podle hodnoty tohoto bitu se rozlišují čísla normalizovaná a nenormalizovaná (tento bit ve skutečnosti nahrazuje implicitně nastavovaný nejvyšší bit mantisy, jak ho známe z předchozích formátů, tedy single a double). Matematický koprocesor Intel 8087 a jeho následovníci sice dokáže pracovat s čísly nenormalizovanými, výsledkem jeho aritmetických operací jsou však vždy hodnoty normalizované. Všechny možnosti, které mohou při ukládání extended FP formátu nastat, jsou přehledně vypsány v následující tabulce:

Pro normalizované i nenormalizované hodnoty je možné uloženou hodnotu vyjádřit pomocí vzorce (všimněte si, že bit 63 je umístěn před binární tečkou):

X_extended=(-1)^s × 2^exp-16383 × m

m=m₆₃⁰+m₆₂^-1+m₆₁^-2+…+m₀^-63

#include <stdio.h>
 
extern long double addLong(long double x, long double y) {
    return x+y;
}
 
extern long double subLong(long double x, long double y) {
    return x-y;
}
 
extern long double mulLong(long double x, long double y) {
    return x*y;
}
 
extern long double divLong(long double x, long double y) {
    return x/y;
}
 
int main(void) {
    printf("%ld\n", sizeof(long double));
    return 0;
}

Po překladu do assembleru snadno zjistíme, že se zde používají původní instrukce matematického koprocesoru:

addLong:
        fld     TBYTE PTR [rsp+8]
        fld     TBYTE PTR [rsp+24]
        faddp   st(1), st
        ret
 
subLong:
        fld     TBYTE PTR [rsp+8]
        fld     TBYTE PTR [rsp+24]
        fsubp   st(1), st
        ret
 
mulLong:
        fld     TBYTE PTR [rsp+8]
        fld     TBYTE PTR [rsp+24]
        fmulp   st(1), st
        ret
 
divLong:
        fld     TBYTE PTR [rsp+8]
        fld     TBYTE PTR [rsp+24]
        fdivp   st(1), st
        ret

3. Varianta quadruple – hodnoty se čtyřnásobnou přesností

„For now the 10-byte Extended format is a tolerable compromise between the value of extra-precise arithmetic and the price of implementing it to run fast; very soon two more bytes of precision will become tolerable, and ultimately a 16-byte format… That kind of gradual evolution towards wider precision was already in view when IEEE Standard 754 for Floating-Point Arithmetic was framed.“
William Kahan, architekt IEEE 754

Rozšiřování bitové šířky mantisy i možných rozsahů exponentu (a tedy i počtu bitů exponentu) se nezastavilo u výše popsaného formátu extended. Ostatně formát extended byl již v době svého vzniku chápán spíše jako pouze dočasné řešení, protože tehdejší technologie neumožňovala, aby matematický koprocesor Intel 8087 obsahoval příliš mnoho tranzistorů (jednalo by se o drahé řešení, navíc zvýšení počtu tranzistorů vedlo k menší výtěžnosti výroby, což opět vede ke zdražení). V novější normě IEEE 754–2008 je kromě single, double a volitelného formátu extended specifikován další nepovinný formát nazvaný binary128, který se ovšem běžně označuje quadruple precision či jen quad precision. Tento formát je založen na slovech širokých celých 128 bitů (16 bajtů), která jsou rozdělena následujícím způsobem:

1. 1 bit pro znaménko
2. 15 bitů pro exponent
3. 112 bitů pro mantisu

Bitově vypadá rozdělení 128 bitů do jednotlivých skupin (bitových polí) následovně:

Exponent je v tomto případě posunutý o hodnotu bias=16383. Dekadická přesnost u tohoto formátu dosahuje 33 až 36 platných číslic! To je skutečně více než dvakrát tolik, než u formátu double (15–17 platných číslic) a více než čtyřikrát tolik v porovnání s formátem single (6–9 platných číslic). Název nového formátu je tedy dosti přesný.

Podpora formátu quadruple v RISC-V

Formát quadruple je podporován těmi mikroprocesory s architekturou RISC-V, které implementují rozšíření instrukční sady „Q“. Toto rozšíření je platné pouze pro 64bitové varianty procesorů, které obsahují sadu RV64I (u 32bitových variant CPU nebude „Q“ podporována). V rámci podpory rozšíření instrukční sady „Q“ se zvětšila šířka FP registrů, a to na 128 bitů. Většina původních instrukcí pro zpracování hodnot s plovoucí řádovou čárkou zůstala zachována, pouze se změnil jejich postfix určující typ zpracovávaných numerických hodnot. Jediné instrukce, které byly odstraněny, jsou instrukce pro bitovou kopii dat mezi celočíselným registrem a FP registrem (to je vlastně pochopitelné), což znamená nutnost načítání a ukládání operandů přímo do operační paměti. V instrukční sadě „Q“ tedy neexistují tyto dvě instrukce: FMV.X.Q a FMV.Q.X.

Pro zajímavost si vypišme nové instrukce pro načítání a ukládání operandů do operační paměti:

Základní aritmetické operace jsou (opět) prakticky stejné, pouze mají odlišný postfix:

Konverze mezi typy single, float a quadruple:

Operace se znaménkem:

4. Varianta octuple – 256 bitů pro jedno číslo

Jen v krátkosti se pro úplnost zmiňme o poslední variantě standardního binárního formátu uložení numerických hodnot s plovoucí řádovou čárkou, který se nazývá binary256 či méně formálně octuple precision („osminásobná přesnost“). Tento formát využívá slova o šířce plných 256 bitů (32 bajtů) s následujícím rozdělením bitů do skupin (bitových polí):

1. 1 bit pro znaménko
2. 19 bitů pro exponent
3. 236 bitů pro mantisu

(povšimněte si, že šířka exponentu nemusí růst nijak radikálně, ostatně již tak je rozsah hodnot více než dostatečný pro všechny myslitelné typy výpočtů)

Bitově vypadá rozdělení následovně:

Exponent je v tomto případě posunutý o hodnotu bias=262143. Dekadická přesnost u tohoto formátu dosahuje 71 platných číslic. Nejmenší (nenormalizovaná) reprezentovatelná hodnota rozdílná od nuly je přibližně 10⁻⁷⁸⁹⁸⁴, maximální hodnota pak 1,611 ×10⁷⁸⁹¹³ (těžko říct, zda je takový rozsah vůbec reálně využitelný).

5. Podpora hodnot typu quadruple v překladači a knihovnách GCC

Překladač GCC do značné míry podporuje práci s hodnotami typu quadruple. Tento formát je reprezentován typem __float128. Podpora tohoto formátu ovšem není univerzální. Plnou podporu nalezneme na platformách i386, x86–64, PowerPC a IA-64. Ovšem i přes relativně dobrou podporu je nutné si dát pozor například při tisku hodnot atd. Podrobnosti budou uvedeny v praktické části dnešního článku.

Většina hlaviček funkcí (ale i definicí konstant), které nějakým způsobem souvisejí s datovými typy __float128 a __complex128, je uložena v hlavičkovém souboru nazvaném quadmath.h. Zaměřme se nejdříve na funkce, které je zde možné nalézt. Tyto funkce lze rozdělit do několika skupin:

1. Funkce akceptující i vracející hodnoty typu __float128
2. Funkce, které na základě hodnoty __float128 vrací znaménko nebo nějaký příznak (pravdivostní hodnotu)
3. Konverzní funkce z typu __float128 na typ long int nebo long long int
4. Funkce akceptující komplexní hodnotu typu __complex128 a vracející __float128 (absolutní hodnota, reálná složka, úhel atd.)
5. Funkce akceptující i vracející komplexní hodnoty typu __complex128 (komplexní funkce)
6. Převody řetězce na __float128 a naopak převody __float128 na řetězec

6. Prototypy všech funkcí definovaných v quadmath.h

Jaké konkrétní hlavičky funkcí v hlavičkovém souboru quadmath.h nalezneme? Všechny zde definované funkce budou vypsány v této kapitole.

Začneme „obyčejnými“ funkcemi, které akceptují i vrací hodnoty typu __float128:

extern __float128 acosq (__float128) __quadmath_throw;
extern __float128 acoshq (__float128) __quadmath_throw;
extern __float128 asinq (__float128) __quadmath_throw;
extern __float128 asinhq (__float128) __quadmath_throw;
extern __float128 atanq (__float128) __quadmath_throw;
extern __float128 atanhq (__float128) __quadmath_throw;
extern __float128 atan2q (__float128, __float128) __quadmath_throw;
extern __float128 cbrtq (__float128) __quadmath_throw;
extern __float128 ceilq (__float128) __quadmath_throw;
extern __float128 copysignq (__float128, __float128) __quadmath_throw;
extern __float128 coshq (__float128) __quadmath_throw;
extern __float128 cosq (__float128) __quadmath_throw;
extern __float128 erfq (__float128) __quadmath_throw;
extern __float128 erfcq (__float128) __quadmath_throw;
extern __float128 exp2q (__float128) __quadmath_throw;
extern __float128 expq (__float128) __quadmath_throw;
extern __float128 expm1q (__float128) __quadmath_throw;
extern __float128 fabsq (__float128) __quadmath_throw;
extern __float128 fdimq (__float128, __float128) __quadmath_throw;
extern __float128 floorq (__float128) __quadmath_throw;
extern __float128 fmaq (__float128, __float128, __float128) __quadmath_throw;
extern __float128 fmaxq (__float128, __float128) __quadmath_throw;
extern __float128 fminq (__float128, __float128) __quadmath_throw;
extern __float128 fmodq (__float128, __float128) __quadmath_throw;
extern __float128 frexpq (__float128, int *) __quadmath_throw;
extern __float128 hypotq (__float128, __float128) __quadmath_throw;
extern __float128 j0q (__float128) __quadmath_throw;
extern __float128 j1q (__float128) __quadmath_throw;
extern __float128 jnq (int, __float128) __quadmath_throw;
extern __float128 ldexpq (__float128, int) __quadmath_throw;
extern __float128 lgammaq (__float128) __quadmath_throw;
extern __float128 logbq (__float128) __quadmath_throw;
extern __float128 logq (__float128) __quadmath_throw;
extern __float128 log10q (__float128) __quadmath_throw;
extern __float128 log2q (__float128) __quadmath_throw;
extern __float128 log1pq (__float128) __quadmath_throw;
extern __float128 modfq (__float128, __float128 *) __quadmath_throw;
extern __float128 nanq (const char *) __quadmath_throw;
extern __float128 nearbyintq (__float128) __quadmath_throw;
extern __float128 nextafterq (__float128, __float128) __quadmath_throw;
extern __float128 powq (__float128, __float128) __quadmath_throw;
extern __float128 remainderq (__float128, __float128) __quadmath_throw;
extern __float128 remquoq (__float128, __float128, int *) __quadmath_throw;
extern __float128 rintq (__float128) __quadmath_throw;
extern __float128 roundq (__float128) __quadmath_throw;
extern __float128 scalblnq (__float128, long int) __quadmath_throw;
extern __float128 scalbnq (__float128, int) __quadmath_throw;
extern __float128 sinhq (__float128) __quadmath_throw;
extern __float128 sinq (__float128) __quadmath_throw;
extern __float128 sqrtq (__float128) __quadmath_throw;
extern __float128 tanq (__float128) __quadmath_throw;
extern __float128 tanhq (__float128) __quadmath_throw;
extern __float128 tgammaq (__float128) __quadmath_throw;
extern __float128 truncq (__float128) __quadmath_throw;
extern __float128 y0q (__float128) __quadmath_throw;
extern __float128 y1q (__float128) __quadmath_throw;
extern __float128 ynq (int, __float128) __quadmath_throw;

Existuje i funkce počítající současně sinus i kosinus. Ta musí mít (pochopitelně) odlišnou hlavičku:

extern void sincosq (__float128, __float128 *, __float128 *) __quadmath_throw;

Dále jsou vypsány funkce, které na základě hodnoty __float128 vrací znaménko nebo nějaký příznak (pravdivostní hodnotu). Tyto funkce mají většinou obdobu pro typy float i double:

extern int signbitq (__float128) __quadmath_throw;
extern int finiteq (__float128) __quadmath_throw;
extern int isinfq (__float128) __quadmath_throw;
extern int ilogbq (__float128) __quadmath_throw;
extern int isnanq (__float128) __quadmath_throw;
extern int issignalingq (__float128) __quadmath_throw;

Následují konverzní funkce z typu __float128 na typ long int nebo long long int:

extern long long int llrintq (__float128) __quadmath_throw;
extern long long int llroundq (__float128) __quadmath_throw;
extern long int lrintq (__float128) __quadmath_throw;
extern long int lroundq (__float128) __quadmath_throw;

Čtveřice funkcí akceptuje komplexní hodnoty typu __complex128 a vrací hodnoty typu __float128. Najdeme zde podporu pro získání reálné i imaginární části komplexního čísla, výpočet jeho absolutní hodnoty nebo úhlu (argumentu):

extern __float128 cabsq (__complex128) __quadmath_throw;
extern __float128 cargq (__complex128) __quadmath_throw;
extern __float128 cimagq (__complex128) __quadmath_throw;
extern __float128 crealq (__complex128) __quadmath_throw;

Funkce začínající na „c“ a končící na „q“ akceptují i vrací hodnoty typu __complex128. Jedná se o standardní funkce rozšířené na obor komplexních čísel:

extern __complex128 cacosq (__complex128) __quadmath_throw;
extern __complex128 cacoshq (__complex128) __quadmath_throw;
extern __complex128 casinq (__complex128) __quadmath_throw;
extern __complex128 casinhq (__complex128) __quadmath_throw;
extern __complex128 catanq (__complex128) __quadmath_throw;
extern __complex128 catanhq (__complex128) __quadmath_throw;
extern __complex128 ccosq (__complex128) __quadmath_throw;
extern __complex128 ccoshq (__complex128) __quadmath_throw;
extern __complex128 cexpq (__complex128) __quadmath_throw;
extern __complex128 cexpiq (__float128) __quadmath_throw;
extern __complex128 clogq (__complex128) __quadmath_throw;
extern __complex128 clog10q (__complex128) __quadmath_throw;
extern __complex128 conjq (__complex128) __quadmath_throw;
extern __complex128 cpowq (__complex128, __complex128) __quadmath_throw;
extern __complex128 cprojq (__complex128) __quadmath_throw;
extern __complex128 csinq (__complex128) __quadmath_throw;
extern __complex128 csinhq (__complex128) __quadmath_throw;
extern __complex128 csqrtq (__complex128) __quadmath_throw;
extern __complex128 ctanq (__complex128) __quadmath_throw;
extern __complex128 ctanhq (__complex128) __quadmath_throw;

Poslední dvě funkce slouží pro převody mezi řetězcem a hodnotou typu __float128 nebo naopak. Obě tyto funkce si ukážeme v demonstračních příkladech:

extern __float128 strtoflt128 (const char *, char **) __quadmath_throw;
extern int quadmath_snprintf (char *str, size_t size, const char *format, ...) __quadmath_throw;

7. Tisk hodnot typu __float128

Pro tisk hodnot typu __float128, resp. pro jejich převod do tisknutelné textové (řetězcové) podoby slouží funkce nazvaná quadmath_snprintf. Teoreticky by bylo možné použít i standardní funkce printf, snprintf atd., ovšem (i když se to může po technologické stránce zdát triviální) by tyto funkce neměly být rozšiřovány o další datové typy nebo specifikátory velikosti (tedy Qf atd.).

V dnešním prvním demonstračním příkladu je ukázáno, jakým způsobem se hodnota typu __float128 převede na řetězcové podoby a následně vytiskne standardní funkcí puts. Povšimněte si, že vyžadujeme tisk třiceti cifer za desetinnou tečkou, což ovšem není přehnané, protože __float128 má 33 až 36 platných číslic:

/*
 * Datový typ _float128 v programovacím jazyku C
 *
 * Realizace tisku hodnot typu _float128 do řetězce i na terminál.
 */
 
#include <stdio.h>
#include <quadmath.h>
 
int main(void) {
#define N 100
    char buffer[N];
    __float128 x = 1.0/3.0;
 
    quadmath_snprintf(buffer, sizeof buffer, "%.30Qf", x);
    puts(buffer);
 
    return 0;
}

Při překladu nesmíme zapomenout na přilinkování knihovny quadmath:

$ gcc float_128_printing.c -lquadmath

Příklad by měl být nyní bez problému spustitelný:

$ ./a.out 
 
0.333333333333333314829616256247

Výsledek je nepřesný z toho důvodu, že 1.0/3.0 je vypočten v přesnosti double. Musíme tedy provést nepatrnou úpravu a použít konstanty s modifikátorem q:

/*
 * Datový typ _float128 v programovacím jazyku C
 *
 * Realizace tisku hodnot typu _float128 do řetězce i na terminál.
 */
 
#include <stdio.h>
#include <quadmath.h>
 
int main(void) {
#define N 100
    char buffer[N];
    __float128 x = 1.0q/3.0q;
 
    quadmath_snprintf(buffer, sizeof buffer, "%.30Qf", x);
    puts(buffer);
 
    return 0;
}

Nyní již bude výsledek přesný na celých třicet desetinných míst:

$ ./a.out 
 
0.333333333333333333333333333333

8. Konstanty definované v hlavičkových souborech math.h a quadmath.h

V hlavičkovém souboru quadmath.h nalezneme i několik numerických konstant s hodnotami typu __float128. Tyto konstanty mají své obdoby ve standardní hlavičkovém souboru math.h, jediný rozdíl spočívá v přidaném znaku q na jejich konci:

V dnešním druhém demonstračním příkladu vytiskneme druhou odmocninu dvojky, přičemž tuto hodnotu uložíme do proměnných typu _Float16, float/single, double i __float128. Počet desetinných míst při tisku záměrně „přeženeme“:

/*
 * Datový typ __float128 v programovacím jazyku C
 *
 * Konstanty definované v hlavičkových souborech math.h a quadmath.h
 */
 
#include <stdio.h>
#include <math.h>
#include <quadmath.h>
 
int main(void) {
#define N 100
    char buffer[N];
     
    _Float16 x16 = M_SQRT2;
    float x32 = M_SQRT2;
    double x64 = M_SQRT2;
    __float128 x128 = M_SQRT2q;
 
    printf("%.*f\n", N-10, (double)x16);
    printf("%.*f\n", N-10, x32);
    printf("%.*f\n", N-10, x64);
    quadmath_snprintf(buffer, sizeof buffer, "%.*Qf", N-10, x128);
    puts(buffer);
 
    return 0;
}

Následují výsledky vytištěné po překladu a spuštění tohoto demonstračního příkladu. Z těchto výsledků je patrné, kolik číslic jednotlivé formáty dokážou zachovat (i když nepřesně – ne všechny číslice jsou plně platné!):

1.414062500000000000000000000000000000000000000000000000000000000000000000000000000000000000
1.414213538169860839843750000000000000000000000000000000000000000000000000000000000000000000
1.414213562373095145474621858738828450441360473632812500000000000000000000000000000000000000
1.414213562373095048801688724209697984347246389158627416255512374049674740456028487756157119

9. Aritmetické operace a speciální hodnoty typu __float128

V dalším demonstračním příkladu je ukázáno, že GCC plně podporuje provádění aritmetických operací s operandy, které jsou typu __float128. Nemusí se tedy volat knihovní funkce, což by nebylo příliš přehledné:

/*
 * Datový typ __float128 v programovacím jazyku C
 *
 * Základní aritmetické operace s hodnotami typu __float128
 */
 
#include <stdio.h>
#include <quadmath.h>
 
extern __float128 add128(__float128 x, __float128 y) {
    return x+y;
}
 
extern __float128 sub128(__float128 x, __float128 y) {
    return x-y;
}
 
extern __float128 mul128(__float128 x, __float128 y) {
    return x*y;
}
 
extern __float128 div128(__float128 x, __float128 y) {
    return x/y;
}
 
extern void print_float128(__float128 x) {
#define N 100
    char buffer[N];
 
    quadmath_snprintf(buffer, sizeof buffer, "%.*Qf", N-10, x);
    puts(buffer);
}
 
int main(int argc, char **argv) {
    __float128 x = 1.0q/3.0q;
    __float128 y = 2.0q/3.0q;
    __float128 z;
 
    z = add128(x, y);
    print_float128(z);
 
    z = sub128(x, y);
    print_float128(z);
 
    z = mul128(x, y);
    print_float128(z);
 
    z = div128(x, y);
    print_float128(z);
 
    return 0;
}

Výsledky všech čtyř operací s hodnotami 1/3 a 2/3 (tyto hodnoty pochopitelně není možné reprezentovat zcela přesně):

1.000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
-0.333333333333333333333333333333333317283917130106367891200183811792272345515819598205098373
0.222222222222222222222222222222222211522611420070911927466789207861514897010546398803398915
0.500000000000000000000000000000000000000000000000000000000000000000000000000000000000000000

Z výpisu výsledku překladu do assembleru je patrné, že výpočty všech čtyř operací jsou realizovány nějakou subrutinou, nikoli přímo instrukcí (překlad pro platformu x86–64):

add128:
        sub     rsp, 8
        call    __addtf3
        add     rsp, 8
        ret
 
sub128:
        sub     rsp, 8
        call    __subtf3
        add     rsp, 8
        ret
 
mul128:
        sub     rsp, 8
        call    __multf3
        add     rsp, 8
        ret
 
div128:
        sub     rsp, 8
        call    __divtf3
        add     rsp, 8
        ret

I pouze při použití základních aritmetických operací mohou být výsledkem speciální hodnoty typu kladné nekonečno, záporné nekonečno nebo dokonce NaN (Not a Number). Toto – zcela korektní a normou popsané – chování si ověříme v dalším demonstračním příkladu, ve kterém budeme provádět dělení nulou, součet nekonečen, rozdíl nekonečen, podíl 0/0 a dokonce i podíl nekonečen. Nakonec se pokusíme o výpočet odmocniny ze záporné hodnoty:

/*
 * Datový typ __float128 v programovacím jazyku C
 *
 * Speciální hodnoty typu __float128.
 */
 
#include <stdio.h>
#include <quadmath.h>
 
extern void print_float128(__float128 x) {
#define N 100
    char buffer[N];
 
    quadmath_snprintf(buffer, sizeof buffer, "%.*Qf", N-10, x);
    puts(buffer);
}
 
int main(int argc, char **argv) {
    __float128 x = 1.0q;
    __float128 y = 0.0q;
    __float128 z;
    __float128 w;
 
    z = x / y;
    print_float128(z);
 
    z = -x / y;
    print_float128(z);
 
    w = x/y + x/y;
    print_float128(w);
 
    w = x/y - x/y;
    print_float128(w);
 
    w = y/y;
    print_float128(w);
 
    w = w/w;
    print_float128(w);
 
    w = sqrtq(-2.0);
    print_float128(w);
 
    return 0;
}

Výsledky budou odpovídat normě IEEE 754 (poznámky jsou dopsány ručně):

inf   (1/0)
-inf  (-1/0)
inf   (inf+inf)
-nan  (inf-inf - zde může být výsledek i kladné NaN)
-nan  (0/0)
-nan  (inf/inf)
-nan  (√-2)

10. Volání funkce definované v hlavičkovém souboru quadmath.h

V šesté kapitole byly vypsány hlavičky všech funkcí, které jsou dostupné v hlavičkovém souboru quadmath.h. Kvůli omezenému rozsahu článku si ukážeme volání jen jediné této funkce, a to konkrétně funkce nazvané sqrtq, která dokáže vypočítat druhou odmocninu svého parametru. Pro zajímavost vypočtenou hodnotu porovnáme s odmocninou hodnoty typu float, double a float16 (tyto výsledky by měly být méně přesné):

/*
 * Datový typ __float128 v programovacím jazyku C
 *
 * Volání funkce definované v hlavičkovém souboru quadmath.h
 */
 
#include <stdio.h>
#include <math.h>
#include <quadmath.h>
 
int main(void) {
#define N 100
    char buffer[N];
 
    _Float16 x16 = sqrtf(2.0f);
    float x32 = sqrtf(2.0f);
    double x64 = sqrt(2.0);
    __float128 x128 = sqrtq(2.0q);
 
    printf("%.*f\n", N-10, (double)x16);
    printf("%.*f\n", N-10, x32);
    printf("%.*f\n", N-10, x64);
    quadmath_snprintf(buffer, sizeof buffer, "%.*Qf", N-10, x128);
    puts(buffer);
 
    return 0;
}

A takto vypadají vypočtené výsledky:

1.414062500000000000000000000000000000000000000000000000000000000000000000000000000000000000
1.414213538169860839843750000000000000000000000000000000000000000000000000000000000000000000
1.414213562373095145474621858738828450441360473632812500000000000000000000000000000000000000
1.414213562373095048801688724209697984347246389158627416255512374049674740456028487756157119

Je pro porovnání – hodnota druhé odmocniny dvou s přesností jednoho sta číslic vypadá takto:

1.41421356237309504880168872420969807856967187537694807317667973799073247846210703885038753432764157

Porovnejme tuto hodnotu číslici po číslici s posledním vypočteným výsledkem:

1.414213562373095048801688724209697984347246389158627416255512374049674740456028487756157119
                                  ^
1.41421356237309504880168872420969807856967187537694807317667973799073247846210703885038753432764157
                                  ^

Číslice se začínají odlišovat na třicátém čtvrtém místě – a to přesně odpovídá přesnosti tohoto datového typu!

11. Rozdílné výsledky operací s hodnotami s plovoucí řádovou čárkou

Již v demonstračním příkladu, který na standardní výstup prováděl tisk hodnot typu float128, jsme narazili na problém zápisu konstantních výrazů, které se mají vyhodnotit právě na hodnotu float128. Numerické konstanty v těchto výrazech mohou mít postfix „q“. V takovém případě překladač provede korektně výpočet s touto přesností.

V dalším demonstračním příkladu se ve zdrojových kódech vyhodnotí (už v době překladu) hodnota 1/3 a výsledek vyhodnocení se uloží do proměnných různých typů. Navíc se v případě hodnot float128 explicitně provede výpočet konstantního výrazu s vyšší přesností:

/*
 * Datový typ __float128 v programovacím jazyku C
 *
 * Rozdílné výsledky operací s hodnotami s plovoucí řádovou čárkou.
 */
 
#include <stdio.h>
#include <quadmath.h>
 
int main(void) {
#define N 100
    char buffer[N];
    _Float16 x16 = 1.0/3.0;
    float x32 = 1.0/3.0;
    double x64 = 1.0/3.0;
    __float128 x128 = 1.0q/3.0q;
 
    printf("%.*f\n", N-10, (double)x16);
    printf("%.*f\n", N-10, x32);
    printf("%.*f\n", N-10, x64);
    quadmath_snprintf(buffer, sizeof buffer, "%.*Qf", N-10, x128);
    puts(buffer);
 
    return 0;
}

Výsledky budou vypadat následovně:

0.333251953125000000000000000000000000000000000000000000000000000000000000000000000000000000
0.333333343267440795898437500000000000000000000000000000000000000000000000000000000000000000
0.333333333333333314829616256247390992939472198486328125000000000000000000000000000000000000
0.333333333333333333333333333333333317283917130106367891200183811792272345515819598205098373

Zajímat nás bude především poslední výsledek. Ten je přesný až na třicáté třetí desetinné místo, což odpovídá počtu platných cifer formátu float128.

12. Iterativní výpočet odmocniny realizovaný s typy float, double i quadruple

Společné znaky, ale i rozdíly mezi výpočty s hodnotami typu float, double a quadruple, si ukážeme na relativně jednoduchém algoritmu. Bude se jednat o výpočet druhé odmocniny nějakého nezáporného reálného čísla. Pro výpočet bude použit iterativní algoritmus, který je známý pod jménem Heronova metoda.

První varianta výpočtu je založena na použití hodnot typu float (jednoduchá přesnost):

/*
 * Datový typ __float128 v programovacím jazyku C
 *
 * Iterativní výpočet odmocniny: realizace pro hodnoty typu float.
 */
 
#include <stdio.h>
 
int main(void) {
    float s = 2.0f;
    float x = 1.0f;
    int i;
 
    for (i=0; i<10; i++) {
        printf("%3d  %.40f\n", i, x);
        x = 1.0f/2.0f*(x+s/x);
    }
 
    return 0;
}

Jak je z dalšího výpisu patrné, dostaneme již po pěti iteracích stabilní (i když pochopitelně nepřesný) výsledek:

  0  1.0000000000000000000000000000000000000000
  1  1.5000000000000000000000000000000000000000
  2  1.4166667461395263671875000000000000000000
  3  1.4142156839370727539062500000000000000000
  4  1.4142135381698608398437500000000000000000
  5  1.4142135381698608398437500000000000000000
  6  1.4142135381698608398437500000000000000000
  7  1.4142135381698608398437500000000000000000
  8  1.4142135381698608398437500000000000000000
  9  1.4142135381698608398437500000000000000000

Druhá varianta při výpočtu využívá hodnoty typu double (dvojnásobná přesnost):

/*
 * Datový typ __float128 v programovacím jazyku C
 *
 * Iterativní výpočet odmocniny: realizace pro hodnoty typu double.
 */
 
#include <stdio.h>
 
int main(void) {
    double s = 2.0;
    double x = 1.0;
    int i;
 
    for (i=0; i<10; i++) {
        printf("%3d  %.40f\n", i, x);
        x = 1.0/2.0*(x+s/x);
    }
 
    return 0;
}

V tomto případě je stabilního výsledku dosaženo až po šesti iteracích:

  0  1.0000000000000000000000000000000000000000
  1  1.5000000000000000000000000000000000000000
  2  1.4166666666666665186369300499791279435158
  3  1.4142156862745096645994635764509439468384
  4  1.4142135623746898698271934335934929549694
  5  1.4142135623730949234300169337075203657150
  6  1.4142135623730949234300169337075203657150
  7  1.4142135623730949234300169337075203657150
  8  1.4142135623730949234300169337075203657150
  9  1.4142135623730949234300169337075203657150

A konečně poslední varianta výpočtu je založena na hodnotách typu __float128 (čtyřnásobná přesnost):

/*
 * Datový typ __float128 v programovacím jazyku C
 *
 * Iterativní výpočet odmocniny: realizace pro hodnoty typu __float128.
 */
 
#include <stdio.h>
#include <quadmath.h>
 
int main(void) {
#define N 50
    char buffer[N];
 
    __float128 s = 2.0q;
    __float128 x = 1.0q;
 
    int i;
 
    for (i=0; i<10; i++) {
        quadmath_snprintf(buffer, sizeof buffer, "%.40Qf", x);
        printf("%3d %s\n", i, buffer);
        x = 1.0q/2.0q*(x+s/x);
    }
 
    return 0;
}

Při tomto výpočtu dojdeme je stabilnímu výsledku až po sedmi iteracích:

  0 1.0000000000000000000000000000000000000000
  1 1.5000000000000000000000000000000000000000
  2 1.4166666666666666666666666666666665382713
  3 1.4142156862745098039215686274509803846042
  4 1.4142135623746899106262955788901348045363
  5 1.4142135623730950488016896235025302739130
  6 1.4142135623730950488016887242096981769402
  7 1.4142135623730950488016887242096981769402
  8 1.4142135623730950488016887242096981769402
  9 1.4142135623730950488016887242096981769402

1.41421356237309504880168872420969807856967187537694807317667973799073247846210703885038753432764157
                                   ^
1.4142135623730950488016887242096981769402

13. Převod řetězce na hodnotu typu __float128

Ještě nám zbývá popis poslední funkce, která přímo souvisí s hodnotami typu __float128. Jedná se o funkci nazvanou strtoflt128, která umožňuje převod (parsing) řetězce právě na hodnotu typu __float128. Tato funkce akceptuje jako své parametry dva ukazatele, ovšem většinou je nutné předat jen první ukazatel – a to konkrétně ukazatel na řetězec se zapsanou numerickou hodnotou. Použití této funkce je tedy v praxi velmi snadné:

/*
 * Datový typ __float128 v programovacím jazyku C
 *
 * Převod z řetězce na hodnoty typu __float128.
 */
 
#include <stdio.h>
#include <quadmath.h>
 
void print_fp128_value(__float128 x) {
#define N 50
    char buffer[N];
    quadmath_snprintf(buffer, sizeof buffer, "%.*Qf", N-5, x);
    puts(buffer);
}
 
int main(void) {
 
    __float128 x;
 
    x = strtoflt128 ("1.0", NULL);
    print_fp128_value(x);
 
    x = strtoflt128 ("0.1", NULL);
    print_fp128_value(x);
 
    x = strtoflt128 ("0.00000000000000000000000000000000001", NULL);
    print_fp128_value(x);
 
    x = strtoflt128 ("0.000000000000000000000000000000000000000000001", NULL);
    print_fp128_value(x);
 
    x = strtoflt128 ("0.0000000000000000000000000000000000000000000001", NULL);
    print_fp128_value(x);
 
    return 0;
}

Z vytištěných výsledků je patrné to, že poslední řetězec obsahoval tak malou hodnotu, že byla převedena na nulu. A pochopitelně nebylo možné korektně reprezentovat hodnotu 0,1, což je typický problém naprosto všech FP formátů založených na exponentu se základem 2. Výsledky ostatních konverzí odpovídají očekávání:

1.000000000000000000000000000000000000000000000
0.100000000000000000000000000000000004814824861
0.000000000000000000000000000000000010000000000
0.000000000000000000000000000000000000000000001
0.000000000000000000000000000000000000000000000

14. Podpora pro práci s komplexními čísly se 128bitovými reálnými a imaginárními složkami

Prozatím jsme si popsali funkce z hlavičkového souboru quadmath.h a knihovny quadmath.so, které slouží ke zpracování reálných čísel (resp. podmnožiny reálných čísel) reprezentovaných ve formátu float128. Ve skutečnosti ovšem tato knihovna podporuje i práci s čísly komplexními. Typ těchto komplexních čísel se jmenuje __complex128 a interně se jedná o dvojici hodnot typu __complex128. V dalších čtyřech kapitolách si ukážeme práci s těmito komplexními hodnotami. Připomeňme si pouze, že v jazyku C je možné zapisovat komplexní číslo formou reálné složky připočtené ke složce imaginární (ta má postfix i nebo j).

15. Tisk numerických hodnot typu __complex128

Přímý tisk komplexních hodnot typu __complex128 sice není možný, ovšem díky tomu, že dokážeme vytisknout hodnoty typu __float128, je možné postupně vytisknout reálnou složku následovanou složkou imaginární. Pro získání reálné složky se používá funkce crealq a pro získání složky imaginární funkce cimagq:

/*
 * Datový typ _complex128 v programovacím jazyku C
 *
 * Realizace tisku hodnot typu _complex128 do řetězce i na terminál.
 */
 
#include <stdio.h>
#include <quadmath.h>
 
int main(void) {
#define N 100
    char buffer[N];
    __complex128 z = 1.0q/3.0q + 1.0iq/9.0q;
 
    quadmath_snprintf(buffer, sizeof buffer, "%.30Qf", crealq(z));
    puts(buffer);
    quadmath_snprintf(buffer, sizeof buffer, "%.30Qf", cimagq(z));
    puts(buffer);
 
    return 0;
}

Výsledkem činnosti tohoto demonstračního příkladu je řádek s reálnou složkou následovaný řádkem se složkou imaginární:

0.333333333333333333333333333333
0.111111111111111111111111111111

16. Základní aritmetické operace s hodnotami typu __complex128

S komplexními proměnnými a hodnotami je možné provádět všechny základní aritmetické operace. To si ostatně ověříme překladem a spuštěním dalšího demonstračního příkladu, který všechny čtyři základní aritmetické operace provede a následně vypíše jejich výsledky (komplexní čísla):

/*
 * Datový typ __complex128 v programovacím jazyku C
 *
 * Základní aritmetické operace s hodnotami typu __complex128
 */
 
#include <stdio.h>
#include <quadmath.h>
 
extern __complex128 add128(__complex128 x, __complex128 y) {
    return x+y;
}
 
extern __complex128 sub128(__complex128 x, __complex128 y) {
    return x-y;
}
 
extern __complex128 mul128(__complex128 x, __complex128 y) {
    return x*y;
}
 
extern __complex128 div128(__complex128 x, __complex128 y) {
    return x/y;
}
 
extern void print_complex128(__complex128 x) {
#define N 100
    char buffer[N];
 
    quadmath_snprintf(buffer, sizeof buffer, "%.*Qf", N-10, crealq(x));
    puts(buffer);
    quadmath_snprintf(buffer, sizeof buffer, "%.*Qf", N-10, cimagq(x));
    puts(buffer);
 
    putchar('\n');
}
 
int main(int argc, char **argv) {
    __complex128 x = 1.0q/3.0q + 0.5iq;
    __complex128 y = 2.0q/3.0q + 0.25iq;
    __complex128 z;
 
    puts("x:");
    print_complex128(x);
 
    puts("y:");
    print_complex128(y);
 
    z = add128(x, y);
    puts("x+y:");
    print_complex128(z);
 
    z = sub128(x, y);
    puts("x-y:");
    print_complex128(z);
 
    z = mul128(x, y);
    puts("x*y:");
    print_complex128(z);
 
    z = div128(x, y);
    puts("x/y:");
    print_complex128(z);
 
    return 0;
}

Po překladu tohoto příkladu se vypíše hodnota x, y a výsledek všech čtyř operací. První řádek obsahuje reálnou složku, druhý řádek složku imaginární:

x:
0.333333333333333333333333333333333317283917130106367891200183811792272345515819598205098373
0.500000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
 
y:
0.666666666666666666666666666666666634567834260212735782400367623584544691031639196410196746
0.250000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
 
x+y:
1.000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0.750000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
 
x-y:
-0.333333333333333333333333333333333317283917130106367891200183811792272345515819598205098373
0.250000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
 
x*y:
0.097222222222222222222222222222222211522611420070911927466789207861514897010546398803398915
0.416666666666666666666666666666666634567834260212735782400367623584544691031639196410196746
 
x/y:
0.684931506849315068493150684931506760933351867161094414554436881376623053388486006554103369
0.493150684931506849315068493150684952612930897394365512942224028327970614116182446196035016

17. Funkce akceptující parametr typu __complex128 a vracející __float128

V knihovně quadmath je definována čtveřice funkcí, které jako svůj parametr akceptují hodnoty typu __complex128 (tedy komplexní číslo) a vrací hodnotu typu __float128, tedy číslo reálné. Použití těchto funkcí je snadné:

/*
 * Datový typ __complex128 v programovacím jazyku C
 *
 * Konverze hodnot typu __complex128 na hodnoty typu __float128.
 */
 
#include <stdio.h>
#include <quadmath.h>
 
extern void print_float128(__float128 x) {
#define N 20
    char buffer[N];
 
    quadmath_snprintf(buffer, sizeof buffer, "%.*Qf", N-10, x);
    puts(buffer);
 
    putchar('\n');
}
 
extern void print_complex128(__complex128 x) {
#define N 20
    char buffer[N];
 
    quadmath_snprintf(buffer, sizeof buffer, "%.*Qf", N-10, crealq(x));
    printf("%s + i", buffer);
    quadmath_snprintf(buffer, sizeof buffer, "%.*Qf", N-10, cimagq(x));
    puts(buffer);
 
    putchar('\n');
}
 
int main(int argc, char **argv) {
    __complex128 x = 3.0q + 4.0iq;
    __float128 y;
 
    puts("x:");
    print_complex128(x);
 
    puts("crealq:");
    y = crealq(x);
    print_float128(y);
 
    puts("cimagq:");
    y = cimagq(x);
    print_float128(y);
 
    puts("cabsq:");
    y = cabsq(x);
    print_float128(y);
 
    puts("cargq:");
    y = cargq(x);
    print_float128(y);
 
    return 0;
}

Z výpisu je patrné, že tyto funkce postupně vrací reálnou složku, imaginární složku, absolutní hodnotu (velikost) komplexního čísla a nakonec jeho úhel:

x:
3.0000000000 + i4.0000000000
 
crealq:
3.0000000000
 
cimagq:
4.0000000000
 
cabsq:
5.0000000000
 
cargq:
0.9272952180

18. Funkce akceptující parametr typu __complex128 a vracející __complex128

V poslední skupině funkcí nabízených knihovnou quadmath nalezneme funkce, které jako svůj parametr akceptují komplexní číslo a současně vrací odlišné komplexní číslo jako svůj výsledek. Těchto funkcí existuje celá řada. My se v demonstračním příkladu zaměříme na funkci určenou pro výpočet odmocniny komplexního čísla (a tedy i záporného čísla reálného – to má nulovou imaginární složku):

/*
 * Datový typ __complex128 v programovacím jazyku C
 *
 * Volání funkcí akceptujících i vracejících hodnoty typu __complex128.
 */
 
#include <stdio.h>
#include <quadmath.h>
 
extern void print_complex128(__complex128 x) {
#define N 20
    char buffer[N];
 
    quadmath_snprintf(buffer, sizeof buffer, "%.*Qf", N-10, crealq(x));
    printf("%s + i", buffer);
    quadmath_snprintf(buffer, sizeof buffer, "%.*Qf", N-10, cimagq(x));
    puts(buffer);
 
    putchar('\n');
}
 
int main(int argc, char **argv) {
    __complex128 x = 2.0q + 0.0iq;
    __complex128 y = -2.0q + 0.0iq;
    __complex128 z = 0.0q + 2.0iq;
    __complex128 w;
 
    puts("x:");
    print_complex128(x);
 
    puts("y:");
    print_complex128(y);
 
    puts("z:");
    print_complex128(z);
 
    w = csqrtq(x);
    puts("√x:");
    print_complex128(w);
 
    w = csqrtq(y);
    puts("√y:");
    print_complex128(w);
 
    w = csqrtq(z);
    puts("√z:");
    print_complex128(w);
 
    return 0;
}

Výsledky ukazují, že výpočty jsou skutečně korektní, a to i pro záporná reálná čísla (v reálném oboru nemá tento vstup řešení):

x:
2.0000000000 + i0.0000000000
 
y:
-2.0000000000 + i0.0000000000
 
z:
0.0000000000 + i2.0000000000
 
√x:
1.4142135624 + i0.0000000000
 
√y:
0.0000000000 + i1.4142135624
 
√z:
1.0000000000 + i1.0000000000

19. Repositář s demonstračními příklady

Demonstrační příklady, s nimiž jsme se v dnešním článku seznámili a které jsou určeny pro překlad s využitím překladače GCC, jsou dostupné, jak je v tomto seriálu zvykem, na GitHubu:

20. Odkazy na Internetu

1. Quadruple-precision floating-point format
   https://en.wikipedia.org/wi­ki/Quadruple-precision_floating-point_format
2. Exploring Quadruple Precision Floating Point Numbers in GCC and ICC
   https://www.nsc.liu.se/~pla/blog/2013/10/17­/quadruple-precision/
3. GCC: Additional Floating Types
   https://gcc.gnu.org/online­docs/gcc/Floating-Types.html
4. GCC Quad-Precision Math Library Application Programming Interface
   https://gcc.gnu.org/online­docs/libquadmath/index.html#SEC_Con­tents
5. Mathematical constants in C provided by math.h
   https://www.lemoda.net/c/maths-constants/
6. Square root of 2
   https://en.wikipedia.org/wi­ki/Square_root_of₂
7. On-Line Encyclopedia of Integer Sequences
   https://en.wikipedia.org/wiki/On-Line_Encyclopedia_of_Integer_Sequences
8. Hero of Alexandria
   https://en.wikipedia.org/wi­ki/Hero_of_Alexandria
9. Norma IEEE 754 a příbuzní: formáty plovoucí řádové tečky
   https://www.root.cz/clanky/norma-ieee-754-a-pribuzni-formaty-plovouci-radove-tecky/
10. IEEE-754 Floating-Point Conversion
    http://babbage.cs.qc.cuny.edu/IEEE-754.old/32bit.html
11. Small Float Formats
    https://www.khronos.org/o­pengl/wiki/Small_Float_For­mats
12. Binary-coded decimal
    https://en.wikipedia.org/wiki/Binary-coded_decimal
13. Why Intel is betting on BFLOAT16 to be a game changer for deep learning training? Hint: Range trumps Precision
    https://hub.packtpub.com/why-intel-is-betting-on-bfloat16-to-be-a-game-changer-for-deep-learning-training-hint-range-trumps-precision/
14. Art of Assembly language programming: The 80×87 Floating Point Coprocessors
    https://courses.engr.illi­nois.edu/ece390/books/arto­fasm/CH14/CH14–3.html
15. Art of Assembly language programming: The FPU Instruction Set
    https://courses.engr.illi­nois.edu/ece390/books/arto­fasm/CH14/CH14–4.html
16. INTEL 80387 PROGRAMMER'S REFERENCE MANUAL
    http://www.ragestorm.net/dow­nloads/387intel.txt
17. Floating-Point Formats
    http://www.quadibloc.com/com­p/cp0201.htm
18. Print __float128, without using quadmath_snprintf
    https://stackoverflow.com/qu­estions/23654693/print-float128-without-using-quadmath-snprintf
19. quadmath and argument type for scanf
    https://stackoverflow.com/qu­estions/21847735/quadmath-and-argument-type-for-scanf
20. Octuple-precision floating point on Apple G4 (archivováno)
    http://images.apple.com/ca/ac­g/pdf/oct3a.pdf
21. Arbitrary-precision arithmetic
    https://en.wikipedia.org/wi­ki/Arbitrary-precision_arithmetic
22. A Guide to GCC's Additional Floating Types and Portable Alternatives
    https://runebook.dev/en/doc­s/gcc/floating-types
23. long double
    https://en.wikipedia.org/wi­ki/Long_double
24. Rust: floating point types
    https://doc.rust-lang.org/book/ch03–02-data-types.html#floating-point-types
25. Go: numeric types
    https://go.dev/ref/spec#Numeric_types