#define pi 3.14159265358979323846264338327950288419716939937510
#include <stdio.h>
#include <stdlib.h>
#include <math.h>


void quasi_random_vector(int dim,double* a,double* b, double* x);
void pseudo_random_vector(int dim,double* a,double* b, double* x);
int plainmc(
	int d, double a[], double b[], double (*func)(double*),
	void (*rndx)(int,double*,double*,double*),
	int N, double* result, double* error);


double fun(double x[]){return 1/(1-cos(x[0])*cos(x[1])*cos(x[2]))/pi/pi/pi;}

int main(){

int N; double exact=1.3932039296856768591842462603255;
int d=3; double a[] ={0,0,0}; double b[]={pi,pi,pi};

double result,error;

N=10000; plainmc(d,a,b,&fun,pseudo_random_vector,N,&result,&error);
printf("-------------- pseudo-random: -------------\n");
printf("N\t= %d\n",N);
printf("result\t= %g\n",result);
printf("exact\t= %g\n",exact);
printf("estimated error\t= %g\n",error);
printf("actual error\t= %g\n",fabs(result-exact));

N=100000; plainmc(d,a,b,&fun,quasi_random_vector,N,&result,&error);
printf("-------------- quasi-random: -------------\n");
printf("N\t= %d\n",N);
printf("result\t= %g\n",result);
printf("exact\t= %g\n",exact);
printf("estimated error\t= %g\n",error);
printf("actual error\t= %g\n",fabs(result-exact));

return 0;
}
