parameter(n=3,pi=3.1415927); dimension a(n),b(n),x(n)
external fun,pseudo

a(1)=0; b(1)=pi; a(2)=0; b(2)=pi; a(3)=0; b(3)=pi;
acc=0.01; eps=0.01

m=9000; call plainmc(n,a,b,x,fun,m,res,acc,eps,pseudo)
print *,'----pseudo----- N\t=',m
print *,'res=\t',res,' acc=\t',acc,' eps=\t',eps
exact=1.3932039296856768591842462603255
print *,'exact=\t',exact,' err=\t',exact-res,' relerr=\t',(exact-res)/exact

end

function fun(n,p)
parameter(pi=3.1415927); dimension p(n)

x=p(1); y=p(2); z=p(3); fun=1./(1-cos(x)*cos(y)*cos(z))/pi**3

return; end