资源描述:
《第九讲随机系统基础》由会员上传分享,免费在线阅读,更多相关内容在学术论文-天天文库。
1、第九讲随机系统与随机微分方程DesmondJ.Higham†,AnAlgorithmicIntroductiontoNumericalSimulationofStochasticDifferentialEquations,SIAMREVIEW2001,Vol.43,No.3,pp.525–546C.W.Gardiner,HandbookofStochasticMethods,Springer,1996,(中文版,上海科技出版社)9.1布朗运动十九世纪末之前,人们都认为微分方程求解解释自然现象,只要收集数据,则我们可以预测未来。但是
2、,量子力学基于纯粹的统计原理,二是混沌的出现发现微小的误差可以被迅速放大,这使得人们对于预测未来的信心开始怀疑。1827年RoberBrown研究了悬浮在液体表面的颗粒无规则运动以证实分子的存在,这种运动的数学解释是1905年Einstein给出的。假设个粒子悬浮在液体中,时间间隔远小于观测时间内,单个粒子坐标增加一个量,位移在和之间的粒子数目为,满足这里就是概率密度,满足,设为单位体积内粒子数目,那么经过时间间隔,时间通过点和之间的粒子数目为按照小参数展开得到令,上式化为这就是著名的扩散方程,偏微分方程的一种。其解为可以看出,沿
3、x轴的方差为是随着时间增加而增加的。法国物理学家Langevin提出了认为比爱因斯坦更加简单的办法描述布朗运动受力分成二部分,一為粘滯阻力-αv,一為隨机作用力(实际就是噪声),粘滯阻力仍來自介質分子對顆粒的碰撞,將顆粒看作半徑為a的小球,在粘滯系數為η的流體中運動,則有上式稱為斯托克斯公式(stokes’slaw)。將上式對大量顆粒子求平均,即把大群顆的運動方程相加然後用顆粒數,Langevin推导出这里是波尔兹曼常数,而由扩散方程(diffusionequation)已经得到,那么比较可得这是和实际相符合的。正是由于对于布朗运
4、动的研究,诞生了第一个随机微分方程经过Ornstein,Langevin,和Winer1918年严格的数学描述,此后Brown运动专指Winer过程。标准的Winer过程是在时间间隔内的连续随机变量,满足下面三个条件:A:,以概率1B:对于时刻,变量的增量-为正态分布,均值为0,方差为,即-,C对于时刻,增量-和增量-是相互独立的,又称为独立增量过程。从计算的角度,我们把Winer过程离散化,采样时间,这里整数N,我们记,那么,按照上面条件有:A:B:,C,且不同下标下相互独立例子:Matlab仿真布朗运动%BPATH1Brown
5、ianpathsimulationrandn('state',100)%setthestateofrandnT=1;N=50;dt=T/N;dW=zeros(1,N);%preallocatearrays预先分配内存W=zeros(1,N);%forefficiencydW(1)=sqrt(dt)*randn;%firstapproximationoutsidetheloop...W(1)=dW(1);%sinceW(0)=0isnotallowedforj=2:NdW(j)=sqrt(dt)*randn;%generalincr
6、ementW(j)=W(j-1)+dW(j);end>>plot([0:dt:T],[0,W],'r-')%plotWagainstt>>xlabel('t','FontSize',16)>>ylabel('W(t)','FontSize',16,'Rotation',0)------------------Listing1M-filebpath1.m.Fig.1布朗运动的模拟解释:不同种子能够得到不同的运动轨迹,由于Matlab向量起始下标为1,所以应该表示为,但是画的过程中将其标记为时间0。bpath1.m.这个程序效率不高,
7、因为没有利用Matlab向量化的特点,将其优化,我们得到程序%BPATH2Brownianpathsimulation:vectorizedrandn('state',100)%setthestateofrandnT=1;N=50;dt=T/N;dW=sqrt(dt)*randn(1,N);%incrementsW=cumsum(dW);%cumulativesumplot([0:dt:T],[0,W],'r-')%plotWagainsttxlabel('t','FontSize',16)ylabel(’W(t)’,’FontS
8、ize’,16,’Rotation’,0)Listing2M-filebpath2.m.我们不用循环产生dW,由于事先知道向量长度N,所以我们一次产生1*N个dW,都是相同方差,相同分布,cumsum命令的含义就是j个元素为前面元素累加和,正是这个