概述编辑本段回目录
以概率和统计的理论、方法为基础的一种计算方法,将所求解的问题同一定的概率模型相联系,用电子计算机实现统计模拟或抽样,以获得问题的近似解,故又称统计模拟法或统计试验法。
蒙特卡罗是摩纳哥的一个城市,以赌博闻名于世界。蒙特卡罗法借用这一城市的名称是为了象征性地表明该方法的概率统计的特点。
蒙特卡罗法作为一种计算方法,是由S.M.乌拉姆和J.冯·诺伊曼在20世纪40年代中叶为研制核武器的需要而首先提出来的。在此之前,该方法的基本思想实际上早已被统计学家所采用了。例如,早在17世纪,人们就知道了依频数来决定概率的方法。
20世纪40年代中叶,出现了电子计算机,使得用数学方法模拟大量的试验成为可能。另外,随着科学技术的不断发展,出现了越来越多的复杂而困难的问题,用通常的解析方法或数值方法都很难加以解决。蒙特卡罗法就是在这些情况下,作为一种可行的而且是不可缺少的计算方法被提出和迅速发展起来的。
应用领域编辑本段回目录
蒙特卡罗方法有很强的适应性,问题的几何形状的复杂性对它的影响不大。该方法的收敛性是指概率意义下的收敛,因此问题维数的增加不会影响它的收敛速度,而且存贮单元也很省,这些是用该方法处理大型复杂问题时的优势。因此,随着电子计算机的发展和科学技术问题的日趋复杂,蒙特卡罗方法的应用也越来越广泛。它不仅较好地解决了多重积分计算、微分方程求解、积分方程求解、特征值计算和非线性方程组求解等高难度和复杂的数学计算问题,而且在统计物理、核物理、真空技术、系统科学 、信息科学 、公用事业、地质、医学,可靠性及计算机科学等广泛的领域都得到成功的应用。
应用程序---Monte Carlo法( 蒙特卡罗)积分编辑本段回目录
#include<time.h>
#include<math.h>
#include<stdlib.h>
using namespace std ;
double FuncOne(double x);
/*被积一元函数*/
double MonteCarloOne(double(*p)(double),double a,double b,double L,double H);
/*一元函数积分调用接口*/
double FuncTwo(double x,double y);
/*被积二元函数*/
double Func1(double x);
/* 积分下限y=Func1(x) */
double Func2(double x);
/* 积分上限y=Func2(x) */
double MonteCarloTwo(double(*func1)(double),double(*func2)(double),double(*p)(double,double),double a,double b,double l,double h,double L,double H);
/*二元函数积分调用接口.在a
double FuncThree(double x,double y,double z);
/*被积分三元函数*/
double Func3(double x,double y);
/* 积分下限z=Func3(x,y) */
double Func4(double x,double y);
/* 积分上限z=Func4(x,y) */
double MonteCarloThree(double(*func1)(double),double(*func2)(double),double(*func3)(double,double),double(*func4)(double,double),double(*p)(double,double,double),double a,double b,double c,double d,double l,double h,double L,double H);
/*三元函数积分调用接口.在a
double FuncOne(double x)
{
return sin(x);
}
double MonteCarloOne(double(*p)(double),double a,double b,double L,double H)
{
unsigned N=100000000 ;
/* 1e8 */
double x,y ;
/*随机点*/
double t ;
/*换元*/
unsigned M=0 ;
/*处于阴影部分内的点的个数*/
srand((unsigned)time(NULL));
/*产生随机数种子*/
for(unsigned int i=0;i
{
t=double(rand())/RAND_MAX ;
/* 产生0~1的随机浮点数 */
y=double(rand())/RAND_MAX ;
x=a+(b-a)*t ;
/* 得到a~b的随机浮点数 */
/*随机点处于阴影部分内*/
if(y*(H-L)<(p(x)-L))M++;
}
double result=(H-L)*(b-a)*M/N+(b-a)*L ;
/* double型乘以unsigned型为double型 */
return result ;
/*积分结果*/
}
double FuncTwo(double x,double y)
{
return x+y ;
}
double Func1(double x)
{
return x ;
}
double Func2(double x)
{
return x*x ;
}
double MonteCarloTwo(double(*func1)(double),double(*func2)(double),double(*p)(double,double),double a,double b,double l,double h,double L,double H)
{
unsigned N=100000000 ;
/* 1e8 */
double x,y,z ;
/*随机点*/
double t,s ;
/*换元*/
unsigned m=0,M=0 ;
/*处于阴影部分内的点的个数*/
srand((unsigned)time(NULL));
/*产生随机数种子*/
for(unsigned int i=0;i
{
t=double(rand())/RAND_MAX ;
/* 产生0~1的随机浮点数 */
s=double(rand())/RAND_MAX ;
x=a+(b-a)*t ;
/* 得到a~b的随机浮点数 */
y=l+(h-l)*s ;
/* 得到l~h的随机浮点数 */
if(y>func1(x)&&y
{
m++;
/* 注意,m是与一元积分相对应的 */
z=double(rand())/RAND_MAX ;
/*随机点处于阴影部分内*/
if(z*(H-L)<(p(x,y)-L))M++;
}
}
double result=(H-L)*(b-a)*(h-l)*M/N+(b-a)*(h-l)*L*m/N ;
/* double型乘以unsigned型为double型 */
return result ;
}
double FuncThree(double x,double y,double z)
{
return x+y+z ;
}
double Func3(double x,double y)
{
return 0 ;
}
double Func4(double x,double y)
{
return FuncTwo(x,y);
}
double MonteCarloThree(double(*func1)(double),double(*func2)(double),double(*func3)(double,double),double(*func4)(double,double),double(*p)(double,double,double),double a,double b,double c,double d,double l,double h,double L,double H)
{
unsigned N=100000000 ;
/* 1e8 */
double x,y,z,u ;
/*随机点*/
double t,s,k ;
/*换元*/
unsigned m=0,M=0 ;
/*处于阴影部分内的点的个数*/
srand((unsigned)time(NULL));
/*产生随机数种子*/
for(unsigned int i=0;i
{
t=double(rand())/RAND_MAX ;
/* 产生0~1的随机浮点数 */
s=double(rand())/RAND_MAX ;
k=double(rand())/RAND_MAX ;
x=a+(b-a)*t ;
/* 得到a~b的随机浮点数 */
y=c+(d-c)*s ;
/* 得到c~d的随机浮点数 */
if(y>func1(x)&&y
{
k=double(rand())/RAND_MAX ;
z=l+(h-l)*k ;
/* 得到l~h的随机浮点数 */
/*随机点处于阴影部分内*/
if(z>func3(x,y)&&z
{
m++;
/* 注意,m是与二元积分对应的 */
u=double(rand())/RAND_MAX ;
if(u*(H-L)<(p(x,y,z)-L))
M++;
}
}
}
double result=(H-L)*(b-a)*(d-c)*(h-l)*M/N+(b-a)*(d-c)*(h-l)*L*m/N ;
/* double型乘以unsigned型为double型 */
return result ;
}
int main()
{
double a,b ;
/* x积分限 */
double l,h ;
/* y的值域或y积分限 */
double result ;
/*积分结果*/
a=0,b=3.1415926 ;
l=0,h=1 ;
result=MonteCarloOne(FuncOne,a,b,l,h);
cout<
a=1.0,b=2.0 ;
l=1.0,h=4.0 ;
double L=2.0,H=6.0 ;
/* z的值域 */
result=MonteCarloTwo(Func1,Func2,FuncTwo,a,b,l,h,L,H);
cout<
a=1.0,b=2.0 ;
double c=1.0,d=4.0 ;
l=0.0,h=6.0 ;
L=2.0,H=12.0 ;
result=MonteCarloThree(Func1,Func2,Func3,Func4,FuncThree,a,b,c,d,l,h,L,H);
cout<
return 0 ;
}
2.00011
3.35045
20.994
准确值分别为:
3.35
20.9964286