棄却サンプリングのプログラム
#include <oxstd.h> #include <oxfloat.h> #include <oxprob.h> main(){ decl x,y,c,n,i; decl n_accept,n_count; decl file; n=1000; y=zeros(n,1); n_accept=n_count=0; for(i=0;i<n;++i){ do{ n_count+=1; x=ranchi(1,1,1); c=sqrt(M_2PI)*3^(1.5)*exp(-1.5); } while(ranu(1,1) > x*exp(-x)/(c*1.0/sqrt(M_2PI*x)*exp(-0.5*x))); n_accept+=1;y[i]=x; } // Output format(600); println("Acceptance Probability= ",n_accept/n_count); file = fopen("gamma.dat", "w"); fprint(file,"%10.5f",y);fclose(file); }