Thread: MonteCarlo neutron transport View Single Post
 P: 1 Hello friends, I am doing simple mc simulation for neutron transport. The principle of random simulation is L=-(1/C)*ln(u), where u is a random number from the uniform distribution [0,1), and c is the total cross section (absorbed + scattered). Could you explain me what is the meaning of u, and how can I use it in my simulation?. I wrote a simple code in c, which I have to parallel it in MPI: #include "mpi.h" #include #include #include #include #include #include #define PI 3.14159265 int main (int argc, char *argv[]) { /* Simple "srand()" seed: just use "time()" */ unsigned int iseed = (unsigned int)time(NULL); srand (iseed); /*c = mean distance between neutron/atom interaction*/ /*cs = scattering component of c*/ /*cc = absorbing component of c*/ float c,cs,cc; /*mean distance between neutron/atom interaction*/ double h,l; /*thickness of the slab and length travelled*/ double u; /*uniform random number*/ double dd,d; /*direction of the neutron [0,pi]*/ double x; /*position particle in the plate*/ int n=10; /*number of samples*/ int bounce; /*boolean value*/ int i,r,b,t; r=0; b=0; t=0; /*reflected,absorbed and transmited neutrons*/ /*generates randon values between 0 and 1*/ /*dd = (float) rand()/RAND_MAX;*/ for (i = 1; i <= n; i++) { x = 0.0; d = 0.0; dd = 0.0; c = 5.5545; cc = 0.0035; cs = 5.551; h = 10.00; bounce = 1; do { dd = (float) rand()/RAND_MAX; d=PI*dd; /*generates randon values between 0 and 1*/ u = (float)rand()/RAND_MAX; l = -(1/c)*log(u); /*lenght travelled by the particle*/ x = x + l*cos(d); //printf("value of x=%2.4f\n",x); if (x < 0.0) /*reflected*/ { r = r + 1; bounce = 0; // printf("reflec r= %d\n",r); } else { if (x >= h) /*transmited*/ { t = t+1; bounce = 0; // printf("trans t=%d\n",t); } else{ if (u < cc/c) { b = b+1; bounce = 0; // printf("absorbed b=%d\n",b); } else d = u * PI; } }/*end first else*/ } while (bounce = 1); } /*end for*/ printf("Particles reflected = %d\n",dd); printf("Particles absorbed = %d\n",b); printf("Particles transmittedd = %d\n",b); }/*end main*/ Please any help is more than welcome