用c语言编程蒙特卡罗法,跪求能用c语言生成蒙特卡罗方法中的sobol序列的代码...

);//read mdeg[]

for(i=1;i<=MAXDIM;i++)

fscanf(frl,"%d",&ip);//read ip[]

fclose(frl);

//set the values of iv[]//

srand((unsigned)time(NULL));//give the rand seed

for(i=1;i<=mdeg[MAXDIM];i++)

{

for(j=1;j<=MAXDIM;j++)

{temp=2*(int((1L<

iv[(i-1)*MAXDIM+j]=temp;

}

}

/*Initialize; don't return a vector*/

for(k=1;k<=MAXDIM;k++)

ix[k]=0;

in=0;

if(iv[1]!=1)

return;

fac=1.0/(1L<

for(j=1,k=0;j<=MAXBIT;j++,k+=MAXDIM)

iu[j]=&iv[k];

for(k=1;k<=MAXDIM;k++)

{

for(j=1;j<=mdeg[k];j++)

iu[j][k]<<=(MAXBIT-j);

/*Stored values only require normalization*/

for(j=mdeg[k]+1;j<=MAXBIT;j++)

{  /*Use the recurrence to get other values.*/

ipp=ip[k];

i=iu[j-mdeg[k]][k];

i ^=(i>>mdeg[k]);

for(l=mdeg[k]-1;l>=1;l--)

{

if(ipp&1)

i ^=iu[j-1][k];

ipp>>=1;

}

iu[j][k]=i;

}

}

}

else

{

/*Calculate the next vector in the sequence*/

im=in++;

for(j=1;j<=MAXBIT;j++)

{

/*Find the rightmost=ero bit*/

if(!(im&1))

break;

im>>=1;

}

if(j>MAXBIT)

printf("\nMAXBIT_too_small_in_sobseq\n");

im=(j-1)*MAXDIM;

for(k=1;k<=IMIN(*n,MAXDIM);k++)

{

ix[k] ^=iv[im+k];

x[k]=ix[k]*fac;

}

}

}

int main()

{FILE *fwl;

static double xsob[S+1];

int i,j;

int ini,run;

ini=-1;

sobseq(&ini,xsob);

if(( fwl=fopen("sobol_points.txt","w"))==NULL)

{ printf( "The_file_'sobol_points'_was_not_opened\n");

exit(0);

}

for(i=1;i<=1024;i++)

{ run=S;

sobseq(&run,xsob);

for(j=1;j<=S;j++)

fprintf(fwl,"%f\t",xsob[j]);

fprintf(fwl,"\n");

}//end of i

fclose(fwl);

}