This site is updated infrequently. For up-to-date information, please visit the new OCaml website at ocaml.org.

Multiplication of matrix in C and OCaml
[ Home ] [ Index: by date | by threads ]
[ Search: ]

[ Message by date: previous | next ] [ Message in thread: previous | next ] [ Thread: previous | next ]
 Date: 2007-02-07 (23:42) From: Frédéric_Gava Subject: Multiplication of matrix in C and OCaml
```Dear All,

I compared multiplication of matrix in C and OCaml and I was a little
surprise to see that the following C code (using -O2) is 8 time faster
than the OCaml one (even with -unsafe).

Anybody have an idea to optimize my OCaml code or know why is there this
"big" difference ?

many thanks

Frédéric Gava

ps: the difference is the same even when I use a non-polymorphic
multiplication

*************** C Code ******************

typedef struct Complexe {double re; double img;} complexe;

__inline__ void* add(void *x, void *y)
{
complexe a=((complexe*) x);
complexe b=((complexe*) y);
return (void *) &((complexe){a.re+b.re, a.img+b.img});
}

__inline__ void* mult(void *x, void *y)
{
complexe a=((complexe*) x);
complexe b=((complexe*) y);
return (void *) &((complexe){a.re*b.re-a.img*b.img, a.re*b.img+a.img*b.re});
}

void multiply_generic(int n, void *zero, void* (*add)(void *a, void *b),
void* (*mult)(void *a, void *b), void *a, void *b, void *c)
{
register int i,j,k;
void *tmp;
void *cc;

for (i=0; i<n; i++)
{
for (j=0;j<n;j++)
{
tmp=zero;
cc=c+(i*n+j);
cc=tmp;
}
}
}

complexe random_complexe(void)
{
complexe c = {c.re=(double)(rand()%1000), c.img=(double)(rand()%1000)};
return c;
}

Call this with (and with random matrix).

complexe zero;
zero.re=(double)0;
zero.img=(double)0;

********************* OCaml Code *******************************

let random_complex n = {Complex.re=Random.float n; im=Random.float n}

let matrixRandom n = Array.init (n*n) (fun _ -> random_complex 1000.)

let multiplication n zero add mul a b c =
let tmp=ref zero in
for i=0 to n-1 do
for j=0 to n-1 do
tmp:=zero;
for k=0 to n-1 do
tmp:= add (!tmp) (mul b.(k*n+j) a.(i*n+k))
done;
c.(i*n+j)<-(!tmp)
done
done

```