随机微分方程数值解在泄洪风险分析中的应用_徐敏.pdf

《随机微分方程数值解在泄洪风险分析中的应用_徐敏.pdf》由会员分享，可在线阅读，更多相关《随机微分方程数值解在泄洪风险分析中的应用_徐敏.pdf（5页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、第36 卷第9 期2006 年9 月数学的实践与认识MATHEMATICS IN PRACTICE AND THEORYVol.36No.9Sep.,2006随机微分方程数值解在泄洪风险分析中的应用徐敏2,胡良剑3,丁永生1,胡盈3,周林峰1(1.东华大学信息科学与技术学院,上海201620)(2.东华大学旭日工商管理学院,上海200051)(3.东华大学应用数学系,上海201620)摘要:根据泄洪过程中库水位过程的随机微分方程,利用数值解方法,模拟了随机干扰下的库水位及其波动状况.采用相应公式计算了洪水漫越坝顶事件的概率以及库水位过程在不同时刻的样本均值.并通过比较在同样强度的随机干扰下库水

2、位的高低状况,确定出各种泄洪方案的优劣,从而对防洪工作具有重要的指导意义.关键词:随机微分方程;数值解;欧拉法;泄洪风险1引言收稿日期:2005-06-27基金项目:国家自然科学基金(60474037);教育部新世纪优秀人才支持计划(NCET-04-415)对于洪水,风暴潮等自然灾害事件,风险分析是一种极为有效的工具 1.由于洪水过程具有很多种不确定性因素,随机性便很自然地被引入到防洪过程的分析.近年来,这方面的很多研究工作都认为洪水过程是一随机点过程 24;Sen 以一阶马尔科夫过程为工具对具有线性相关结构的水文系列风险进行计算 5.特别地,随机微分方程被引入防洪风险分析,由此建立了水库调洪

3、演算的随机数学模型 6,7.由于随机微分方程本身的复杂性,除了一些线性的或者特殊结构的方程以外,可求出显示解的随机微分方程很少 8,9.本文中讨论的随机微分方程不具有上述性质,因此无法求出显示解.姜树海根据其解过程的一阶概率密度函数满足Fokker-Plank 向前方程,而这一方程又是一偏微分方程,从而利用偏微分方程的有限差分法求出其数值解 6,但这种方法不能求得概率特征,于是JC 计算方法被用于近似地算出洪水漫越坝顶的概率 7.不难看出,这种方法由于采用多次转化,误差比较大.本文利用随机微分方程数值解方法,结合实际例子,分析总结了库水位在布朗运动干扰下的随机波动状况;直接求出了洪水漫坝的风险

4、概率和库水位过程在不同时刻的数学期望.并且还对不同的方案进行分析比较,以确定哪种方案的效果更好,从而可对防洪决策过程提供一定的依据.2调洪过程的随机微分方程调洪过程中入库洪水和出库泄量是随机过程,其库容水位满足随机微分方程 6:dH(t)=Q-(t)-q-(H,c)G(H)dt+dB(t)G(H)H(t0)=H0(1)其中 H(t)为库水位过程;H0为初始库水位,它是一个随机变量;Q(t)为任意时刻入库洪水量;q(h,c)为相应时刻的泄洪流量;Q-,q-分别为来流和泄洪的均值过程线;c为流量系数等水利参数.G(H)=dW(H)dH,W(H)是水库的库容量,B(t)是一均值为零的Wiener 过

5、程,dB(t)/dt 是一正态白噪声,B(t)的一维概率密度函数 f(B)为:f(B)=12P t Rexp-B22R2t.由上式可以看出,E B(t)=0,D B(t)=R2t.洪水漫越坝顶的泄洪风险率定义为 Pf=Pf H E Z,其中,Z 为相应的坝高.3计算方法由于随机微分方程很少可求出显示解,故其数值解方法得到广泛的研究和应用.相对于常微分方程数值法而言,随机微分方程数值解方法引入了随机增量,它将所考虑的时间区间做有限划分,一步一步地在节点处生成样本轨道的逼近值,其数值解方法主要有:Eu-ler 法、Milstein 法、Runge-Kutta 法等.这里采用Euler 法.3.1随

6、机微分方程解的欧拉逼近法考虑一般随机微分方程:dXt=a(t,Xt)dt+b(t,Xt)dWt(2)其中,t0F t F T,初始条件是 Xt0=X0.我们对时间区间 t0,T 进行离散化:t0=S0 S1 Sn 0)阶强收敛于一个Itod过程 X,如果存在一不依赖于 D的正常数 C 及一 D0(0),使得下式成立:E(D)=E(XT-YD(T)F CDC,PD(0,D0)实际上,就Euler 法而言,如果方程(2)的两个系数 a(t,Xt)和b(t,Xt)满足Lipschitz 条件和线性增长条件,则它为 C=0.5 阶强收敛 8.3.2计算公式对调洪过程中的随机微分方程(1),运用Eule

7、r 法,有下面的迭代公式:H(j,k)=H(j-1,k)+Q(j-1)-q(j-1)G(H(j-1,k)(tj-tj-1)+1G(H(j-1,k)B(j)-B(j-1),j=1,2,M;k=1,2,K(3)154数学的实践与认识36 卷其中,j 和k 分别代表时间节点和轨道,M 是时间节点的个数,K 是轨道数.对每一条固定轨道 k,根据以上公式及已知条件,可求出不同时间节点上的库水位值 H(1,k),H(2,k),H(M,k).对每一条轨道 k 及每一个时间节点 j,定义函数 f(j,k):f(j,k)=1,当 H(j,k)E Z;0,当 H(j,k)Z.其中,Z 表示坝高,即对同一条轨道下的

8、每一个时间节点处的库水位是否超过坝高进行判断,若超过则定义为1,否则定义为0.接下来,对每一条轨道 k,再定义函数 g(k):g(k)=max1F jF Mf(j,k)即只要在一个时间节点处出现了洪水漫越坝顶,那么就算漫越坝顶.洪水漫坝的概率为:pf=Kk=1g(k)K另外,库水位过程 H(t)的均值为:H-(j)=1KMk=1H(j,k),1 F j F M.4仿真结果以湖北某水库的调洪实际 10作为仿真实例,在(3)中取M=33,K=100,运用随机微分方程数值解方法.通过实际测量,给出该水库的w-h 关系曲线,并拟合经验关系式:w=1.4594+0.0475h+0.0175h2 108(

9、m3)该水库大坝按千年一遇洪水标准设计,坝顶高程 Z=19.5m(以堰顶高程为基准),调洪方案共有三种 7:方案一:水库的泄洪建筑物有主坝5 孔,副坝1 孔溢流堰,每孔宽 9.4m,堰顶高程为 0,流量系数m-1=0.49;另设一底孔,宽4m,高6m,孔中心线高程为-5.5m,流量系数m-2=0.90,该库的防洪限制水位 h0=12.5m.防洪调度规则为:当入库流量Q 3000m3/s时不泄流,q=0;当 3000m3/s Q 6300m3/s 时,泄洪设施全部敞泄.方案二:方案一的泄洪建筑物规模和调洪规则保持不变,仅以主坝5 孔溢流堰泄洪.方案三:方案一的泄洪设施规模保持不变,改变调洪规则,

10、将敞泄的控制限值 Q 降低至5500m3/s,敞泄的时间提前到 t=27h.表1洪水漫坝的风险概率过程强度R方案方案一方案二方案三310401091040.040.8701.51050.210.790.02采 用方案一,当Wiener 过程的 R 取 3104时,对每一条固定轨道,运用以上计算方法,可得出对应的库水位轨线,图1中以第 100 条轨道为例.相应地可得出100 条库水位线的随机分布,如图 2 所示.由此可计算出洪水漫坝的风险概率,如表1所示.同时为了对三种方案进行比较,以 k=100 为例,得到其对应的1559 期徐敏,等:随机微分方程数值解在泄洪风险分析中的应用库水位线,如图3

11、所示.从图3 中可以看出,在同样强度的随机干扰下,在泄洪开始的一段时间,三种方案的库水位大体相同.随着时间的推移,方案一和方案二的库水位仍维持基本相同,但方案三的库水位明显下降;此后,方案一的库水位高度低于方案二的高度,方案三的库水位高度低于方案一的高度.同时,在同样强度(即Wiener 过程的 R 取相同的值)的随机干扰下,可得出三种方案的洪水漫坝的风险概率,如表 1 所示.由此可知三种方案的优劣依次为:方案三优于方案一,方案一又优于方案二.图1库水位轨线(方案一,R=3104,k=100)图2 库水位随机分布(方案一,取R=3104,共100 条轨道)图 3三种方案库水位轨线(R=3104

12、,k=100)图 4方案一的库水位均值(R=3104)在表1 中,值得注意的是方案二,随着随机干扰强度的增大,相应的风险概率反而减小,这是由于随机干扰强度的增大触发了泄洪机制,使得当水位到达一定高度时泄洪设施全部敞泄,水位回落很快,导致了漫坝的概率有所降低.进一步地,还可以求出库水位的均值,图4为方案一在 R=3 104时的库水位均值.5结论本文利用随机微分方程的数值解方法对调洪过程中的库水位进行分析,得出其在不同强度的布朗运动干扰下的随机波动状况,并可求出相应的洪水漫坝的风险概率及库水位的样本均值,本文采用的方法更直观更直接,且精确度较高.156数学的实践与认识36 卷参考文献:1 王栋,潘

13、少明,吴吉春,朱庆平.洪水的风险分析 J.应用基础与工程科学学报,2004,(增刊):134140.2 徐宗学,邓永录.洪水风险率HSPPB模型及其应用 J.水利水电学报,1989,(1):4655.3 徐宗学,曾光明.洪水频率分析HSPPC 模型应用研究 J.水科学进展,1992,3(3):175180.4 Xu Z X,Li J Y,Ito K.Clustering stochastic point process model for flood risk analysisJ.Stochastic Hydrologyand Hydraulics,1998,12(1):5364.5 Z Se

14、n.Simple risk calculations in dependent hydrological series J.Hydrological Science,1999,44(6):871878.6 姜树海.水库调洪演算的随机数学模型 J.水科学进展,1993,4(4):294300.7 姜树海.随机微分方程在泄洪风险分析中的应用J.水利学报,1994,(3):19.8 Kloeden P E,Platen E.N umerical Solution of Stochastic Differential Equations M .Spring-Verlag,1992.9 Kouritzi

15、n M A,Deli L.On explicit solution to stochastic differential equations J.Stochastic Analysis and Appli-cations,2000,18(4):571580.10 水电规划设计院,长江流域办公室主编.水利动能设计手册(防洪分册)M.水利电力出版社,1988.Numerical Solution of Stochastic DifferentialEquations with Application to RiskAnalysis of Releasing FloodwaterXU Min2,HU

16、 Liang-jian3,DING Yong-sheng1HU Ying3,ZHOU Lin-feng1(1.College of Information Sciences and Technology,Donghua University,Shanghai 201620,China)(2.Glorious Sun School of Business and Management,Donghua University,Shanghai 200051,China)(3.Department of Applied Marthematics,Donghua University,Shanghai

17、201620,China)Abstract:Based on the process of the water level of a reservoir satisfying stochastic differen-tial equations during the process of releasing floodwater,this paper simulates the water leveland its fluctuation under stochastic disturbance of Wiener Process by means of numerical solu-tion

18、 of stochastic differential equations.It also calculates the mean value at different time of theprocess of water level as well as the probability of overflowing with corresponding formula.Furthermore,it compares the water level under stochastic disturbance with the same strengthin order to ascertain which one is more effective.As a result,it is helpful for the flood practice.Keywords:stochastic differential equation;numerical solution;Euler method;risk of releas-ing floodwater1579 期徐敏,等:随机微分方程数值解在泄洪风险分析中的应用