Probability Distribution

#include<stdio.h>
#include<stdlib.h>
#include<math.h>

#define e 2.7182818284

double xy( int, double );

void media_fraz( double [][1001], double [], int );

double xy( int n, double a ){
  double x=0., y=0., C=0.;
  while(y>=C){
    x=((double)lrand48()/RAND_MAX*14)-4;
    y=(double)lrand48()/RAND_MAX;
    if(x>(-a)){
      C=pow(e, (-pow(x,2.)/2.));
    }
    else{
      C=(pow((n/a),n)*pow(e, pow(a,2.)/2.))*pow((n/a - a)-x, -n);
    }
  }
  return x;
}

void media_fraz( double x[][1001], double mf[], int n){
  double media=0, frac=0;
  int i;
  for( i=1; i<1001; i++){
    media+=x[n-1][i]/1000.;
    if(fabs(x[n-1][i])<1){
    frac++;
    }
  }
  mf[0]=media;
  mf[1]=frac/1000.;
}

main(){
  double a, x[4][1001]={0}, mf[2]={0}, frac;
  int i, n;
  srand48(time(0));
  printf("Questo programma calcola la media e la frazione degli eventi con |x| < 1 di un campione di dati distribuito secondo la distribuzione di probabilita` C(x) per vari valori del suo parametro n.\n");
  printf("Inserisci un valore di alfa strettamente maggiore di zero che sarà utile a calcolare la distribuzione. \n a = ");
  scanf("%lf", &a);
  while(a<=0){
    printf("Il valore di alfa è sbagliato! Reinseriscilo. \n a = ");
    scanf("%lf", &a);
  }
  for(n=1; n<5; n++){
    x[n-1][0]=1;
    for(i=1; i<1001; i++){
      x[n-1][i]=xy( n, a );
    }
    media_fraz( x, mf, n);
    printf("I risultati per n = %d sono : \n < x > = %lf \n frac = %lf \n", n, mf[0], mf[1]);
    if(mf[1]>0.68){
      printf("Il valore di frac = %lf per n = %d è maggiore di 0.68.\n", mf[1], n);
    }
  }
}