|
2023年全国大学生数学建模竞赛E题“黄河水沙监测数据分析”（见全国大学生数学建模竞赛网站：http:// www.mcm.edu.cn/）中数据来源于黄河下游的孙口水文站监测的实际数据，为了适合竞赛的要求略作修改。对黄河水沙监测数据分析的主要目的是结合当前国家战略的迫切需求，探究黄河水沙通量变化规律及其对气候变化和人类活动的响应，了解黄河泥沙的变化情况，为优化黄河流域水资源分配、协调人地关系、调水调沙、防洪减灾等提供理论指导。命题人在《数学建模及其应用》杂志2023年第4期中刊发论文“黄河水沙监测数据分析中的问题探讨”较为详细地介绍了赛题的背景。这里只做简要地介绍，内容摘自于吕教授等人的论文，详细情况可参阅：
1 i" l2 k* J1 B9 F" t$ b 吕良军,李雪婷,郝炯驹. 黄河水沙监测数据分析中的问题探讨[J]. 数学建模及其应用，2023，12(4),44-448.黄河水沙监测数据分析中的问题探讨_吕良军.pdf
% k0 T0 x3 T |2 W, a9 H \ 1 孙口水文站的基本情况 {1 T8 v1 f3 w: x
孙口水文站始建于1949年8月21
" F6 V9 [1 ?4 K% d! ^ 日，是国家重点报汛站、东平湖分滞洪决策控制站，具有常年向国家防总、黄河防总上报水情的任务。 " g# ~ A* v# L; Q) `9 i
& V/ \' y, E, w# G# U9 Y2 a2 H 测站地址位于山东省梁山县赵堌堆乡蔡楼村，集水面积为734824km2，距河源距离为5105km，距入海口距离为449km- w) Y( ]' c! E4 V' Y {
# C" U0 L s! p% t
，详见下图：
: f. ~. O7 K/ i
) k1 Q, m$ o( n 
, S4 N, X2 o* _8 F5 U5 R 2 河段特性及过流能力情况 6 S; Q+ Z/ \, j& |
2.1测验河段 测验河段长约1000m，河道比降约为千分之1.1，测流断面位于弯道下500m处，属于两滩一槽的复式河槽，沙质河床，滩槽流向变异较大。测验河段呈“S”形，断面左岸上游2km处为梁路口护滩工程，断面右岸上游500m处为蔡楼控导工程，左岸下游2.5km处有影堂险工，河道断面较稳定。
4 ^ n1 X( l* b7 C+ T 2.2断面布设与测量 断面布设呈“Z”字形折线，断面总宽度4630m，其中主槽即两个生产堤之间，宽690m其主槽内常年过水宽度500m左右，左右两岸皆有嫩滩；左边滩宽1484m，右边滩宽2456m。该站冲淤变化较大，主流摆动较频繁。上比降断面位于基本水尺断面上游1000m，下比降断面与基本水尺断面重合，测流断面位于基本水尺断面上游145m。黄河右岸大堤堤顶高程53.25m，左岸大堤堤顶高程53.65m，设防流量为17500m³/s。左右岸滩地各有一道生产堤，左岸生产堤高程为48.54m；右岸生产堤高程为48.23m。自2003年以后，主槽趋于平坦，主流均匀化。随着水位抬高，流量增大河床两边冲刷，形成两边深槽，中间较平坦的河床。: R3 l; L* N; n- z" ]5 V: |
2.3过流能力 # I! A1 z0 p* { U
2002年调水调沙前，孙口河段平滩流量仅为1800m³/s。经过历年调水调沙以后，2022年平滩流量5100m³/s,相应水位46.68m，河道行洪能力明显增强，如下图所示: 7 i, I. U1 K- q/ W! T
4 N; i/ [. m. `' ~0 j
/ J- u/ _1 Q* S& I
黄河孙口站历年调水调沙主槽过流能力趋势 由孙口站“调水调沙”后过流能力趋势图可以看出：水位～流量关系曲线逐渐右移，同水位流量增加，2022年较2002年同水位流量偏大2720m³/s，同水位45.00m对应的2002年相应流量556m³/s，2022年相应流量3280m³/s，过流能力增加明显。
0 K D& ]6 Z9 q 2.4 断面历年冲积淤情况 自1976年汛前至2002年汛前，孙口水文站所处的河道主槽河底高程提高显著，2002年开始“调水调沙”至2022年6月24日，经数次“调水调沙”及较大洪水后，孙口水文站在这多年断面测验发现,主槽断面冲刷深度为1.69m,冲刷面积约为1018.61,达到了黄河下游河底不抬高的要求，见下图:, q9 K2 x6 m) }! l
2 X5 B6 u6 q# m" c! w0 E 孙口水文站历年调水调沙主槽断面冲淤变化
3 d7 ~0 k2 a1 `1 b" {6 P 3 近年来水文站径流量与输沙量特征分析 / U! l- r' q- `0 ? v
3.1水沙年内分配 孙口水文站近年（2016-2021）的年平均径流量为330.01亿, 其中最大月平均径流量出现在7月, 占全年径流量的15.52%, 最小月平均径流量出现在１月,占全年径流量的2.73%.年平均输沙量为20213.70万, 其中最大月平均输沙量出现在7月, 占全年输沙量的31.71%, 最小月平均输沙量出现在1月, 仅占全年输沙量的0.43%.孙口水文站输沙量和径流量年内变化趋势基本一致, 但输沙量和径流量年内分配不均匀, 输沙量主要集中在6-10月的汛期, 占全年输沙量的82.40%,同期径流量占全年径流量的61.04%。
9 E* {, X" F! r9 Z 3.2趋势与周期分析 对近6年的水沙通量按月统计，通过可视化分析、滑动平均和曼-肯德尔（Mann-Kendall）等方法分析研究水文站的径流量、输沙量的变化趋势和周期性等水文特征。利用曼-肯德尔法对2016-2021年月径流量和月输沙量做趋势检验。月径流量趋势检验整体呈现显著增加趋势，输沙量的趋势整体也为增加趋势, 但增加趋势不显著。可以显著发现月径流量和月输沙量具有周期性，且周期为12个月，即以年为一个周期。
9 q/ l- L5 c0 \/ U3 P: j/ } 4 测验时遇到的问题
. q8 D$ M/ u& Q0 G 4.1 测验方案优化 如下图所示，目前水文站在流量与含沙量监测采样时先布线，再在线上按比例设采样点，即“多线多点”法采样，此法的施测方案误差小，精度高，但是在实际流量测验工作中，如果施测时间过长，往往水位会发生变化，也会受到其他不稳定因素影响，反而降低了测验精度。因此，需要优选误差小、精度高、施测时间短、减少工作量、经济适用的流量施测方案。
# I3 z* C, k% T4 E8 ] 
断面测验布线: ~; V, C5 f9 o3 p2 L" I/ P. A
4.1.1实测流量计算模型 水文站最常用的流量测验方法，其将河道断面划分为若干部分，测算出各部分断面的面积，并用流速仪测量、计算出各部分面积的平均流速，两者的乘积，即为部分流量，各部分流量的和为全断面的流量。
4 T0 S7 \$ @7 Y- j2 C 
由于河水流速与水深、离岸距离函数式不易确定，实际测量时采用多线多点测水深和流速，再按如下方法进行计算断面流量：
; n; ^2 F3 Z) s 
 对2021年7月9日的测验进行分析计算。首先是测点的布设，孙口水文站采用缆道流速仪常没法，利用23条固定垂线、每条垂线上在水面、相对水深0.2、0.6、0.8和河底处布设5点法测速，如下图所示。测验结果数据详见2023全国大学生数学建模E题附件3
- S6 V/ ^, ]6 V; n" c8 E! C6 X 。
' L& v+ Z/ S$ t) K6 @! a' {' F
& p7 {( h+ m. z* j+ a- m! T 
' N# [1 `' P, O$ m. _5 v: G 孙口水文站测验测点布设 $ w. b" l1 g" [: t& a0 R9 X
然后用实测流量计算模型，可对上面的测验数据进行计算。
& q9 l( Y1 w0 E 4.1.2实测流量优化问题 问题：每次布线与测验都需较长时间，由于历时较长，水位会发生变化，测验过程也会受到其他不稳定因素影响，往往会产生系统误差，反而降低了测验精度。因此，需要优选误差小、精度高、施测时间短、减少工作量、最经济适用的流量施测方案。4.1.3实测含沙量计算模型
. V" _9 q: \7 W* ]8 a7 X 
& N9 |: f; {$ m 
0 Z& p7 s3 I/ g Q! O! F' P/ ?
由于含沙量与河水流速、水深、离岸距离函数式不易确定，实际测量时选用全断面采集水样，即采用多线多点测水深、流速、含沙量，再按如下方法进行计算断面平均含沙量。 2 [0 q) C, T7 F1 g
6 v% a$ A4 B' M0 Z! x5 r5 ? 
, {) Z$ k" ?( `8 z9 \ Y8 L5 {( D
含沙量测验垂线布设，如下图所示
- Q4 U, u9 H: \# F : 0 G/ v j# S n( X$ J7 a' H+ U Y
1 F% t; I! E! k8 j
这里布设12条垂线，每条垂线上分别在水面、相对水深0.2、0.6、0.8和河底处布设5点取水，经24小时沉淀后测深水的质量，再用比对法测出水样的含沙量。利用上述计算含沙量计算公式，对2021年7月9日的测含沙量相关数据进行计算。
& o& j; J8 o! n2 \ 4.1.4 实测含沙量优化问题 测验含沙量与测验流量问题类似，如果布线历时较长，测点水位会发生变化，取样点位置就会出现偏差，另在比对时受气温等其他不稳定因素影响，会产生系统误差，从而降低了测验精度。
/ X2 S1 ~7 Y. f1 {, K 4.2 含沙量与流速关系 由于近年ADCP-声学多普勒流速剖面仪在孙口水文站投入使用，使测量流速和流量工作量减小，它比传统的河流流量测量方法提高效率几十倍，它标志着河流流量测量的现代化。ADCP是用于采集流速大小和方向的声学设备，用ADCP进行河流流速测量是近十年才发展和应用的新的流速和测量方法。如果能依据孙口水文站的水文特征，得到流速与含水沙量之间的关系，同样可以大大减少含沙量监测工作量。含沙量与流速、水深、沉降等诸多因素有关。为了简单起见，可以探讨断面平均水流速度与断面平均含沙量的关系。
8 g9 Z4 R4 x+ m5 }6 T* K P 4.2.1水流断面面积计算 对河底高程进行三次样条插值,假设在调水调沙前，河床没发生大的变化，河底高程基本不变下，对河底高程进行三次样条插值。在对河底高程插值的基础上计算过水断面面积，可以采用离散梯形面积求和法：& H0 a1 z0 o" o; I
8 W- i0 V4 L( v( p5 l: a- g0 h
4.2.2水流断面平均流速计算
: u ^/ U/ U- J* p; B7 o% C 
. Z5 _) b1 {: \! J+ M" ]$ T+ X
% s5 j1 N# d( ^. \
4.2.3拟合含沙与平均流速关系 由实际测量数据计算平均流速，画出平均流速与断面平均含沙量的散点图，分析建立流速与含沙量数学模型：
& ?; c+ q) Z; \/ ~/ F 
经机理分析或对数据作非线性拟合，得出断面平均流速与平均含沙之间的关系。如果验证通过，就可以通过ADCP-声学多普勒流速剖面仪测水流速，再计算含断面水沙量就大大减少工作量，同时还减少了测量时产生的系统误差。
& o9 l# Z; K A( @6 k 4.3“调水调沙”对河床演变的影响 以2022年为例，孙口站6月25日至7月1日，共施测流量6次。实测最大流量3970m3/s，相应水位45.60m，最大含沙量为6.17kg/m3，最大水深9.3m，最大流速3.46m/s,最大断面平均流速2.65m/s。对断面进行了实时套绘，本次洪水前至今，整体呈冲刷形态，冲刷位置由右岸逐渐向左移动，主流开始左移，如下图所示，图中
2 Q% Y3 K$ v2 O 是指第42次测量。 2 {/ H' `2 L% C* Y! t8 c
' p3 F. G" Q7 `1 Z1 [" I3 X. A 
" o- W) H" y$ T9 ]0 b8 F: }* h) |
与洪水前实测流量断面进行了比较分析，从比较的结果来看，截止6月30日，主槽冲刷明显；与6月18日洪水前标准水位下平均河底高程冲刷1.13m，面积淤积冲刷350，具体数据见下表。 9 [( T0 x' F0 I8 |3 _
6 v# l0 y! {) Z 2022年6月18日-6月30日监测数据
5 [( l# a6 G: O' V s
; p5 y* k( n$ {; p, e A9 D) r& Z 
现在的问题是：在依据孙口水文站现有的水文特征和“调水调沙”期间的测验数据，能否较准确的预测出“调水调沙”对河床的冲涮效果以及河床会如何演变？6 _) J$ V0 m! J4 N" S! U
5 结语 % f8 C4 }; f8 a0 W6 g* M
由于含沙量与诸多因素有关，含沙量与水流速回归模型精度不够高，有待进一步改进。因该水文站河底高程数据度跨度不一，加上泥沙淤积机理复杂，预测“调水调沙”对河床冲涮演变的影响有等进一步研究分析。
5 e* }7 a9 m- L# ~% j $ V5 d+ @2 q6 Y" {1 M# T" K9 K+ H
主要参考文献 + q9 I3 U) ]# s$ {# c7 T4 X
[1]魏玲娜,金宇斌,王文等.水文测验中两种测速垂线精简算法比较研究[J].人民长江,2020,51(03):99-102+173.
9 d( v8 w+ U' y }" s, \ [2]张莉,陶武亮,张孝军等.新江口水文站悬移质泥沙测验方案优化研究[J].水利水电快报,2019,40(01):18-22. & \3 n4 A0 P- I7 @
[3]程瑞修.流量测验方法简化分析及应用——以沙里寨水文站为例[J].黑龙江水利科技,2020,48(04):185-187.
' S5 b; {$ P3 W2 o [4]刘欣,刘远征.小浪底水库调水调沙以来黄河下游游荡河段河床演变研究[J].泥沙研究,2019,44(05):56-60.
5 f T' S: o1 a, U [5]许永占,李红薇,李迎春.单断关系法推求断面平均含沙量[J].吉林农业,2013(09):90.
& @& u6 w9 ^$ H# ? J7 ^ [6]李韶星.双桥水文站流量测验测速垂线精简研究[J].黑龙江水利科技,2018,46(01):39-41. 7 U2 f' m: j$ I3 G8 @$ q
[7]全国大学生数学建模竞赛组委会.2023全国大学生数学建模竞赛E题[Z/OL]. [2023-09-26]. http://www.mcm.edu.cn/html_cn/node/c74d72127066f510a5723a94b5323a26.html.
9 A9 `' c3 i0 w3 Y [8]赵兰.韩江潮安水文站水沙特征分析[J].广东水利水电,2019(01):22-25. 8 F* p- X9 h7 I/ Q) o. ]4 T$ u
论文作者介绍 吕良军，男，黄河水利职业技术学院教授，多次为全国大学生数学建模竞赛命题.李雪婷,女,工程师,主要研究方向为水文监测.郝炯驹,男,博士,主要研究方向为人工智能.
* L1 h( ^5 b+ u2 a$ p2 N
1 M7 c* L$ Y) T) c
6 p3 @' j% }+ U7 n$ e
* [ k- w6 W) T; m' T" D3 c: {( D! ]- P$ y+ v
|
