CN109345777A - 基于陡坡汇流和断面流量计算的山洪泥石流预警方法和*** - Google Patents
基于陡坡汇流和断面流量计算的山洪泥石流预警方法和*** Download PDFInfo
- Publication number
- CN109345777A CN109345777A CN201811178161.4A CN201811178161A CN109345777A CN 109345777 A CN109345777 A CN 109345777A CN 201811178161 A CN201811178161 A CN 201811178161A CN 109345777 A CN109345777 A CN 109345777A
- Authority
- CN
- China
- Prior art keywords
- flow
- section
- water
- early
- warning
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Classifications
-
- G—PHYSICS
- G08—SIGNALLING
- G08B—SIGNALLING OR CALLING SYSTEMS; ORDER TELEGRAPHS; ALARM SYSTEMS
- G08B21/00—Alarms responsive to a single specified undesired or abnormal condition and not otherwise provided for
- G08B21/02—Alarms for ensuring the safety of persons
- G08B21/10—Alarms for ensuring the safety of persons responsive to calamitous events, e.g. tornados or earthquakes
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A50/00—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE in human health protection, e.g. against extreme weather
Landscapes
- Life Sciences & Earth Sciences (AREA)
- Environmental & Geological Engineering (AREA)
- General Life Sciences & Earth Sciences (AREA)
- Geology (AREA)
- Business, Economics & Management (AREA)
- Emergency Management (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Alarm Systems (AREA)
Abstract
本发明提供了一种基于陡坡汇流和断面流量计算的山洪泥石流预警方法,包括步骤:1、通过对泥石流单沟及其所在流域的高精度DEM地形数据和遥感数据进行处理,得到各网格单元的散点高程与坡度;2、根据不同泥石流沟的实际地质条件,结合连续长历时降雨数据,得到各网格单元的雨水下渗量;3、采用基于陡坡汇流的计算方法,计算泥石流沟及其所在流域的地表汇流;4、通过不同时间段内降雨量所对应的山洪泥石流淹没范围进行预警。还提供了一种山洪泥石流预警***。本发明应用于山洪、沟谷型泥石流的预警,对于及早形成山洪泥石流预警成果,提升预警效率效果明显。对保障人员安全和应急避险、降低人民生命财产损失有较大的实际意义。
Description
技术领域
本发明涉及灾害预警技术领域,特别涉及一种基于陡坡汇流和断面流量计算的山洪泥石流预警方法和***。
背景技术
泥石流是北京地区造成人员伤亡最大的突发性地质灾害,自1949年以来,北京地区各类突发地质灾害造成的死亡人数已超过600人,直接经济损失达数亿元,其中泥石流灾害共致死500余人。沟谷型泥石流在北京地区分布最为广泛,其形成大体可分为三个阶段:气象诱因(主要是降雨)、形成洪水或坡面泥石流、形成沟谷泥石流。而从第一阶段雨量的累积、土壤饱和发展到第三阶段形成泥石流灾害,是有一段时间过程的,在这段时间内,越早发出预警,就给应急避险留下越多的宝贵时间。
现阶段我国地质灾害的预警预报精度水平依然有限,发布的预警信息主要从气象风险诱发地质灾害风险的角度来考虑,并不具体区分地质灾害类型,且精度一般只到区县级。虽然大范围的区域预警可以有效降低漏报率,防患于未然,但在实际过程中对启动应急预案、人员疏散和应急避险的指导效果并不理想。
而针对沟域的预警,现有技术的一些预警做法是采用设定临界雨量的方法,即在雨量计测得在达到某固定预警值时即报警并发布预警信息。但此方法也存在不合理性,山洪泥石流的形成虽然主要由降雨激发,但前期累积降雨量、实时雨量和雨强都是重要的激发因素,不同沟域地质条件和植被条件都不相同,前期累积降雨量对于土壤饱和程度和稳定程度产生影响,而短时雨强对灾害的启动影响也至关重要,单纯的采用临界雨量值作为预警指标其精准性是有待商榷的。
发明内容
本发明的预警模型其研究范围精度是针对沟域的,用于解决区域预警的精度太低的问题;运用陡坡汇流计算方法,来计算坡面洪水的渗流和汇流过程,实现基于降水量和降水时间的过流断面流量值计算,并计算在该流量值下的下游淹没范围,根据淹没范围的大小和水深高度来设定该断面流量的分级预警阈值,并将其对应的实时雨量或预报雨量设为分级预警的阀值,用以解决山洪泥石流预警指标精准确定和时间效率的问题。为达到上述目的,本发明采取如下技术方案。
基于陡坡汇流和断面流量计算的山洪泥石流预警方法,包括步骤:
S101、通过对泥石流单沟及其所在流域的高精度DEM地形数据和遥感数据进行处理,依据典型地形地物特征对研究区域进行网格划分,得到各网格单元的散点高程与坡度;
S102、根据不同泥石流沟的地形地貌特征,计算不同地质条件下水流的下渗以及不同植被条件下地表的汇流,结合长历时降雨数据,得出各网格单元的下渗量及汇流量;
S103、采用基于陡坡汇流的计算方法,对研究区域内的地表汇流进行水动力过程模拟,基于所述的步骤S101和S102,计算泥石流沟及其所在流域的地表汇流;
S104、通过在泥石流沟的流通区范围内选取一个或多个过流断面,作为预警计算的标准参考断面,计算所述标准参考断面在降雨过程中的流量值,并形成所述标准参考断面流量随降雨及预报过程中时间变化的流量曲线;依据村庄、公路和/或桥梁设施位置,设置参考断面下游的不同预警级别的淹没范围,当所述标准参考断面流量计算值达到一个特定值,使得其下游的水流淹没范围达到预警淹没范围且水深达到某一预先设定值时进行相应级别的预警提醒。
进一步地,所述的典型地形地物特征包括:分水岭、沟道和/或河道。
进一步地,通过所述步骤S101获得各网格单元的散点高程与坡度的具体步骤如下:
S01:勾勒泥石流单沟汇流流域,主要依据分水岭对泥石流单沟进行模型范围确定;
S02:在模型范围基础上,按照等尺度划分任意不规则网格,网格划分依据边界点间距,寻找最优化网格结点;
S03:模型单元网格地形根据实测资料插值获取。
进一步地,所述步骤S102中,采用干旱地区超渗产流的格林-安普顿公式,忽略湿润锋毛细管水压力计算渗流:
f=K(H/Z+1)
式中,f为下渗速率,K为饱和水力传导度,H为地面滞水深度,Z为地面以下入渗饱和水舌厚度;
饱和水力传导度K将地质条件以下渗系数的形式添加至表达式中,在模型的计算过程中,通过系数的率定来对不同的地质条件进行调整。
进一步地,所述步骤S103中,采用基于网格的水动力陡坡汇流计算方法,在溃坝大比降陡变水面水流计算模式的基础上,采用有限体积Riemann通量计算格式,建立基于网格的水动力陡坡汇流计算模型:
其中,h为水深;u为局部坐标x方向的流速;qe为源汇项;α为地形坡度;φ为降雨角度;g为重力加速度;S0为底坡比降;Sf为阻力项;σ为风应力;ρ为空气密度;p为降雨强度;vm为雨滴降落速度;d为下渗强度;i为蒸发强度。
进一步地,所述步骤S104中,所选取的该过流断面在整个降雨(及预报)过程中通过的山洪流量Q表示为以降雨量和时间为输入条件的函数Q=f(q,t;k),q为降雨量,t为时间,k为调整系数,用以调整不同的地质条件和植被条件带来的变化;对于确定的泥石流沟具有确定的地形地貌和地质条件,在降雨时间段内得出该过流断面的流量Q随时间变化的连续数值;根据断面流量Q计算出在断面下游山洪泥石流推进过程和淹没范围A=f(Q,t;z),其中Q为所选过流断面的流量值,z为该过流断面的位置;当A达到了淹没范围的分级预警值Af,且淹没水深h达到了预设值hf,则将其所对应的断面流量Qf设定为该分级预警的阈值,而Qf所对应的当时的雨量qf作为预警阀值。
进一步地,所述基于网格的水动力陡坡汇流计算方法,从三维流体力学的N-S方程出发,沿水深方向积分后获得二维浅水方程,所述二维浅水方程包括连续方程和运动方程,其中:
连续方程:
运动方程:
其中,h为水深,z为水位,z=z0+h,z0为底高程,u,v分别为x和y方向的平均流速,q为源汇项,g为重力加速度,n为糙率;
令H=h,Qx=hu,Qy=hv;选择任意形状的网格,在网格中心处计算水深H;网格周边的边线相当于容器的边壁,称为通道,水体通过通道流入单元,在网格通道处计算流量Q;时间上采取水位-流量交错的计算方式;
所述连续方程,用矢量表示流动速度分量:
其中,表示单位宽度流量;
所述连续方程也可以表示为:
其中,Ai为第i个控制体的面积,m为控制体的通道个数,Qij为控制体i的第j个通道上的单宽流量;Lij为控制体i的第j个通道长度;
所述连续方程的离散形式为:
其中,Hpk为第k个管道计算断面的水深或压力水头;Apk为第k个管道计算断面的过流面积;A′pk第k个管道计算断面的等效底面积,即管段中自由水面的面积;
所述运动方程,用矢量表示流动速度分量:
或者,所述运动方程表示为:
所述运动方程的离散形式为:
其中,Qpk为第k个管道计算断面的单宽流量;Zk2与Zk1分别为计算断面k两侧网格中的水深;dxk为计算断面k两侧网格形心到断面中心的距离之和;
所述运动方程的离散形式适用于水深较大的单元通道。
进一步地,针对某特定泥石流沟,包括:地形处理模块、土壤下渗模块、地表汇流模块和预警模块,其中,
所述地形处理模块用于对泥石流单沟及其所在流域的高精度DEM地形数据和遥感数据进行处理,依据典型地形地物特征对研究区域进行网格划分,得到各网格单元的散点高程与坡度;
所述土壤下渗模块用于根据不同泥石流沟的地形地貌特征,计算不同地质条件下水流的下渗以及不同植被条件下地表的汇流,结合长历时降雨数据,得出各网格单元的下渗量及汇流量;
所述地表汇流模块用于基于陡坡汇流的计算方法,对研究区域内的地表汇流进行水动力过程模拟,并基于所述的地形处理模块与土壤下渗模块,计算泥石流沟及其所在流域的地表汇流;
所述预警模块用于通过在泥石流沟的流通区范围内选取一个或多个过流断面,作为预警计算的标准参考断面,计算参考断面在降雨及预报过程中的流量值,并形成该断面流量随降雨过程中时间变化的流量曲线;依据村庄、公路、桥梁等设施位置,设置参考断面下游的不同预警级别的淹没范围,当参考断面流量计算值达到一个特定值,使得其下游的水流淹没范围达到预警淹没范围且水深达到某一预先设定值时进行相应级别的预警提醒。
附图说明
图1是本发明的基于陡坡汇流和断面流量计算的山洪泥石流预警方法流程图。
图2是本发明的基于陡坡汇流和断面流量计算的山洪泥石流预警方法示意图。
图3是本发明的基于陡坡汇流和断面流量计算的山洪泥石流预警方法中水位H与流量Q的空间布置示意图。
图4是本发明的基于陡坡汇流和断面流量计算的山洪泥石流预警方法中水位H与流量Q的时间交错计算示意图。
具体实施方式
下面结合附图对本发明的内容作进一步详细说明。
针对某特定泥石流沟,结合附图1所示的本发明一种基于陡坡汇流和断面流量计算的山洪泥石流预警方法流程图,以满足模型中各计算部分的功能需求。
汛期时,泥石流沟出现连续性降雨或短时强降雨时,沟道的淹没范围及水流流速是衡量其是否诱发突发性地质灾害的一个重要因素。泥石流在正式形成之前首先形成山区洪水或坡面泥石流,综合考虑地质条件影响参数,根据实时降雨量和预报雨量计算山区某一泥石流沟洪水演进范围及淹没水深,实现对不同时间点上的山洪演进过程变化对比,这对后续形成泥石流的判断具有极大的参考意义。
实现山区泥石流沟的某过流断面随时间变化的流量模拟。泥石流的形成主要由降雨激发,前期降雨量、实时雨量和雨强都是重要的激发因素,单纯的采用临界雨量作为预警指标是不精准的。泥石流沟域内某断面的流量可以综合全面的反映整个降雨过程中所产生的地表渗流、汇流等情况,因此具有更实际的意义。
实现基于降雨量的山洪泥石流预警。根据实时雨量和预报雨量,可计算出实时和未来若干时间内随时间连续变化的泥石流沟域的过流断面流量Q、山洪淹没范围A、流速v、水深h,流向等参数,结合泥石流危险程度等级和淹没范围内是否有居民点、景区、学校、交通干线、重大工程等因素,在充分参考历史发灾资料、人员撤离时间、避险路线等条件,由专业技术人员确定分级预警流量阈值Qf,确定降雨量预警阀值qf,并由地质灾害主管部门依据实际情况生成泥石流灾害分级预警产品并发布。
如图1-4所示,本发明提供一种基于陡坡汇流和断面流量计算的山洪泥石流预警方法,包括如下步骤:
S101、通过对泥石流单沟及其所在流域的高精度DEM地形数据和遥感数据进行处理,依据典型地形地物特征对研究区域进行网格划分,得到各网格单元的散点高程与坡度;
S102、根据不同泥石流沟的地形地貌特征,计算不同地质条件下水流的下渗以及不同植被条件下地表的汇流,结合长历时降雨数据,得出各网格单元的下渗量及汇流量;
S103、采用基于陡坡汇流的计算方法,对研究区域内的地表汇流进行水动力过程模拟,基于所述的步骤S101和S102,计算泥石流沟及其所在流域的地表汇流;
S104、通过在泥石流沟的流通区范围内选取一个或多个过流断面,作为预警计算的标准参考断面,计算所述标准参考断面在降雨过程中的流量值,并形成所述标准参考断面流量随降雨及预报过程中时间变化的流量曲线;依据村庄、公路和/或桥梁设施位置,设置参考断面下游的不同预警级别的淹没范围,当所述标准参考断面流量计算值达到一个特定值,使得其下游的水流淹没范围达到预警淹没范围且水深达到某一预先设定值时进行相应级别的预警提醒。
其中,优选的,步骤S101对泥石流单沟及其所在流域的高精度DEM地形数据和遥感数据进行处理,依据典型地形地物特征(分水岭、沟道、河道等)对研究区域进行网格划分,得到各网格单元的散点高程与坡度。通过该模块计算得到的网格单元地形数据,是之后水动力模拟计算中不可或缺的部分。
所述步骤S101首先勾勒泥石流单沟汇流流域,主要依据分水岭对泥石流单沟进行模型范围确定;在模型范围基础上,按照等尺度划分任意不规则网格,网格划分依据边界点间距,寻找最优化网格结点;模型单元网格地形根据实测资料插值获取。
优选的,步骤S102根据不同泥石流沟的地形地貌特征,主要用来计算不同地质条件下水流的下渗以及不同植被条件下地表的汇流,结合长历时降雨数据,得出各网格单元的下渗量及汇流量。通过该步骤S102计算得到的网格单元下渗情况,是之后汇流计算中不可或缺的部分。
对于渗流的计算,采用干旱地区超渗产流的格林-安普顿公式,忽略湿润锋毛细管水压力:
f=K(H/Z+1),
式中,f为下渗速率,K为饱和水力传导度,H为地面滞水深度,Z为地面以下入渗饱和水舌厚度。
式中系数K为饱和水力传导度,将地质条件以下渗系数添加至公式中,在模型的实际计算过程中,可以通过系数的率定来对不同的地质条件进行调整。
所述方法中步骤S103对研究区域内的地表汇流进行水动力过程模拟计算,采用基于陡坡汇流的计算方法。基于所述步骤S101和步骤S102,计算得到泥石流沟及其所在流域的地表汇流。通过该模块计算得到的汇流情况,便可用于概化山洪泥石流的淹没范围与程度(如附图2)。
采用基于网格的陡坡汇流计算方法。在溃坝大比降陡变水面水流计算模式的基础上,采用有限体积Riemann通量计算格式,建立基于网格的水动力陡坡汇流计算方法,其基本方程如式(1)所示:
其中,h为水深;u为局部坐标x方向的流速;qe为源汇项;α为地形坡度;φ为降雨角度;g为重力加速度;S0为底坡比降;Sf为阻力项;σ为风应力;ρ为空气密度;p为降雨强度;vm为雨滴降落速度;d为下渗强度;i为蒸发强度。
优选的,步骤S104进行预警时,在泥石流沟的流通区范围内选取一个或多个过流断面,作为预警计算的标准参考断面。该过流断面(标准参考断面)在整个降雨过程中通过的山洪流量Q可以表示为以降雨量和时间为输入条件的函数Q=f(q,t;k),q为降雨量,t为时间,k为调整系数,用以调整不同的地质条件和植被条件带来的变化。对于确定的泥石流沟具有确定的地形地貌和地质条件,故在降雨时间段内可以得出该断面的流量Q随时间变化的连续数值。根据断面流量Q可计算出在断面下游山洪泥石流推进过程和淹没范围A=f(Q,t;z),其中Q为所选断面的流量值,z为该断面的位置。故以断面流量为基准,可以得到不同时间段内降雨量所对应的山洪泥石流淹没范围A。综合考虑经济指标、人身财产安全、撤离时间等因素,如果A达到了淹没范围的分级预警值Af,且淹没水深h达到了预设值hf,则将其所对应的断面流量Qf设定为该分级预警的阈值,而Qf所对应的当时的实时雨量qf作为预警阀值,将其所对应的断面流量Qf设定为该分级预警的阀值,而Qf所对应的当时的雨量qf作为预警阀值。
从三维流体力学的N-S方程出发,沿水深方向积分后得到二维浅水方程,所述的二维浅水方程包括连续方程和运动方程,其中:
连续方程:
运动方程:
其中,h为水深,z为水位,z=z0+h,z0为底高程,u,v分别为x和y方向的平均流速,q为源汇项,g为重力加速度,n为糙率。
令H=h,Qx=hu,Qy=hv选择任意形状的网格,每一个网格相当于一个蓄水容器,称为单元,网格中心的高程表示平底网格的高程,在网格中心处计算水深H;网格周边的边线相当于容器的边壁,称为通道,水体通过通道流入单元,在网格通道处计算流量Q(见附图3)。时间上采取水位-流量交错的计算方式,若在时刻T计算网格内水位,各通道流量不变;则下一个计算时刻T+dt计算网格周边通道的流量,各单元水位、水深不变(见附图4)。
用矢量表示流动速度分量,连续方程(2)为:
其中,表示单位宽度流量。由于水深定义在单元中,沿单元面积积分得:
假设水深和源汇项在同一控制体内均呈均匀分布,并且采用格林公式将方程(5)第2项沿控制体的积分转化为沿其周边通道的积分:
式中,Ai为第i个控制体的面积;Li为第i个控制体的周边通道;为控制体周边通道的外法向单位向量。对任意一个控制体,方程左边第二项可转化为:
式中,m为控制体的通道个数;Qi为第i个控制体周边通道的单宽流量;Qij为控制体i的第j个通道上的单宽流量;Lij为控制体i的第j个通道长度。将(7)代入方程(6),得到:
用矢量表示流动速度分量,运动方程(3)为:
方程(9)表示穿过单元通道的流量矢量运动平衡关系,这些流量分布在每一条通道上,沿单元通道积分可以得到流量矢量和,即:
若在求解方程(10)时,各单元水位不变,单元源汇项在计算通道流量时不直接起作用,则连续方程(5)表现为:
由于通道流量计算时刻选择各单元水位不变,则在通道法线方向为确定量,近似取:
方程(10)式化简为:
因为通道面积任意选取,方程(13)的被积函数恒为零。即:
方程(14)说明在此计算模式下,运动方程中不包含对流项的影响,或对流项影响较小,可以忽略不计。
连续方程(8)的离散形式为:
式中,Hpk为第k个管道计算断面的水深或压力水头;Apk为第k个管道计算断面的过流面积;A′pk第k个管道计算断面的等效底面积,即管段中自由水面的面积。
运动方程(14)的离散形式为:
式中,Qpk为第k个管道计算断面的单宽流量;Zk2与Zk1分别为计算断面k两侧网格中的水深;dxk为计算断面k两侧网格形心到断面中心的距离之和。公式(16)适用于水深较大的单元通道。
当单元通道水深很浅,通道两侧单元为陆地地面时,形成地面型通道,在方程(14)中流量的变化速度很小,阻力项与动力项平衡。则方程(14)为:
或写为:
式中,sign为符号函数,表示Q的正负与ΔZ的正负相同。公式(16)适用于水深较小的单元通道。离散方程形式为:
式中,为通道两侧单元的水深;Hj为通道上的平均水深;dLj取通道两侧单元形心至通道中点距离之和。
当单元通道为溢流堤坝时,采用宽顶堰流公式计算通道流量。单宽流量公式为:
式中,σs为淹没系数,m为流量系数。离散方程形式为:
当单元通道为明渠渠道时,采用明渠流公式计算通道流量。单宽流量公式为:
公式(22)与(18)有相同的表达形式。
当单元通道为闸孔时,采用闸孔出流公式计算通道流量。闸孔单宽流量公式为:
式中,e为闸孔开度,μ为闸孔自由出流流量系数,H0为闸上作用水头。离散方程形式为:
当单元通道为桥涵时,采用桥孔过流公式计算通道流量。桥孔单宽流量公式为:
式中,ε为侧收缩系数,为桥孔自由出流流量系数,hc为桥下游水深。离散方程形式为:
若单元通道为河道纵断面时,不仅有垂直于通道的流量发生,还会有沿通道的流量发生,垂直于通道的流量按堰流公式确定;沿通道的流量可按明渠流公式确定:
式中,Qh为河道流量,hh为河道水深,C为谢齐系数,R为水力半径,J为水力比降,ΔZh为河段间的水位差。
土壤含水量未达到饱和前地表产流强度计算:
在土壤含水量达到土壤持水量前,由于土壤入渗速率低于降水强度,产生的地表径流由下式计算:
式中rs为地表径流,i为降水强度,fm为最大入渗速率,f为透水层平均入渗速率。
土壤入渗速率平稳期地表产流计算:
由于入渗速率平稳,根据霍顿入渗公式原理,时段产流模型采用下式进行计算:
其中,n=2.5;Rp为时段降水量;h为时段地表径流。
在陡坡或薄水层汇流时,因为单元水量较小,坡面比降大,容易出现小水量产生大流量的虚假流量现象。解决这一问题一方面要适当减小时间步长,另一方面要从计算模式上避免出现虚假流量。在模式计算上,将方程(15)中的单元汇流量分解为流入累计量和流出累计量,(15)改写为
其中,用input下标表示流入单元累计量;用output下标表示流出单元累计量。
当为负时说明流出量过大,单元中没有这么多的水量,产生虚假流量,应成比例减小流出水量output部分水量。减小比例应为
但减小output中各量,将会影响相邻单元的入流量,可能导致其他单元产生虚假流量,需要进行重新计算来确定。若重复计算几次后使流入流出水量平衡,而影响区域不大时,则可继续下一时间步长的计算;若重复计算几次后仍不满足平衡条件,或引起大面积单元重复计算,可考虑进一步减小时间步长。
依据国家防汛抗旱总指挥部编写的“洪水风险图编制导则”中的提供的经验系数,以及实验室以往建立水力学模型的经验进行修正,最终得到不同植被类型下的糙率。
A3.8糙率(n)是水力学计算的关键参数,建议依据下垫面情况参照表1选择,风险分析成果经过实际洪水的验证后再修正。
表1不同下垫面情况下的糙率n
下垫面 | 村庄 | 树丛 | 旱田 | 水田 | 道路 | 空地 | 河道 |
糙率(n) | 0.07 | 0.065 | 0.06 | 0.05 | 0.035 | 0.035 | 0.025-0.035 |
针对数据资料提供的不同植被类型,具体赋值情况如下:
当植被类型为有林地时,n=0.065;
当植被类型为灌木林地时,n=0.063;
当植被类型为果园时,n=0.062;
当植被类型为其他草地时,n=0.061;
当植被类型为旱地时,n=0.06;
当植被类型为裸地时,n=0.035;
当植被类型为农村老宅地时,n=0.07;
当植被类型为设施农用地时,n=0.067;
当植被类型为农村道路时,n=0.035;
当植被类型为公路用地时,n=0.032;
当植被类型为采矿用地时,n=0.075。
糙率n主要影响谢才系数C,进而影响通过谢才公式计算得到的通道流速和流量。下面介绍主要计算公式和计算过程:
计算通道上的流速和流量时,采用谢才公式,即:
式中,v为断面平均流速(m/s),R为水力半径(m),A为过流断面面积,Pw为水流与固体边界接触部分的周长,称为湿周,J为水力坡度,C为谢才系数。
本次模型主要采用巴甫洛夫斯基公式来计算C值,即:
式中,n为糙率。
渗透性系数赋值根据《水利水电工程地质勘察规范》(GB50287-99)规定,岩土渗透性可按表2分级。
表2岩土渗透性分级
注:Lu为吕荣单位,是1MPa压力下,每米试段的平均压入流量,以L/min计。
本发明还提供一种基于陡坡汇流和断面流量计算的山洪泥石流预警***,包括:地形处理模块、土壤下渗模块、地表汇流模块和预警模块,其中,
所述地形处理模块对泥石流单沟及其所在流域的高精度DEM地形数据和遥感数据进行处理,依据典型地形地物特征对研究区域进行网格划分,得到各网格单元的散点高程与坡度;
所述土壤下渗模块据不同泥石流沟的地形地貌特征,计算不同地质条件下水流的下渗以及不同植被条件下地表的汇流,结合长历时降雨数据,得出各网格单元的下渗量及汇流量。
所述地表汇流模块采用基于陡坡汇流的计算方法,对研究区域内的地表汇流进行水动力过程模拟,基于所述的地形处理模块与土壤下渗模块,计算泥石流沟及其所在流域的地表汇流;
所述预警模块通过在泥石流沟的流通区范围内选取一个或多个过流断面,作为预警计算的标准参考断面,计算标准参考断面在降雨及预报过程中的流量值,并形成该断面流量随降雨过程中时间变化的流量曲线。依据村庄、公路、桥梁等设施位置,设置参考断面下游的不同预警级别的淹没范围,当标准参考断面流量计算值达到一个特定值,使得其下游的水流淹没范围达到预警淹没范围且水深达到某一预先设定值时进行相应级别的预警提醒。
优选地,所述的典型地形地物特征包括分水岭、沟道、河道等。
与现有技术相比,本发明所具有的优点和有益效果如下:
采用了基于网格陡坡汇流计算山洪泥石流的演进过程:山区地形较陡、坡度较大,将会形成流态较为复杂的坡面汇流,国内外学者虽然提出了一系列模拟陡坡汇流的计算方法,但是其精度及稳定性仍难以满足复杂的山区雨洪模拟要求。
为了全方面、多角度地研究山洪泥石流在沟域内的演进及淹没规律,本发明采用基于网格的陡坡汇流计算方法。相较于一般雨洪数学模型(尤其是商业软件)局限于处理小变坡汇流问题,本发明提出的陡坡汇流计算方法能够为山区陡坡(大变坡)汇流提供一个更精确、更可靠的计算模拟方法。
创新的预警方法模式:采用在泥石流沟内选取过流断面作为预警计算的标准参考断面,根据累积雨量、实时降雨量和预报降雨量,结合地质条件、地形条件等因素计算出在降雨(及预报)过程中给定的任意时间上该断面的流量,并以断面流量对应的下游淹没范围和水深值作为预警依据来进行分级预警阈值的设定,比以往采用固定的临界雨量值来设置分级预警值要更为合理。用这种方法来对山洪泥石流预警具有更加合理和精准的特点,具有创新性。
预警效率的提升:本发明对山区泥石流沟域,在形成沟谷型泥石流之初的洪水阶段,建立水动力数学模型进行了预警过程模拟。地质条件因素对计算的影响主要通过遥感影像、高程数据、糙率取值和土壤下渗计算来体现。该方法既可以应用于山洪预警,又可以应用于沟谷型泥石流的预警,对于及早形成山洪泥石流预警成果,提升预警效率效果明显。对保障人员安全和应急避险、降低人民生命财产损失有较大的实际意义。
修正的计算方法:(1)如图4,在时间上采用了水位-流量交错的计算方式,相较于以往单纯采用流量或者水位的方法,在算法上有所创新,采用这种计算方式模拟的结果会更加精确,更贴近实际情况;(2)采用了多种形式的通道流量计算,相比于以往只利用单一通道流量计算,多种形式更符合实际情况,有利于精确模拟不同类型的山区地形,从而提供更精确的计算模拟结果;(3)汇流的虚假流量修正,在汇流阶段,山区降雨形成的山洪,由于水量小,坡面陡,容易形成小水量,大流量的虚假情况,这里采用了减少时间步长,以及改进计算模式来修正,经过多次运算改进后,可以有效降低虚假流量,取得了比较好的效果。
最后所应说明的是:以上实施例仅以说明而非限制本发明的技术方案,尽管参照上述实施例对本发明进行了详细说明,本领域的普通技术人员应当理解:依然可以对本发明进行修改或者等同替换,而不脱离本发明的精神和范围的任何修改或局部替换,其均应涵盖在本发明的权利要求范围当中。
Claims (8)
1.一种基于陡坡汇流和断面流量计算的山洪泥石流预警方法,其特征在于,包括步骤:
S101、通过对泥石流单沟及其所在流域的高精度DEM地形数据和遥感数据进行处理,依据典型地形地物特征对研究区域进行网格划分,得到各网格单元的散点高程与坡度;
S102、根据不同泥石流沟的地形地貌特征,计算不同地质条件下水流的下渗以及不同植被条件下地表的汇流,结合长历时降雨数据,得出各网格单元的下渗量及汇流量;
S103、采用基于陡坡汇流的计算方法,对研究区域内的地表汇流进行水动力过程模拟,基于所述的步骤S101和S102,计算泥石流沟及其所在流域的地表汇流;
S104、通过在泥石流沟的流通区范围内选取一个或多个过流断面,作为预警计算的标准参考断面,计算所述标准参考断面在降雨过程中的流量值,并形成所述标准参考断面流量随降雨及预报过程中时间变化的流量曲线;依据村庄、公路和/或桥梁设施位置,设置参考断面下游的不同预警级别的淹没范围,当所述标准参考断面流量计算值达到一个特定值,使得其下游的水流淹没范围达到预警淹没范围且水深达到某一预先设定值时进行相应级别的预警提醒。
2.如权利要求1所述的山洪泥石流预警方法,其特征在于,所述的典型地形地物特征包括:分水岭、沟道和/或河道。
3.如权利要求2所述的山洪泥石流预警方法,其特征在于,通过所述步骤S101获得各网格单元的散点高程与坡度的具体步骤如下:
S01:勾勒泥石流单沟汇流流域,主要依据分水岭对泥石流单沟进行模型范围确定;
S02:在模型范围基础上,按照等尺度划分任意不规则网格,网格划分依据边界点间距,寻找最优化网格结点;
S03:模型单元网格地形根据实测资料插值获取。
4.如权利要求1所述的山洪泥石流预警方法,其特征在于,所述步骤S102中,采用干旱地区超渗产流的格林-安普顿公式,忽略湿润锋毛细管水压力计算渗流:
f=K(H/Z+1)
式中,f为下渗速率,K为饱和水力传导度,H为地面滞水深度,Z为地面以下入渗饱和水舌厚度;
饱和水力传导度K将地质条件以下渗系数的形式添加至表达式中,在模型的计算过程中,通过系数的率定来对不同的地质条件进行调整。
5.如权利要求1所述的山洪泥石流预警方法,其特征在于,所述步骤S103中,采用基于网格的水动力陡坡汇流计算方法,在溃坝大比降陡变水面水流计算模式的基础上,采用有限体积Riemann通量计算格式,建立基于网格的水动力陡坡汇流计算模型:
其中,h为水深;u为局部坐标x方向的流速;qe为源汇项;α为地形坡度;φ为降雨角度;g为重力加速度;S0为底坡比降;Sf为阻力项;σ为风应力;ρ为空气密度;p为降雨强度;vm为雨滴降落速度;d为下渗强度;i为蒸发强度。
6.如权利要求1所述的山洪泥石流预警方法,其特征在于,所述步骤S104中,所选取的该过流断面在整个降雨(及预报)过程中通过的山洪流量Q表示为以降雨量和时间为输入条件的函数Q=f(q,t;k),q为降雨量,t为时间,k为调整系数,用以调整不同的地质条件和植被条件带来的变化;对于确定的泥石流沟具有确定的地形地貌和地质条件,在降雨时间段内得出该过流断面的流量Q随时间变化的连续数值;根据断面流量Q计算出在断面下游山洪泥石流推进过程和淹没范围A=f(Q,t;z),其中Q为所选过流断面的流量值,z为该过流断面的位置;当A达到了淹没范围的分级预警值Af,且淹没水深h达到了预设值hf,则将其所对应的断面流量Qf设定为该分级预警的阈值,而Qf所对应的当时的雨量qf作为预警阀值。
7.如权利要求5所述的山洪泥石流预警方法,其特征在于,所述基于网格的水动力陡坡汇流计算方法,从三维流体力学的N-S方程出发,沿水深方向积分后获得二维浅水方程,所述二维浅水方程包括连续方程和运动方程,其中:
连续方程:
运动方程:
其中,h为水深,z为水位,z=z0+h,z0为底高程,u,v分别为x和y方向的平均流速,q为源汇项,g为重力加速度,n为糙率;
令H=h,Qx=hu,Qy=hv;选择任意形状的网格,在网格中心处计算水深H;网格周边的边线相当于容器的边壁,称为通道,水体通过通道流入单元,在网格通道处计算流量Q;时间上采取水位-流量交错的计算方式;
所述连续方程,用矢量表示流动速度分量:
其中,表示单位宽度流量;
所述连续方程也可以表示为:
其中,Ai为第i个控制体的面积,m为控制体的通道个数,Qij为控制体i的第j个通道上的单宽流量;Lij为控制体i的第j个通道长度;
所述连续方程的离散形式为:
其中,Hpk为第k个管道计算断面的水深或压力水头;Apk为第k个管道计算断面的过流面积;A′pk第k个管道计算断面的等效底面积,即管段中自由水面的面积;
所述运动方程,用矢量表示流动速度分量:
或者,所述运动方程表示为:
所述运动方程的离散形式为:
其中,Qpk为第k个管道计算断面的单宽流量;Zk2与Zk1分别为计算断面k两侧网格中的水深;dxk为计算断面k两侧网格形心到断面中心的距离之和;
所述运动方程的离散形式适用于水深较大的单元通道。
8.一种基于陡坡汇流和断面流量计算的山洪泥石流预警***,其特征在于,针对某特定泥石流沟,包括:地形处理模块、土壤下渗模块、地表汇流模块和预警模块,其中,
所述地形处理模块用于对泥石流单沟及其所在流域的高精度DEM地形数据和遥感数据进行处理,依据典型地形地物特征对研究区域进行网格划分,得到各网格单元的散点高程与坡度;
所述土壤下渗模块用于根据不同泥石流沟的地形地貌特征,计算不同地质条件下水流的下渗以及不同植被条件下地表的汇流,结合长历时降雨数据,得出各网格单元的下渗量及汇流量;
所述地表汇流模块用于基于陡坡汇流的计算方法,对研究区域内的地表汇流进行水动力过程模拟,并基于所述的地形处理模块与土壤下渗模块,计算泥石流沟及其所在流域的地表汇流;
所述预警模块用于通过在泥石流沟的流通区范围内选取一个或多个过流断面,作为预警计算的标准参考断面,计算参考断面在降雨及预报过程中的流量值,并形成该断面流量随降雨过程中时间变化的流量曲线;依据村庄、公路、桥梁等设施位置,设置参考断面下游的不同预警级别的淹没范围,当参考断面流量计算值达到一个特定值,使得其下游的水流淹没范围达到预警淹没范围且水深达到某一预先设定值时进行相应级别的预警提醒。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201811178161.4A CN109345777B (zh) | 2018-10-10 | 2018-10-10 | 基于陡坡汇流和断面流量计算的山洪泥石流预警方法和*** |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201811178161.4A CN109345777B (zh) | 2018-10-10 | 2018-10-10 | 基于陡坡汇流和断面流量计算的山洪泥石流预警方法和*** |
Publications (2)
Publication Number | Publication Date |
---|---|
CN109345777A true CN109345777A (zh) | 2019-02-15 |
CN109345777B CN109345777B (zh) | 2021-10-22 |
Family
ID=65309164
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201811178161.4A Active CN109345777B (zh) | 2018-10-10 | 2018-10-10 | 基于陡坡汇流和断面流量计算的山洪泥石流预警方法和*** |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN109345777B (zh) |
Cited By (14)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN109785586A (zh) * | 2019-03-04 | 2019-05-21 | 四川大学 | 山区小流域山洪灾害分区预警方法 |
CN110428586A (zh) * | 2019-07-08 | 2019-11-08 | 大连理工大学 | 基于前期降雨和上下游拓扑关系的农村基层洪涝预警方法 |
CN110485360A (zh) * | 2019-07-01 | 2019-11-22 | 中国水利水电科学研究院 | 一种基于虚拟水库溃决过程的山洪流量计算方法 |
CN110543692A (zh) * | 2019-08-07 | 2019-12-06 | 河海大学 | 一种基于下垫面特征的重配置汇流模拟方法 |
CN110570107A (zh) * | 2019-08-28 | 2019-12-13 | 浙江仁欣环科院有限责任公司 | 一种基于dem的山洪灾害风险评估方法 |
CN112182701A (zh) * | 2020-09-18 | 2021-01-05 | 黄河水利职业技术学院 | 一种气盾坝流量系数及流量的计算方法 |
CN113222283A (zh) * | 2021-05-31 | 2021-08-06 | 中国水利水电科学研究院 | 一种基于数字孪生的山洪预报预警方法及*** |
CN113409550A (zh) * | 2021-06-25 | 2021-09-17 | 西藏林芝市气象局 | 基于径流汇流模拟的泥石流灾害预警方法及*** |
CN114754838A (zh) * | 2022-03-31 | 2022-07-15 | 北京城市气象研究院 | 适用于自然坡面和城市内涝积水深度的预警方法 |
CN115423157A (zh) * | 2022-08-15 | 2022-12-02 | 中国水利水电科学研究院 | 一种山洪灾害风险动态预警方法及装置 |
CN116084337A (zh) * | 2022-10-28 | 2023-05-09 | 四川大学 | 用于产卵场修复中宽级配河床糙率确定的水槽试验方法 |
CN116384279A (zh) * | 2023-04-07 | 2023-07-04 | 中南林业科技大学 | 一种洪水演进过程模拟方法 |
CN117829031A (zh) * | 2024-03-01 | 2024-04-05 | 中国科学院、水利部成都山地灾害与环境研究所 | 考虑径流和泥石流相互作用的动力学模拟方法 |
CN118095139A (zh) * | 2024-04-17 | 2024-05-28 | 中国水利水电科学研究院 | 一种考虑淹没来源的洪涝识别方法 |
Citations (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP2005122562A (ja) * | 2003-10-17 | 2005-05-12 | Nec Corp | 災害時支援システムおよび災害時支援方法 |
JP2006092398A (ja) * | 2004-09-27 | 2006-04-06 | Shinko Electric Co Ltd | 災害事前検知システム |
CN106056851A (zh) * | 2016-05-13 | 2016-10-26 | 国网浙江省电力公司台州供电公司 | 电网设施暴雨预警方法 |
CN106599457A (zh) * | 2016-12-13 | 2017-04-26 | 中国水利水电科学研究院 | 一种基于Godunov格式一、二维耦合技术的山洪数值模拟方法 |
CN106599471A (zh) * | 2016-12-15 | 2017-04-26 | 中国水利水电科学研究院 | 一种分布式流域水沙过程耦合模拟方法 |
CN108446502A (zh) * | 2018-03-22 | 2018-08-24 | 中国水利水电科学研究院 | 一种利用完整二维浅水方程组获得流域单位线的方法 |
-
2018
- 2018-10-10 CN CN201811178161.4A patent/CN109345777B/zh active Active
Patent Citations (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP2005122562A (ja) * | 2003-10-17 | 2005-05-12 | Nec Corp | 災害時支援システムおよび災害時支援方法 |
JP2006092398A (ja) * | 2004-09-27 | 2006-04-06 | Shinko Electric Co Ltd | 災害事前検知システム |
CN106056851A (zh) * | 2016-05-13 | 2016-10-26 | 国网浙江省电力公司台州供电公司 | 电网设施暴雨预警方法 |
CN106599457A (zh) * | 2016-12-13 | 2017-04-26 | 中国水利水电科学研究院 | 一种基于Godunov格式一、二维耦合技术的山洪数值模拟方法 |
CN106599471A (zh) * | 2016-12-15 | 2017-04-26 | 中国水利水电科学研究院 | 一种分布式流域水沙过程耦合模拟方法 |
CN108446502A (zh) * | 2018-03-22 | 2018-08-24 | 中国水利水电科学研究院 | 一种利用完整二维浅水方程组获得流域单位线的方法 |
Cited By (22)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN109785586A (zh) * | 2019-03-04 | 2019-05-21 | 四川大学 | 山区小流域山洪灾害分区预警方法 |
CN110485360A (zh) * | 2019-07-01 | 2019-11-22 | 中国水利水电科学研究院 | 一种基于虚拟水库溃决过程的山洪流量计算方法 |
CN110428586A (zh) * | 2019-07-08 | 2019-11-08 | 大连理工大学 | 基于前期降雨和上下游拓扑关系的农村基层洪涝预警方法 |
CN110428586B (zh) * | 2019-07-08 | 2021-01-05 | 大连理工大学 | 基于前期降雨和上下游拓扑关系的农村基层洪涝预警方法 |
CN110543692A (zh) * | 2019-08-07 | 2019-12-06 | 河海大学 | 一种基于下垫面特征的重配置汇流模拟方法 |
CN110570107A (zh) * | 2019-08-28 | 2019-12-13 | 浙江仁欣环科院有限责任公司 | 一种基于dem的山洪灾害风险评估方法 |
CN112182701A (zh) * | 2020-09-18 | 2021-01-05 | 黄河水利职业技术学院 | 一种气盾坝流量系数及流量的计算方法 |
CN112182701B (zh) * | 2020-09-18 | 2024-06-11 | 黄河水利职业技术学院 | 一种气盾坝流量系数及流量的计算方法 |
CN113222283B (zh) * | 2021-05-31 | 2023-12-26 | 中国水利水电科学研究院 | 一种基于数字孪生的山洪预报预警方法及*** |
CN113222283A (zh) * | 2021-05-31 | 2021-08-06 | 中国水利水电科学研究院 | 一种基于数字孪生的山洪预报预警方法及*** |
CN113409550A (zh) * | 2021-06-25 | 2021-09-17 | 西藏林芝市气象局 | 基于径流汇流模拟的泥石流灾害预警方法及*** |
CN113409550B (zh) * | 2021-06-25 | 2022-05-27 | 西藏林芝市气象局 | 基于径流汇流模拟的泥石流灾害预警方法及*** |
CN114754838A (zh) * | 2022-03-31 | 2022-07-15 | 北京城市气象研究院 | 适用于自然坡面和城市内涝积水深度的预警方法 |
CN114754838B (zh) * | 2022-03-31 | 2022-11-11 | 北京城市气象研究院 | 适用于自然坡面和城市内涝积水深度的预警方法 |
CN115423157B (zh) * | 2022-08-15 | 2023-04-28 | 中国水利水电科学研究院 | 一种山洪灾害风险动态预警方法及装置 |
CN115423157A (zh) * | 2022-08-15 | 2022-12-02 | 中国水利水电科学研究院 | 一种山洪灾害风险动态预警方法及装置 |
CN116084337A (zh) * | 2022-10-28 | 2023-05-09 | 四川大学 | 用于产卵场修复中宽级配河床糙率确定的水槽试验方法 |
CN116384279A (zh) * | 2023-04-07 | 2023-07-04 | 中南林业科技大学 | 一种洪水演进过程模拟方法 |
CN116384279B (zh) * | 2023-04-07 | 2023-10-17 | 中南林业科技大学 | 一种洪水演进过程模拟方法 |
CN117829031A (zh) * | 2024-03-01 | 2024-04-05 | 中国科学院、水利部成都山地灾害与环境研究所 | 考虑径流和泥石流相互作用的动力学模拟方法 |
CN118095139A (zh) * | 2024-04-17 | 2024-05-28 | 中国水利水电科学研究院 | 一种考虑淹没来源的洪涝识别方法 |
CN118095139B (zh) * | 2024-04-17 | 2024-06-21 | 中国水利水电科学研究院 | 一种考虑淹没来源的洪涝识别方法 |
Also Published As
Publication number | Publication date |
---|---|
CN109345777B (zh) | 2021-10-22 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN109345777A (zh) | 基于陡坡汇流和断面流量计算的山洪泥石流预警方法和*** | |
Wang et al. | Urban flood prediction under heavy precipitation | |
Dingman | Fluvial hydraulics | |
CN106884405B (zh) | 一种无资料地区溃堤型山洪灾害分析评价方法 | |
Harrison et al. | Hydraulic and geomorphic processes in an overbank flood along a meandering, gravel‐bed river: Implications for chute formation | |
Jeevan et al. | Estimation of Soil Erosion Using RUSLE Model and GIS Techniques forConservation Planning from Kulekhani Reservoir Catchment, Nepal | |
Laks et al. | Problems with modelling water distribution in open channels with hydraulic engineering structures | |
Psilovikos et al. | An empirical model of sediment deposition processes in Lake Kerkini, Central Macedonia Greece | |
Hrissanthou et al. | Estimate of sediment inflow into Vistonis Lake, Greece | |
Chen et al. | Simulation of runoff and suspended sediment transport rate in a basin with multiple watersheds | |
Benjankar et al. | The impact of different sediment concentrations and sediment transport formulas on the simulated floodplain processes | |
Ljubenkov et al. | Numerical model of stratified flow- case study of the Neretva riverbed salination(2004) | |
Bielecka et al. | A 3D mathematical model of Vistula Lagoon hydrodynamics-general assumptions and results of preliminary calculations | |
Nayar et al. | Quantitative Morphometric analysis of Kosasthalaiyar sub basin (Chennai basin) using remote sensing (SRTM) data and GIS techniques | |
Morales-Marín et al. | Estimating sediment loadings in the South Saskatchewan river catchment | |
Millares et al. | Reservoir sedimentation and erosion processes in a snow-influenced basin in Southern Spain | |
Balouin et al. | Longshore sediment movements from tracers and models, Praia de Faro, South Portugal | |
Hrissanthou | Comparative application of two mathematical models to predict sedimentation in Yermasoyia Reservoir, Cyprus | |
Papalaskaris et al. | Calibration of a bed load transport rate model in streams of NE Greece | |
Hrissanthou et al. | Simulation modeling of runoff, groundwater flow and sediment transport into Kastoria Lake, Greece | |
Kardavani et al. | Efficiency of hydraulic models for flood zoning using GIS (Case Study: Ay-Doghmush River Basin) | |
Sobeih et al. | Evaluating Efficiency of the river nile navigational path | |
Nasution | Determining groundwater recharge from stream flow with seasonal recession method | |
Nama et al. | Review and state of the art for the hydro-morphological modeling of transboundary rivers, Tigris River as a case study | |
Sangale et al. | Variation of flow of water from Rivers Nzoia, Yala and Sio into Lake Victoria |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |