1、随机水文学随机水文学基础部分基础部分 讲述随机讲述随机水文学的水文学的概念、方概念、方法和发展法和发展 第第1 1章章随机模型随机模型 归纳与总结归纳与总结绪论绪论介绍随机过程介绍随机过程的基本理论、的基本理论、水文随机过程水文随机过程分析方法和随分析方法和随机模拟技术,机模拟技术,同时介绍水文同时介绍水文序列组成分析序列组成分析 第第23章章介绍随机介绍随机水文模型水文模型第第47章章介绍随机介绍随机模型的选模型的选择及应用择及应用 第第8章章平稳平稳随机随机模型模型季节性季节性随机模随机模型型尊重弱者 比比“大大”字等级更高的是字等级更高的是“老老”字,一个人先成为字,一个人先成为“大大人
2、人”才能成为才能成为“老人老人”。那么,。那么,我已经做了大半辈子的我已经做了大半辈子的“老老师师”,那就退后一步叫,那就退后一步叫“大师大师”也可以吧。也可以吧。关于设立关于设立“余秋雨大师工余秋雨大师工作室作室”时,余秋雨回应。时,余秋雨回应。3.1 水文序列及其组成水文序列及其组成 3.3 水文序列谱分析水文序列谱分析 3.2 水文序列相关分析水文序列相关分析 3.4 水文序列组成成分识别水文序列组成成分识别第三章第三章 水文序列分析方法及其随机模拟技术水文序列分析方法及其随机模拟技术 3.5 纯随机序列的随机模拟纯随机序列的随机模拟3.1 水文序列及其组成水文序列及其组成 水文现象随时
3、间水文现象随时间变化的过程,称为水变化的过程,称为水文过程。实际获得的文过程。实际获得的水文过程都是离散的,水文过程都是离散的,因此又常称为水文序因此又常称为水文序列。列。Hydrology time series(a)3000400050006000700080009000100001948 1954 1960 1966 1972 1978 1984 1990 1996 2002t/年Xt/m3/s(b)050001000015000200002500030000Oct-95Jul-96Mar-97Nov-97Jul-98Mar-99Dec-99Aug-00Apr-01Dec-01Sep-0
4、2t/月Xt/m3/s一一 水文序列水文序列 水文现象(如河川水文现象(如河川流量、河中溶解氧、水流量、河中溶解氧、水位等)随时间的变化一位等)随时间的变化一般是连续的。然而,连般是连续的。然而,连续过程有一个缺点,续过程有一个缺点,就就是不可能在数字计算机上进行处理。为了研究和计算的方便,常常将连续是不可能在数字计算机上进行处理。为了研究和计算的方便,常常将连续水文过程离散化处理得到水文时间序列(简称水文序列)。水文过程离散化处理得到水文时间序列(简称水文序列)。通常有三种方式获取水文序列:通常有三种方式获取水文序列:取时间区间上的平均值取时间区间上的平均值 日平均流量过程、月平均流量过程日
5、平均流量过程、月平均流量过程 旬平均流量过程、年平均流量过程旬平均流量过程、年平均流量过程 按某种规则选取特征值按某种规则选取特征值 年年/月月/日最大值,如洪峰流量,日最高水位;日最大值,如洪峰流量,日最高水位;年最大过程,如最大年最大过程,如最大3 3日洪量,最大日洪量,最大7 7日洪量等。日洪量等。在离散时间上取样在离散时间上取样 固定时刻取样固定时刻取样 时间轴上固定时间间隔取样时间轴上固定时间间隔取样二二 水文序列的组成水文序列的组成 图图3-1 水文序列组成图水文序列组成图ttttSNPX 从上图可看出,水文序列由多种成分组成。一般假定水文序列是多从上图可看出,水文序列由多种成分组
6、成。一般假定水文序列是多种成分的线性叠加,即种成分的线性叠加,即tPtNtS不规则振荡和随不规则振荡和随机影响造成的机影响造成的水文序列水文序列确定性成分确定性成分随机性成分随机性成分周期的周期的非周期的非周期的平稳的平稳的非平稳的非平稳的趋势趋势突变突变跳跃跳跃简单周期简单周期复合周期复合周期近似周期近似周期相依的相依的独立的独立的tX 周期成分、非周期成分、随机成分是水文序列的主要成分。三者并不周期成分、非周期成分、随机成分是水文序列的主要成分。三者并不一定同时存在。一定同时存在。当当Pt+Nt=0时,时,Xt=St,为随机成分序列;,为随机成分序列;当当St=0时,时,Xt=Pt+Nt,
7、为近似确定序列。,为近似确定序列。1、Xt=Pt+St 2、Xt=Nt+St 水文序列是一定自然条件和气候条件下的产物,当其中某种条件发生水文序列是一定自然条件和气候条件下的产物,当其中某种条件发生了变化,就会出现趋势、跳跃和突变现象。了变化,就会出现趋势、跳跃和突变现象。整体趋势随机成分整体趋势随机成分局部趋势随机成分局部趋势随机成分0100200t100150200X(t)跳跃随机成分跳跃随机成分突变随机成分突变随机成分 3、Xt=Pt+Nt+St趋势周期成分随机成分趋势周期成分随机成分如何识别周期成分、非周期成分和随机成分呢?如何识别周期成分、非周期成分和随机成分呢?非周期成分非周期成分
8、Nt随机成分随机成分St周期成分周期成分Pt谱分析技术谱分析技术统计技术统计技术成因分析成因分析相关分相关分析技术析技术 3.2 水文序列相关分析水文序列相关分析 判断水文序列是否相依判断水文序列是否相依(相关相关),相依程度如何,需要进行相关分析。,相依程度如何,需要进行相关分析。当研究水文序列自身相依性时,使用自相关分析技术;当研究两个当研究水文序列自身相依性时,使用自相关分析技术;当研究两个水水文序列之间的关系时,使用互相关分析技术。文序列之间的关系时,使用互相关分析技术。1 概念概念 描述水文序列自身内部线性相依程度的统计技术。一般用自相关系数描述水文序列自身内部线性相依程度的统计技术
9、。一般用自相关系数或自相关图刻画。或自相关图刻画。2)()(xCovTTdtutxutxTCov0)()(1lim)(TTdttxTtXEu0)(1lim)(TTxdtutxT022)(1lim其中其中为时间间隔为时间间隔,又称滞时又称滞时.一一 自相关分析(自相关分析(Auto-correlation analysis)2 理论自相关系数和理论自相关图理论自相关系数和理论自相关图 (1)对于连续平稳随机过程对于连续平稳随机过程X(t)的一个样本函数的一个样本函数x(t),其自相关系数为,其自相关系数为(2)对于离散平稳随机过程对于离散平稳随机过程Xt的一个样本函数的一个样本函数xt,其自相关
10、系数为,其自相关系数为2)()(xkCovk其中其中,k(=0,1,2,)为时间间隔为时间间隔,又称滞时。又称滞时。knikiinuxuxknkCov1)(1lim)(niinxntXEu11lim)(niinxuxn122)(1lim (3)理论自相关系数和自相关图理论自相关系数和自相关图 称称(k)、()为滞时为滞时k、的的理论自相关系数。满足理论自相关系数。满足)()(-,1)(1)()(-,1)(1kkk(k)、()随随滞时滞时k、而而变化,称这个变化过程图变化,称这个变化过程图为理论自相关图。为理论自相关图。k(k)()kttkkvCokr)()(3 由实测资料估计自相关系数和自相关
11、图由实测资料估计自相关系数和自相关图 设有实测资料设有实测资料xt(t=1,2,n)。样本自相关系数为。样本自相关系数为其中其中k(=0,1,2,)为滞时。为滞时。kntktktttxxxxknkvCo1)(1)(kntttxknx11knttttxxkn122)(1样本协方差样本协方差tt+k时刻样本标准差时刻样本标准差kntktktktxxkn122)(1kntktktxknx11样本自相关系数样本自相关系数nmmkxxxxxxrnttkntkttk;,2,1,0 ,)()(121当当n很大且很大且k较小时,样本自相关系数计算公式可简化为较小时,样本自相关系数计算公式可简化为其中其中m(=
12、0,1,2,)为最大滞时。当为最大滞时。当n50时,时,mn/4;当当ni)中中xiU /2,故该站故该站19832000年年径年年径流有趋势成分。流有趋势成分。设水文序列由线性趋势成分设水文序列由线性趋势成分Tt和随机成分和随机成分 t组成,即组成,即tttTxbtaTt式中式中 。线性趋势识别步骤:。线性趋势识别步骤:(1)回归分析确定参数回归分析确定参数a和和bt bxattxxttbntntt)(/)(121 (2)回归效果显著性检验回归效果显著性检验nttnt11nttxnx11其中,其中,构造统计量:构造统计量:bSbTntbttSS1222)(其中,其中,2)()(112222n
13、ttbxxSntntt统计量统计量T服从自由度为服从自由度为(n-2)的的t分布。分布。假设原序列无线性趋假设原序列无线性趋势。给定显著水平势。给定显著水平 后,后,查算查算T /2。当当 时,时,接受原假设,即线性趋势接受原假设,即线性趋势不显著;反之,线性趋势不显著;反之,线性趋势显著。显著。2/TT 3 线性趋势回归检验法线性趋势回归检验法4 非线性趋势回归检验法非线性趋势回归检验法mmttbtbtbaT221 由最小二乘法(由最小二乘法(Least square metod)估计回归参数估计回归参数b1,b2,bm。趋势趋势成分是否显著,必须对回归系数成分是否显著,必须对回归系数b1,
14、b2,bm和回归方程进行假设检验。和回归方程进行假设检验。T(t)=0.386t2+4.923t+968.9R2=0.99396097098099010001010102010301040105010601070024681012t/年x(t)/万人成都市人口与成都市人口与t的拟合曲线的拟合曲线 设水文序列由非线性趋势成分设水文序列由非线性趋势成分Tt和随机成分和随机成分 t组成,其中组成,其中二二 跳跃成分识别跳跃成分识别 跳跃(跳跃(shift or jump)指水文序列从一种状态过渡到另)指水文序列从一种状态过渡到另一种状态表现出来的急剧变化形式。下图给出跳跃成分示意一种状态表现出来的急
15、剧变化形式。下图给出跳跃成分示意图,其中图,其中 为突变点。例如一个流域若突发大面积的森林火为突变点。例如一个流域若突发大面积的森林火灾,则径流会突然变化,形成跳跃成分。灾,则径流会突然变化,形成跳跃成分。跳跃一般表现在均跳跃一般表现在均值、方差、自相关系数等统计特性上值、方差、自相关系数等统计特性上。统计分析国内外大量已建水库资料发现,河流拦蓄后,多数都引起了统计分析国内外大量已建水库资料发现,河流拦蓄后,多数都引起了库区上游流域降雨量和径流量的改变,绝大部分表现为降雨量和径流量减库区上游流域降雨量和径流量的改变,绝大部分表现为降雨量和径流量减少,且减少程度与河流拦蓄规模有关,拦蓄规模(即水
16、库总库容与年径流少,且减少程度与河流拦蓄规模有关,拦蓄规模(即水库总库容与年径流量的比值)越大,减少的程度也越大。但水库下游和其它不受水库影响的量的比值)越大,减少的程度也越大。但水库下游和其它不受水库影响的干支流径流量并未改变,主要表现为流域降水量与径流量时空分布发生了干支流径流量并未改变,主要表现为流域降水量与径流量时空分布发生了一定变化。如丹江口水库修建后,年均径流量减少了一定变化。如丹江口水库修建后,年均径流量减少了10.3%(统计(统计1935 1990 年资料),尤其自年资料),尤其自1990 年以来,径流量减少幅度加大,年以来,径流量减少幅度加大,19902004 的的14 年径
17、流量平均值只有年径流量平均值只有287.2 亿亿m3,比建库前的,比建库前的387.6 亿亿m3减少了减少了25.9%(见图见图2)水文序列中的跳跃是人为的或自然的原因引起的。水文序列中的跳跃是人为的或自然的原因引起的。突变是跳跃的一种特殊形式,是瞬间的行为,如下图所突变是跳跃的一种特殊形式,是瞬间的行为,如下图所示。突变发生后,水文序列又保持原来的特性,如溃坝、泥示。突变发生后,水文序列又保持原来的特性,如溃坝、泥石流导致河道堵塞等,这将引起流量的突变,但随着临时水石流导致河道堵塞等,这将引起流量的突变,但随着临时水坝的冲毁,又恢复到原来状态。坝的冲毁,又恢复到原来状态。(3)统计推断统计推
18、断(有序聚类分析法有序聚类分析法)1 先确定分割点先确定分割点 2 显著性检验显著性检验(1)成因分析成因分析(2)时序累计值相关曲线法时序累计值相关曲线法(3)均值、方差齐次性检验均值、方差齐次性检验(1)游程检验法游程检验法(2)秩和检验法秩和检验法非参数检验法非参数检验法分割点分割点/干扰点干扰点 流域系统一旦受到干扰(例如流域系统一旦受到干扰(例如人类活动的影响),则其平稳性遭人类活动的影响),则其平稳性遭受破坏,而演变为非平稳性序列。受破坏,而演变为非平稳性序列。对于受到干扰影响的水文序列,推对于受到干扰影响的水文序列,推估出干扰的时刻,即干扰点(估出干扰的时刻,即干扰点(突变突变点
19、点),在水文分析计算中具有实际),在水文分析计算中具有实际的意义。的意义。成因分析成因分析 水文现象受到干扰后,其时序水文现象受到干扰后,其时序变化错综复杂,它是固有的随机变变化错综复杂,它是固有的随机变化和干扰引起的变化综合而得的结化和干扰引起的变化综合而得的结果。一般而言,干扰点的推估是十果。一般而言,干扰点的推估是十分困难的,但是在某些场合下,例分困难的,但是在某些场合下,例如获得如获得人类大规模活动或其他外界人类大规模活动或其他外界因素因素引起干扰变化发生时刻的信息引起干扰变化发生时刻的信息时,干扰点比较明显时,干扰点比较明显,易于确定易于确定。这就是成因分析。这就是成因分析。研究序列
20、研究序列x1,x2,xn,参证序列参证序列y1,y2,yn。分别计算其时序累计值。分别计算其时序累计值时序累计值相关曲线法时序累计值相关曲线法njxgjttj,2,1,1njymjttj,2,1,1点绘点绘mj与与gj的关系图,由图形进行观察。若出现转折,转折点即为跳跃的关系图,由图形进行观察。若出现转折,转折点即为跳跃点点(突变点突变点)。gjmj 当干扰点并不明显,必须统计推估。所谓统计推估指按某种事先规定当干扰点并不明显,必须统计推估。所谓统计推估指按某种事先规定的目标,通过统计分析推估出水文序列最可能的干扰点。在干扰点统计推的目标,通过统计分析推估出水文序列最可能的干扰点。在干扰点统计
21、推断中有两种方法;一是断中有两种方法;一是有序聚类分析法有序聚类分析法;另一是;另一是T检验统计量最大值法。检验统计量最大值法。在分类时若不能打乱次序,这样的分类称为有序分类。以有序分类来在分类时若不能打乱次序,这样的分类称为有序分类。以有序分类来推估最可能的干扰点推估最可能的干扰点 0,其实质是求最优分割点,是同类之间的离差平方其实质是求最优分割点,是同类之间的离差平方和较小而类与类之间的离差平方和较大。对于水文序列和较小而类与类之间的离差平方和较大。对于水文序列x1,x2,xn,最优最优二分割法的要点如下:二分割法的要点如下:设可能干扰点为设可能干扰点为(如右下图如右下图),则干扰前后的离
22、差平方和可表示为,则干扰前后的离差平方和可表示为 iixxV12)(nininxxV12)(xnxn-nVVS)()min11(SSnn 这样总离差平方和为这样总离差平方和为 其中,其中,和和 分别为前后两部分的均值。分别为前后两部分的均值。那么那么,当当 时时,为为最优二分割最优二分割。统计推断统计推断)1()2(2,21(2212121nnNnnnnnnnNK AABBABAAABBBABAABB 将序列将序列x1,x2,x,x+1,x+2,xn的前后两部分用不同字母表示,统计的前后两部分用不同字母表示,统计游程个数游程个数k。当当n1、n220时,时,k趋于正态分布趋于正态分布构造统计量
23、构造统计量游程检验法游程检验法 x1,x2,x-1,x,x+1,x+2,xnF1(x)n1F2(x)n2游程(游程(running):):连续出现同字母的序列连续出现同字母的序列游程长度:每个游程中字母的个数游程长度:每个游程中字母的个数游程个数:游程总个数游程个数:游程总个数)1,0()1()2(2)21(2212121NnnnnnnnnnnkU 假设假设F1(x)=F2(x)。给定显著水平。给定显著水平 后,查算后,查算U /2。当当 时,接受时,接受原假设;反之,原假设;反之,F1(x)不等于不等于 F2(x),有有显著差异。显著差异。2/UU 实例:实例:序列序列xt:11 9 7 1
24、2 14 15 16 10 13 n1=4 n2=5 A B从小到大排序:从小到大排序:7 9 10 11 12 13 14 15 16 A A B A A B B B B 统计游程个数统计游程个数k=4。计算统计量计算统计量U)19(9)9542(542)95421(4)1()2(2)21(22212121nnnnnnnnnnkU 当检验出跳跃成分后,当检验出跳跃成分后,可用适当方法排除。对于可用适当方法排除。对于以年为时间尺度的水文序以年为时间尺度的水文序列,其跳跃一般表现在均列,其跳跃一般表现在均值或方差上。若仅表现在值或方差上。若仅表现在均值上,进行中心化处理均值上,进行中心化处理即可
25、(图即可(图a);若还表现在);若还表现在方差上,再进行标准化处方差上,再进行标准化处理(图理(图c)。对于时间尺度)。对于时间尺度小于年的水文序列,其跳小于年的水文序列,其跳跃除表现在均值、方差上,跃除表现在均值、方差上,还可能表现在自相关系数还可能表现在自相关系数上。均值、方差上的跳跃上。均值、方差上的跳跃可用标准化方式排除,自可用标准化方式排除,自相关系数上的跳跃可以通相关系数上的跳跃可以通过一定模型排除。过一定模型排除。三三 周期成分识别周期成分识别周期图法周期图法方差谱密度图法方差谱密度图法累计解释方差图法累计解释方差图法(一一)周期成分识别周期成分识别 (二二)周期成分的提取周期成
26、分的提取 周期确定后相应的周期成分就识别出来了。例如,对黄河陕县站周期确定后相应的周期成分就识别出来了。例如,对黄河陕县站1919年年-1981年实测年径流量序列年实测年径流量序列Xt进行分析,发现有进行分析,发现有3年的周期成分,即年的周期成分,即 ttPt32sin0.5332cos3.192.502 对于季节性水文序列,其周期成分可能表现在均值、均方差和自相对于季节性水文序列,其周期成分可能表现在均值、均方差和自相关系数上。关系数上。此排除周此排除周期的方法期的方法称为参数称为参数法法 四四 随机成分识别随机成分识别 除去周期成分除去周期成分Pt、非周期成分、非周期成分Nt的剩余部分的剩
27、余部分St=Xt-Nt-Pt一般为平稳随一般为平稳随机序列。机序列。对于平稳随机序列,主要任务是判断其是独立的还是相依的。是独对于平稳随机序列,主要任务是判断其是独立的还是相依的。是独立的,称为独立平稳随机序列(纯随机序列)。例如年最大流量序列,立的,称为独立平稳随机序列(纯随机序列)。例如年最大流量序列,年最大年最大3日暴雨量序列等。在随机水文学中常用正态分布型、对数正态日暴雨量序列等。在随机水文学中常用正态分布型、对数正态分布型、分布型、P-型等概率模型描述纯随机序列。型等概率模型描述纯随机序列。若若St是相依的,称为相依平稳随机序列。水文序列中常存在着一种是相依的,称为相依平稳随机序列。
28、水文序列中常存在着一种持续变化现象,如许多河流径流年际变化就存在持续丰水年组与枯水年持续变化现象,如许多河流径流年际变化就存在持续丰水年组与枯水年组交替出现的现象。这种持续变化现象,表明序列存在着相依性。年径组交替出现的现象。这种持续变化现象,表明序列存在着相依性。年径流序列的相依性,除了因气候因素造成的原因外,主要与流域的地表和流序列的相依性,除了因气候因素造成的原因外,主要与流域的地表和地下水库对径流的调蓄能力有关,如流域有融雪稳定补给,湖沼度和植地下水库对径流的调蓄能力有关,如流域有融雪稳定补给,湖沼度和植被度大,岩性和土壤有利下渗和持水等。因此,在地表和地下径流均有被度大,岩性和土壤有
29、利下渗和持水等。因此,在地表和地下径流均有较强调蓄的情况下,径流年际之间可能有较好的相依关系,自相关系数较强调蓄的情况下,径流年际之间可能有较好的相依关系,自相关系数就较大。就较大。大学是知识的源泉,少大学是知识的源泉,少数的人如饥似渴地畅饮数的人如饥似渴地畅饮,更多的人优哉游哉地,更多的人优哉游哉地品吮,绝大多数的人则品吮,绝大多数的人则只是漱了漱口。只是漱了漱口。3.5 纯随机序列的随机模拟纯随机序列的随机模拟设有一随机模型:设有一随机模型:问:问:如何生成随机序列如何生成随机序列x1,x2,xn(n=200)?ttttxaxax2211 对于纯随机序列,常用正态分布型、对数正态分布型、对
30、于纯随机序列,常用正态分布型、对数正态分布型、P-型等概率模型描述。型等概率模型描述。问:问:如何生成纯随机序列如何生成纯随机序列x1,x2,xn(n=500)?纯随机序列纯随机序列 1)当当 t服从正态分布时,如何随机模拟?服从正态分布时,如何随机模拟?2)当当 t服从对数正态分布时,如何随机模拟服从对数正态分布时,如何随机模拟?3)当当 t服从服从PIII型分布时,如何随机模拟?型分布时,如何随机模拟?一一 概述概述 将生成随机序列将生成随机序列x1,x2,xn的过程,称为随机模拟的过程,称为随机模拟(stochastic simulation,generation),有时又称,有时又称统
31、计试验统计试验或或Monte Carlo法法。PQmPIII型频率曲线型频率曲线(随机模型)(随机模型)(1)如何生成如何生成Pi?Pi均匀随机数均匀随机数Pi(2)如何将如何将Pi转换成指定分布的随机序列转换成指定分布的随机序列?两个问题两个问题:随机模拟这个概念早就有了,并在各个领域得到了大量的应用。随机模拟这个概念早就有了,并在各个领域得到了大量的应用。Qmi二二 均匀随机数的模拟均匀随机数的模拟 在在0,1区间上呈均匀分布的随机序列,以区间上呈均匀分布的随机序列,以u1,u2,un表示,称为均匀表示,称为均匀随机数。随机数。均匀随机数具有如下特征:均匀随机数具有如下特征:独立性独立性、
32、均匀性均匀性。均匀随机数有很多方法可生成,如均匀随机数有很多方法可生成,如随机数表随机数表、物理随机数发生器、物理随机数发生器、乘乘同余法同余法等。等。乘同余法乘同余法生成公式为生成公式为MxuiMxxiiii),2,1(),(mod(1式中,式中,(乘子乘子)、M(模模)、x0(初始值)为选定的常数;(初始值)为选定的常数;ui为为均匀随机数。均匀随机数。xi-1被被M整除整除后的余数后的余数、M、x0的优选的优选 书中三组参数可供选择。书中三组参数可供选择。不同的不同的xi、ui最多只有最多只有M个,因此,产生的均匀随机数不是真正意义上个,因此,产生的均匀随机数不是真正意义上的均匀随机数,
33、称为的均匀随机数,称为伪均匀随机数伪均匀随机数。0,1 上均匀分布的随机数表上均匀分布的随机数表实例实例 1 解解 x0=71=7,x1=7,u1=x1/M=0.007 x1=77=49,x2=49,u2=x2/M=0.049 x2=749=343,x3=343,u3=x3/M=0.343 x3=7343=2401,x4=401,u4=x4/M=0.401 x4=7401=2807,x5=807,u5=x5/M=0.807 x5=7807=5649,x6=649,u6=x6/M=0.649 设设=7、M=103、x0=1。请生成。请生成10个随机数?个随机数?MxuiMxxiiii),2,1(
34、),(mod(1实例实例 2 需要时,一是可以自己编写代码,而是可以使用内部函数。下面给需要时,一是可以自己编写代码,而是可以使用内部函数。下面给出了出了代码。代码。设设 =7、M=10、x0=1。请生成。请生成10个随机数?个随机数?解解 x0=71=7,x1=7,u1=x1/M=0.7 x1=77=49,x2=9,u2=x2/M=0.9 x2=79=63,x3=3,u3=x3/M=0.3 x3=73=21,x4=1,u4=x4/M=0.1 x4=71=7,x5=7,u5=x5/M=0.7 x5=77=49,x6=9,u6=x6/M=0.9 x6=79=63,x7=3,u7=x7/M=0.3
35、 x7=73=21,x8=1,u8=x8/M=0.1 x8=71=7,x9=7,u9=x9/M=0.7 FUNCTION RAND(IX)INTEGER A,P,IX,B15,B16,XHI,XALO,LEFTLO,FHI,K DATA A/16807/,B15/32768/,B16/65536/,P/2147483647/XHI=IX/B16 XALO=(IX-XHI*B16)*A LEFTLO=XALO/B16 FHI=XHI*A+LEFTLO K=FHI/B15 IX=(XALO-LEFTLO*B16)-P)+(FHI-K*B15)*B16)+K IF(IX.LT.0.D0)IX=IX+
36、P RAND=FLOAT(IX)*4.656612875E-10 RETURN END(1)模模M 与随机数的周期有关,愈大,周期愈长。要求计算简便。与随机数的周期有关,愈大,周期愈长。要求计算简便。与计算机性能(字长)有关。与计算机性能(字长)有关。M=2k最大周期为最大周期为2k-2(2)初值初值x0 要保证有有最大周期,要保证有有最大周期,x0需奇数,通常取需奇数,通常取1或或3。、M、x0的选取原则的选取原则(3)乘子乘子 =52s+1其中,其中,s为正整数,满足为正整数,满足 52s+12k 52s+331k34 s=6;35k39 s=740k44 s=8;45k0.5时,宜采用舍选法。时,宜采用舍选法。分母小于分母小于或等于或等于1作作 业业第第3章第章第3题题本章小结本章小结
侵权处理QQ:3464097650--上传资料QQ:3464097650
【声明】本站为“文档C2C交易模式”,即用户上传的文档直接卖给(下载)用户,本站只是网络空间服务平台,本站所有原创文档下载所得归上传人所有,如您发现上传作品侵犯了您的版权,请立刻联系我们并提供证据,我们将在3个工作日内予以改正。