1. Wstęp

W artykule zostanie zaprezentowane wyświetlanie liczb typu single (32-bitowych). Jednak algorytm jest o tyle uniwersalny, iż nie powinno być problemów z jego adaptacją dla innych formatów zmiennoprzecinkowych.

Jak powszechnie wiadomo, najmniej znaczącą cyfrę liczby całkowitej, można otrzymać jako resztę z dzielenia przez podstawę systemu — w naszym przypadku podstawą jest 10. Najbardziej znaczącą cyfrę liczby ułamkowej (mniejszej od 1), uzyskujemy poprzez przemnożenie jej przez podstawę systemu.

Ponieważ wszelkie operacje na liczbach zmiennoprzecinkowych są zawsze obarczone błędem, przeto wspomniane dzielenie i mnożenie przez 10 należy wykonać na postaci stałoprzecinkowej.

Konwersja z postaci zmiennoprzecinkowej na stałoprzecinkową polega na przesunięciu bitowemu mantysy o wartość zapisaną w wykładniku. Ponieważ wykładnik jest rzędu +-128, więc liczba stałoprzecinkowa musi posiadać 2*128+24 bity, co daje 35 bajtów. W przykładowym kodzie, dla równego rachunku, przeznaczono po 20 bajtów na część całkowitą i ułamkową.

2. Wyświetlanie części całkowitej

Jak zostało powiedziane, najmłodszą cyfrę można uzyskać poprzez podzielenie 10 — dokładniej cyfrą jest resztą z tego dzielenia. Proces ten powtarzany jest, aż do chwili, gdy dzielona liczba będzie równa zero; pseudokod przedstawia się następująco:

while x > 0:
        cyfra[i] = x mod 10
               x = x div 10
               i = i + 1

Jak jednak podzielić duże liczby stałoprzecinkowe? Można potraktować je, jako liczby o podstawie 256, wtedy każdy bajt reprezentuje cyfrę. Jeśli dana jest reszta z dzielenia ze starszej pozycji (ozn. r) oraz i-ta cyfra (ozn. a) to dzielenie przeprowadzamy na liczbie:

r ⋅ 256 + a

Dodatkowo, ponieważ dzielnik jest stały, wyniki dzielenia można prekalkulować — rozmiar dwóch tablic: wynik i reszta z dzielenia --- wyniesie 2560 bajtów (r = 0..9, a = 0..255, zatem liczba możliwych kombinacji jest równa 2559).

W przypadku tak dużych liczb stałoprzecinkowych, ilość iteracji może być dość duża. Ale ponieważ początkowo ilość bitów znaczących jest względnie niewielka więc można wyznaczyć przynajmniej najstarszy znaczący bajt. Gdy wynik dzielenia dla tego bajtu będzie równy 0, to zostanie zawężony przedział obliczeń.

3. Konwersja części ułamkowej

Uzyskanie najstarszej cyfry wymaga przemnożenia przez 10, oczywiście możliwe jest prekalkulowanie wyników, ale obecnie mnożenie nie jest kosztowne. Uwagi podobne jak w punkcie 2.

4. Wartości specjalne

Format IEEE-754 definiuje następujące wartości specjalne:

Jeśli liczba jest NaN, to najstarszy — 23 bit mantysy — służy rozróżnieniu pomiędzy QNaN i SNaN. Jeśli jest ustawiony to liczba jest QNaN (Quiet NaN) jej wystąpienie nie powoduje generacji wyjątku; w przeciwnym razie mamy do czynienia z SNaN (Significant NaN) — generuje wyjątek.

Wyjaśnienia wymaga jeszcze liczba zdenormalizowana. Mantysa liczby zdenormalizowanej jest mniejsza niż jeden, wobec czego nie uzupełnia się bitu jedności.

5. Funkcja

// -- float2string.cc --
#include <stdio.h>
#include <string.h>
#include <assert.h>

typedef unsigned long long int  qword; // 64 bity
typedef unsigned int            dword; // 32 bity
typedef unsigned short int       word; // 16 bity
typedef unsigned char            byte; //  8 bitów

/*
 * byte quotient [10*256] = {for (i=0; i<2560; i++) element[i]=i/10};
 * byte remainder[10*256] = {for (i=0; i<2560; i++) element[i]=i%10};
 */
#include "precalc.c"

/*
 * funkcja dzieli liczbę o podstawie 256 przez 10.
 *
 * 'start' i 'end' (zmienna 'end' może zostać zmieniona)
 * określają zakres tablicy 'number' na jakim przeprowadzane
 * są obliczenia
 *
 * zwracana jest reszta z dzielenia
 */
byte div_by_10(byte number[], int& start, int& end)
{
        assert(end >= start);

        int  i   = end;
        word rem = 0;   // reszta

        // dzielenie
        while (i >= start)
                {
                 word digit = number[i];
                 word tmp   = 256*rem + digit;

                 number[i]  = quotient[tmp];  // tmp/10;
                 rem        = remainder[tmp]; // tmp%10;
                 i--;
                }
        // zawężenie przedziału obliczeń
        while (number[end]==0 && start <= end)
                end--;

        return rem;
}

/*
 * funkcja mnoży liczbę o podstawie 256 przez 10.
 *
 * 'start' i 'end' (zmienna 'start' może zostać zmieniona)
 * określają zakres tablicy 'number' na jakim przeprowadzane
 * są obliczenia
 *
 * zwracane jest przeniesienie
 */
byte mul_by_10(byte number[], int& start, int& end)
{
        assert(end >= start);

        int i = start;
        dword carry = 0;        // przeniesienie

        // mnożenie
        while (i <= end)
                {
                 word digit = number[i];
                 word tmp   = 10*digit + carry;

                 number[i]  = tmp &  0xff;
                 carry      = tmp >> 8;
                 i++;
                }
        // zawężenie przedziału obliczeń
        while (number[start]==0 && start <= end)
                start++;

        return carry;
}

char* float2string(char* S, float x)
{
 static byte fixedpoint[40]; //  0..19 - część ułamkowa
                             // 20..39 - część całkowita
#define bias 127
#define dot_position (20*8)  /* pozycja kropki dziesiętnej (w bitach)  */

 // 1. wydzielenie pól bitowych
 dword bin  =  *(dword*)&x;
 dword sign = (bin >> 31);              // znak
 dword exp  = (bin >> 23) & 0xff;       // wykładnik+bias
 dword mantissa = (bin & 0x007fffff);   // mantysa

 // 2. wartości specjalne
 // 2a. zero
 if (exp==0x00 && mantissa==0)
        {
         if (sign) strcpy(S, "-0.0");
              else strcpy(S, "+0.0");
         return S;
        }
 // 2b. nieskończoność
 if (exp==0xff && mantissa==0)
        {
         if (sign) strcpy(S, "-inf");
              else strcpy(S, "+inf");
         return S;
        }
 // 2c. NaN -- rozróżnienie między QNaN i SNaN
 if (exp==0xff && mantissa!=0)
        {
         if (mantissa & 0x800000)
                strcpy(S, "QNaN");
         else
                strcpy(S, "SNaN");
         return S;
        }

 // 3. uzupełnienie części całkowitej
 //    (jeśli liczba nie jest zdenormalizowana)
 if (exp==0x00 && mantissa!=0)
        ;
 else
        mantissa |= 0x800000;

 // 4. wstawienie znaku liczby
 if (sign)
        {
         S[0] = '-';
         S[1] = '\0';
        }
 else
        {
         //S[0] = '+'; S[1]='\0'
         S[0] = '\0';
        }

 // 5. zamiana liczby zmiennoprzecinkowej na fixed-point

 memset(fixedpoint, sizeof(fixedpoint), 0); // zerowanie

 // exp-23, ponieważ mantysa jest przesunięta o 23 bity w prawo
 int position = (exp-23) + dot_position - bias;

 int shift    = position & 0x7; // przesunięcie w obrębie bajtu
 int index    = position >> 3;  // pierwszy bajt

 // "wstawienie" mantysy do liczby stałoprzecinkowej
 // ponieważ shift=0..7, wobec czego liczba bitów znaczących
 // wynosi zaledwie 32
 mantissa <<= shift;
 fixedpoint[index+0] = (mantissa >> 0x00) & 0xff;
 fixedpoint[index+1] = (mantissa >> 0x08) & 0xff;
 fixedpoint[index+2] = (mantissa >> 0x10) & 0xff;
 fixedpoint[index+3] = (mantissa >> 0x18) & 0xff;

 // 6. konwersja
 static char tmp[128];

 if (exp >= 0+bias) // liczba posiada część całkowitą
        {
         memset(tmp, sizeof(tmp), 0);
         // zakres w którym zostaną przeprowadzone obliczenia
         int start = 20;
         int end   = index+3;
         //
         int idx = sizeof(tmp)-1;
         while (1)
                {
                 tmp[--idx] = div_by_10(fixedpoint, start, end) + '0';
                 if (start > end)
                        break;
                }
         strcat(S, &tmp[idx]);
        }
 else
         strcat(S, "0");

 strcat(S, ".");

 if (exp < 23+bias) // liczba posiada część ułamkową
        {
         memset(tmp, sizeof(tmp), 0);
         // zakres w którym zostaną przeprowadzone obliczenia
         int start = index;
         int end   = 19;
         //
         int idx   = 0;
         while (1)
                {
                 tmp[idx++] = mul_by_10(fixedpoint, start, end) + '0';
                 if (start > end)
                        break;
                }
         strcat(S, tmp);
        }
 else
         strcat(S, "0");

 return S;
}

Włączany jest plik precalc.c, który zawiera tablice quotient i remainder (wynik dzielenia i reszta). Oczywiście ich zawartość może być liczona na początku programu, co z oczywistych względów nie zostało pokazane. Są one wyznaczane wg algorytmu podanego w komentarzu; poniżej skrypt w bashu który generuje ten plik.

#!/bin/bash
file="precalc.c"

echo "byte quotient[2560] = {" > $file
for ((i=0 ; i<256 ; i++))
do
    if [ $i -ne 255 ]
    then
        echo $i, $i, $i, $i, $i, $i, $i, $i, $i, $i, >> $file
    else
        echo $i, $i, $i, $i, $i, $i, $i, $i, $i, $i  >> $file
    fi
done
echo "};" >> $file
echo ""   >> $file

echo "byte remainder[2560] = {" >> $file
for ((i=0 ; i<256 ; i++))
do
    if [ $i -ne 255 ]
    then
        echo "0, 1, 2, 3, 4, 5, 6, 7, 8, 9," >> $file
    else
        echo "0, 1, 2, 3, 4, 5, 6, 7, 8, 9"  >> $file
    fi
done
echo "};" >> $file