一种未来情景下河口设计高潮位定量研究方法

摘要

本发明公开的一种未来情景下河口设计高潮位定量研究方法，包括以下步骤：步骤S10，利用Mann‑Kendall法对经过修正的河口潮位的长系列观测数据进行处理，以判断河口潮位的趋势性；步骤S20，采用Markov Chain Monte Carlo算法并结合Metropolis‑Hastings抽样构建随机模型，并通过随机模型模拟潮位序列；步骤S30，对历史潮位序列和模拟潮位序列进行组合，得到一个新潮位序列，并对新潮位序列进行检验；步骤S40，根据新潮位序列计算初值，并利用P‑III型曲线进行目估适线，最终得到设计高潮位。本发明为河口地区设计高潮位的预测研究短板提供了一种定量解决方案，其计算结果合理可信，在防汛防潮工作相关及相应工程设计方面均具有较强的实用价值。

说明书

技术领域

本发明涉及水利工程技术领域，尤其涉及一种未来情景下河口设计高潮位定量研究方法。

背景技术

设计潮位是海岸防灾工程设计中的一个重要水文数据，它不仅直接影响着堤防工程高程的确定，而且影响到建筑物类型的选择以及结构计算等。海岸防灾工程的规模、等级和使用情况不同，选用的设计潮位也不同。过去，我国一些单位在设计中，堤防工程的设计高水位采用历年最高潮位，在实测资料的年限较短的情况下，历年最高潮位则根据调查和论证确定。我国有关单位经过大量潮位资料分析比较后，建议采用P-III等年频率统计的方法来确定设计潮位。

极端高潮位是影响设计高潮位的关键因素。受海平面上升和地面沉降等诸多因素影响，现状的设计高潮位极有可能满足不了未来的防洪要求，将给堤防和人民生命财产安全带来极大威胁。未来情景下河口的设计高潮位受到研究人员关注，并已取得了诸多研究成果，提出了海平面上升和潮位变化的回归关系等方式。但大多数研究仅考虑海平面上升等因素带来的影响，但并没有进行定量分析，特别是未针对极端高潮位出现频次对设计高潮位影响做定量分析，给实际应用带来了诸多不便。

为此，本申请人经过有益的探索和研究，找到了解决上述问题的方法，下面将要介绍的技术方案便是在这种背景下产生的。

发明内容

本发明的目的在于：提供一种方便实际应用的未来情景下河口设计高潮位定量研究方法。

为了实现上述目的，本发明可以采用如下技术方案：

一种未来情景下河口设计高潮位定量研究方法，包括以下步骤：

步骤S10，利用Mann-Kendall法对经过修正的河口潮位的长系列观测数据进行处理，以判断河口潮位的趋势性；

步骤S20，采用Markov Chain Monte Carlo算法并结合Metropolis-Hastings抽样构建随机模型，并通过随机模型模拟潮位序列；

步骤S30，对历史潮位序列和模拟潮位序列进行组合，得到一个新潮位序列，并对新潮位序列进行检验；

步骤S40，根据新潮位序列计算初值，并利用P-III型曲线进行目估适线，最终得到设计高潮位。

在本发明的一个优选实施例中，在步骤S10中，所述利用Mann-Kendall法对经过修正的河口潮位的长系列观测数据进行处理，包括以下步骤：

步骤S11，定义Mann-Kendall检验统计量S：

式中，sign0为符号函数。当X

步骤S12，统计量公式S大于、等于、小于零时分别为：

Z为正值表示增加趋势，负值表示减少趋势，Z的绝对值在大于等于1.28、1.64、2.32时表示分别通过了信度90％、95％、99％显著性检验。

在本发明的一个优选实施例中，在步骤S20中，所述采用Markov Chain MonteCarlo算法并结合Metropolis-Hastings抽样构建随机模型，并通过随机模型模拟潮位序列，包括以下步骤：

步骤S21，依据序列Y

步骤S22，马尔科夫链的状态转移的概率只依赖于前一个状态，其数学表达为：

P(X

步骤S23，由于马氏链的收敛的行为和初始概率分布无关，仅由概率转移矩阵P决定。当n足够大的时候，矩阵P

通常情况下P

记Q

步骤S24，为避免因为接受率偏小可能导致的MCMC大量产生拒绝跳转，缩减马氏链遍历所有状态空间所花费的时间，加快收敛速度，本发明在MCMC基础之上，将接受率α

取新的接受率：

在本发明的一个优选实施例中，在步骤S30中，所述对历史潮位序列和模拟潮位序列进行组合，得到一个新潮位序列，并对新潮位序列进行检验，包括以下步骤：

步骤S31，累计偏差计算方式

步骤S32，计算调整部分：

式中，D

在本发明的一个优选实施例中，在步骤S40中，根据新潮位序列计算初值，并利用P-III型曲线进行目估适线，最终得到设计高潮位，包括以下步骤：

步骤S41，根据新潮位序列，利用无偏估计计算样本均值、变异系数和偏态系数作为目估适线的初值；

步骤S42，利用P-III型曲线不断调整参数使得曲线拟合较好：

P-III型曲线为伽马分布，其数学表述为：

式中，Γ(α)为α的伽马函数，α、β、α

式中，

步骤S43，根据指定的设计保证率求得对应的设计高潮位和重现期。

由于采用了如上技术方案，本发明的有益效果在于：

1.本发明提出了基于MCMC和Metropolis-Hastings抽样的随机模型，生成了理论上的未来年最高潮位序列，便于对序列进行分析，按工程规范要求计算设计潮位。

2.本发明对随机生成的潮位进行一致性、趋势性等检验，筛选淘汰不符合水文规律的序列，保留与原实测序列具备相同趋势和一致性的结果。

3.本发明得到的设计潮位，可以多次批量计算，在大量计算下可以看出设计潮位的趋势，削减未来潮位模拟的人为干扰因素，并能针对极端高潮位等一些安全条件，进行相应的模拟预测。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1是本发明的未来情景下河口设计高潮位定量研究方法的流程图。

图2是上海市吴淞口区位示意图。

图3是上海市吴淞站经修订后的长系列观测资料。

图4是上海市吴淞站1997年之前P-III曲线目估适线拟合图。

图5是上海市吴淞站2018年P-III曲线目估适线拟合图。。

图6是保守性物质团单一方向纯扩散问题下扩散系数与物质团0.5mg/L等值线包络范围的关系。

具体实施方式

为了使本发明实现的技术手段、创作特征、达成目的与功效易于明白了解，下面结合具体图示，进一步阐述本发明。

参见图1，图中给出的是一种未来情景下河口设计高潮位定量研究方法，包括以下步骤：

步骤S10，利用Mann-Kendall法对经过修正的河口潮位的长系列观测数据进行处理，以判断河口潮位的趋势性。

Mann-Kendall检验是世界气象组织推荐并已被广泛使用的非参数检验方法，最初由Mann和Kendall提出，可用来分析潮位时间序列的趋势变化。

在步骤S10中，利用Mann-Kendall法对经过修正的河口潮位的长系列观测数据进行处理，包括以下步骤：

步骤S11，定义Mann-Kendall检验统计量S：

式中，sign()为符号函数。当X

步骤S12，统计量公式S大于、等于、小于零时分别为：

Z为正值表示增加趋势，负值表示减少趋势，Z的绝对值在大于等于1.28、1.64、2.32时表示分别通过了信度90％、95％、99％显著性检验。

步骤S20，采用Markov Chain Monte Carlo算法并结合Metropolis-Hastings抽样构建随机模型，并通过随机模型模拟潮位序列。

马尔科夫链蒙特卡洛方法(Markov Chain Monte Carlo)，简称MCMC，该方法是将马尔科夫(Markov)过程引入到Monte Carlo模拟中，实现抽样分布随模拟的进行而改变的动态模拟。

在步骤S20中，所述采用Markov Chain Monte Carlo算法并结合Metropolis-Hastings抽样构建随机模型，并通过随机模型模拟潮位序列，包括以下步骤：

步骤S21，依据序列Y

步骤S22，马尔科夫链的状态转移的概率只依赖于前一个状态，其数学表达为：

P(X

步骤S23，由于马氏链的收敛的行为和初始概率分布无关，仅由概率转移矩阵P决定。当n足够大的时候，矩阵P

通常情况下P

记Q

步骤S24，为避免因为接受率偏小可能导致的MCMC大量产生拒绝跳转，缩减马氏链遍历所有状态空间所花费的时间，加快收敛速度，本发明在MCMC基础之上，将接受率α

取新的接受率：

步骤S30，对历史潮位序列和模拟潮位序列进行组合，得到一个新潮位序列，并对新潮位序列进行检验，其中检验包括一致性检验、趋势性检验等。

随机模型生成的序列不一定同原序列具有一致的趋势，也不一定与原序列直接连续。为此，需要对新的组合序列进行分析，舍去不符合一致性和相同趋势的序列。趋势性仍用Mann-Kendall检验法，一致性采用累计偏差法分析。

在步骤S30中，对历史潮位序列和模拟潮位序列进行组合，得到一个新潮位序列，并对新潮位序列进行一致性、趋势性等检验，包括以下步骤：

步骤S31，累计偏差计算方式

步骤S32，计算调整部分：

式中，D

步骤S40，根据新潮位序列计算初值，并利用P-III型曲线进行目估适线，最终得到设计高潮位。

在步骤S40中，根据新潮位序列计算初值，并利用P-III型曲线进行目估适线，最终得到设计高潮位，包括以下步骤：

步骤S41，根据新潮位新序列，利用无偏估计计算样本均值、变异系数和偏态系数作为目估适线的初值；

步骤S42，利用P-III型曲线不断调整参数使得曲线拟合较好：

P-III型曲线为伽马分布，其数学表述为：

式中，Γ(α)为α的伽马函数，α、β、α

式中，

步骤S43，根据指定的设计保证率求得对应的设计高潮位和重现期。

为了更好地说明本发明的各个步骤，本发明针对极端高潮位情景对上海吴淞口设计高潮位影响研究进行分析说明，具体说明如下：

已知吴淞口历史上最高潮位为9711台风期间的5.99m，根据未来30年吴淞口可能再出现5.99m极端高潮位频次(0～2次)的设计情景，分别分析3种未来尺度下的极端高潮位情景对吴淞口设计高潮位(千年一遇)的影响。

参见图1，本发明的未来情景下河口设计高潮位定量研究方法，包括以下步骤：

1.利用Mann-Kendall法对经过修正的河口潮位的长系列观测数据进行处理，以判断河口潮位的趋势性。

吴淞站地理区位如图2所示，其经修正的多年观测潮位序列S1如图3所示。M-K检验统计值为-1.997，通过95％置信度的显著性检验，序列具有下降的趋势。

2.采用Markov Chain Monte Carlo算法并结合Metropolis-Hastings抽样构建随机模型，并通过随机模型模拟潮位序列。

吴淞口历史上出现过3次5.80m以上的极端高潮位，分别发生在8114号台风、9711号台风和0012号台风期间。极端高潮位平均40年发生一次，且近20年未发生，未来30年极端高潮位存在再发生0～2次的可能性。吴淞口历史上出现的最高潮位为5.99m。据此，设计未来30年吴淞口可能再出现5.99m极端高潮位频次(0～2次)，分别对应3种未来尺度下的极端高潮位情景，详见表1所示。

表1极端高潮位设计情景

对历史情景，在95％置信水平下用Mann-Kendall法检验趋势性，结果分别为-2.203；对于历史和现状两种情景，分别采用公式(12)计算初值，经目估后的P-III拟合曲线及参数如图4和图5所示，两种情景下的千年一遇设计高潮位分别为6.389m和6.582m。

依据序列S1，计算Markov Chain单步转移矩阵为：

单步转移频率矩阵

由此迭代可以计算得到平稳分布，进而可以构建随机模型，进行得到模拟序列S2。

3.对历史潮位序列和模拟潮位序列进行组合，得到一个新潮位序列，并对新潮位序列进行检验。

对表1中所列的3种未来情景，分别对其S1+S2序列进行检验，结果如表2所示。

表2各情景检验及潮位设计成果

4.根据新潮位序列计算初值，并利用P-III型曲线进行目估适线，最终得到设计高潮位。

初值计算如表2所示，采用Cs与Cv比值(25.0)固定不变、无偏估计计算均值、变异系数同倍比放大的方式分别计算其千年一遇设计高潮位，如表2最后一列所示。

为消除随机模拟带来的误差，保证定量研究的代表性，对三种未来情景，分别采用相同方法，随机生成满足相同置信水平下的一致性、趋势性和马尔科夫性的1000组随机序列，计算生成的设计高潮位情况如图6所示。三种情景序列的均值分别为6.605m、6.717m、6.823m；标准差分别为0.034m、0.033m、0.031m，基本可以认为吴淞口未来30年内发生0次、1次、2次极端高潮位情景下，吴淞口2050年千年一遇设计高潮位约为6.60m、6.72m、6.82m。三种情景的1000次计算结果分别集中在一个误差较小的范围，可见采用本发明，吴淞口未来情景模式下的设计潮位预测问题得到了很好的解决。

以上显示和描述了本发明的基本原理和主要特征和本发明的优点。本行业的技术人员应该了解，本发明不受上述实施例的限制，上述实施例和说明书中描述的只是说明本发明的原理，在不脱离本发明精神和范围的前提下，本发明还会有各种变化和改进，这些变化和改进都落入要求保护的本发明范围内。本发明要求保护范围由所附的权利要求书及其等效物界定。