Fixed point math e Mandelbrot

Un esempio di utilizzo della matematica a virgola fissa, un metodo di rappresentazione dei numeri reali con un numero prefissato di cifre decimali utilizzato in particolar modo per calcolatori privi di coprocessore matematico (per esempio il software Fractint usava questo tipo di approccio per il calcolo di frattali su processori 386 e 486 all'inizio degli anni 90)

Per effettuare un calcolo dell'insieme di Mandelbrot si puo' usare attualmente la libreria FIXEDPTC (Fixed Point Math Library for C)  (rispetto ad altre librerie, come FPMLib o LibFoixMath , e' possibile settare il proprio formato dati secondo le proprie necessita')

Visto che il calcolo e' tutto compreso per valori inferiori a 8 puo' essere utile, per massimizzare la precisione, usare una rappresentazione 4.28 (4 digit per la parte intera e 28 per la parte frazionaria) su un calcolatore a 32 bit. In questo modo il piu' piccolo numero calcolabile e' 0.00000000372529. Per la rappresentazione grafica e' stata usata la libreria PngWriter gia' usata qui
L'uso della libreria e' abbastanza autoesplicativo...l'unico aspetto un po' confuso e' che per effettuare dei confronti in un ciclo if...then si deve convertire il dato da formato a virgola fissa in float.

mand.c
------------------------------------------------------
#include <pngwriter.h>

#define FIXEDPT_WBITS 4
#define WIDTH 640
#define    HEIGHT 480

#include "fixedptc.h"

int main() {
int     j,k,i;
float    test;
fixedpt Re,Im,DRe,DIm;
fixedpt X,Y,DUE;
fixedpt XN,YN,YNN;
fixedpt    A,B;

pngwriter png(WIDTH,HEIGHT,0,"mandelbrot.png");

Re = fixedpt_rconst(-2.00390625); //finestra reale tra -2 e +0.5
Im = fixedpt_rconst(-1.205); //finestra immaginaria tra -1.2 e 1.2
DRe = fixedpt_rconst(0.00390625); //2.5/640
DIm = fixedpt_rconst(0.005); // 2.4/480
DUE = fixedpt_rconst(2.0);

A = Re;

for (j=0;j<WIDTH;j++)
    {
    A = fixedpt_add(A,DRe);
    B = Im;
    for (k=0;k<HEIGHT;k++)
    {
    B = fixedpt_add(B,DIm);

    X = fixedpt_rconst(0.0);
    Y = fixedpt_rconst(0.0);

    for (i=0;i<=255;i++)
        {
        XN=fixedpt_sub(fixedpt_mul(X,X),fixedpt_mul(Y,Y))+A; // (x*x) - (y*y) + A
        YN=fixedpt_mul(X,Y); // x*y
        YNN=fixedpt_mul(DUE,YN); // 2*x*y
        YN=YNN + B; // 2*x*y*+B
        test = fixedpt_tofloat(fixedpt_mul(XN,XN) + fixedpt_mul(YN,YN)); //(XN*XN)+(YN*YN)
        if (test > 4.0)
            {
            png.plot(j,k,1.0,1.0,1.0);
            break;
            }
        X = XN;
        Y = YN;
        }
    }
    }
png.close();

}
------------------------------------------------------

Per compilare il sorgente si puo' usare i seguenti comandi
------------------------------------------------------
rm mandelbrot.png
g++ mand.c -o mand  -I/usr/local/include  -L/usr/local/lib -lpng -lpngwriter
./mand
display mandelbrot.png
------------------------------------------------------