正在加载图片...
7.8 Adaptive and Recursive Monte Carlo Methods 327 free_vector(fminr,1,ndim); free_vector(fminl,1,ndim); free_vector(fmaxr,1,ndim); free_vector(fmaxl,1,ndim); if (!jb)jb=1+(ndim*iran)/175000; MNPT may be too small. rgl-regn[jb]; Apportion the remaining points between rgm=rmid[jb]; left and right. rgr=regn [ndim+jb]; fracl=fabs((rgm-rgl)/(rgr-rgl)); nptl=(unsigned long)(MNPT+(npts-npre-2*MNPT)*fracl*siglb /(fracl*siglb+(1.0-fracl)*sigrb)); Equation (7.8.23). nptr=npts-npre-nptl; regn_temp=vector(1,2*ndim); Now allocate and integrate the two sub- for (j=1;j<mndim;j++){ regions. regn_temp[j]=regn[j]; 鱼 nted for 19881992 regn_temp[ndim+j]=regn[ndim+j]; 1.00 regn_temp [ndim+jb]=rmid[jb]; miser(func,regn_temp,ndim,nptl,dith,&avel,&varl); regn_temp[jb]=rmid[jb]; Dispatch recursive call;will return back Cambridge from NUMERICAL RECIPESI regn_temp[ndim+jb]=regn[ndim+jb]; here eventually. miser(func,regn_temp,ndim,nptr,dith,ave,var); free_vector(regn_temp,1,2*ndim); *ave=fracl*avel+(1-fracl)*(*ave); server computer, (Nor *var=fracl*fracl*varl+(1-fracl)*(1-fracl)*(*var); Combine left and right regions by equation (7.8.11)(1st line) America UnN电.t THE free_vector(rmid,1,ndim); one paper ART free_vector(pt,1,ndim); 是 strictly prohibited Programs The miser routine calls a short function ranpt to get a random point within a specified d-dimensional region.The following version of ranpt makes consecutive calls to a uniform to dir random number generator and does the obvious scaling.One can easily modify ranpt to generate its points via the quasi-random routine sobseq(87.7).We find that miser with sobseq can be considerably more accurate than miser with uniform random deviates.Since the use of RSS and the use of quasi-random numbers are completely separable,however,we OF SCIENTIFIC COMPUTING (ISBN 0-521 have not made the code given here dependent on sobseg.A similar remark might be made regarding importance sampling,which could in principle be combined with RSS.(One could v@cam in principle combine vegas and miser,although the programming would be intricate.) .Further reproduction, 1988-1992 by Numerical Recipes extern long idum; -43108 void ranpt(float pt[],float regn[],int n) Returns a uniformly random point pt in an n-dimensional rectangular region.Used by miser; (outside calls ran1 for uniform deviates.Your main program should initialize the global variable idum Software. to a negative seed integer. float ran1(long *idum); Ame int ji for (j=1;j<=n;j++) pt [j]=regn [j]+(regn[n+j]-regn[j])*ran1(&idum);7.8 Adaptive and Recursive Monte Carlo Methods 327 Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copyin Copyright (C) 1988-1992 by Cambridge University Press. Programs Copyright (C) 1988-1992 by Numerical Recipes Software. Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5) g of machine￾readable files (including this one) to any server computer, is strictly prohibited. To order Numerical Recipes books or CDROMs, visit website http://www.nr.com or call 1-800-872-7423 (North America only), or send email to directcustserv@cambridge.org (outside North America). } free_vector(fminr,1,ndim); free_vector(fminl,1,ndim); free_vector(fmaxr,1,ndim); free_vector(fmaxl,1,ndim); if (!jb) jb=1+(ndim*iran)/175000; MNPT may be too small. rgl=regn[jb]; Apportion the remaining points between rgm=rmid[jb]; left and right. rgr=regn[ndim+jb]; fracl=fabs((rgm-rgl)/(rgr-rgl)); nptl=(unsigned long)(MNPT+(npts-npre-2*MNPT)*fracl*siglb /(fracl*siglb+(1.0-fracl)*sigrb)); Equation (7.8.23). nptr=npts-npre-nptl; regn_temp=vector(1,2*ndim); Now allocate and integrate the two sub￾for (j=1;j<=ndim;j++) { regions. regn_temp[j]=regn[j]; regn_temp[ndim+j]=regn[ndim+j]; } regn_temp[ndim+jb]=rmid[jb]; miser(func,regn_temp,ndim,nptl,dith,&avel,&varl); regn_temp[jb]=rmid[jb]; Dispatch recursive call; will return back regn_temp[ndim+jb]=regn[ndim+jb]; here eventually. miser(func,regn_temp,ndim,nptr,dith,ave,var); free_vector(regn_temp,1,2*ndim); *ave=fracl*avel+(1-fracl)*(*ave); *var=fracl*fracl*varl+(1-fracl)*(1-fracl)*(*var); Combine left and right regions by equation (7.8.11) (1st line). free_vector(rmid,1,ndim); } free_vector(pt,1,ndim); } The miser routine calls a short function ranpt to get a random point within a specified d-dimensional region. The following version of ranpt makes consecutive calls to a uniform random number generator and does the obvious scaling. One can easily modify ranpt to generate its points via the quasi-random routine sobseq (§7.7). We find that miser with sobseq can be considerably more accurate than miser with uniform random deviates. Since the use of RSS and the use of quasi-random numbers are completely separable, however, we have not made the code given here dependent on sobseq. A similar remark might be made regarding importance sampling, which could in principle be combined with RSS. (One could in principle combine vegas and miser, although the programming would be intricate.) extern long idum; void ranpt(float pt[], float regn[], int n) Returns a uniformly random point pt in an n-dimensional rectangular region. Used by miser; calls ran1 for uniform deviates. Your main program should initialize the global variable idum to a negative seed integer. { float ran1(long *idum); int j; for (j=1;j<=n;j++) pt[j]=regn[j]+(regn[n+j]-regn[j])*ran1(&idum); }
<<向上翻页向下翻页>>
©2008-现在 cucdc.com 高等教育资讯网 版权所有