CN110599001B - 一种泥石流流域系统的泥石流危险性评估方法及装置 - Google Patents

一种泥石流流域系统的泥石流危险性评估方法及装置 Download PDF

Info

Publication number
CN110599001B
CN110599001B CN201910777738.1A CN201910777738A CN110599001B CN 110599001 B CN110599001 B CN 110599001B CN 201910777738 A CN201910777738 A CN 201910777738A CN 110599001 B CN110599001 B CN 110599001B
Authority
CN
China
Prior art keywords
debris flow
subsystem
state variable
basin
flow basin
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.)
Active
Application number
CN201910777738.1A
Other languages
English (en)
Other versions
CN110599001A (zh
Inventor
王钧
宇岩
袁少雄
陈军
黄光庆
宫清华
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Guangzhou Institute of Geography of GDAS
Southern Marine Science and Engineering Guangdong Laboratory Guangzhou
Original Assignee
Guangzhou Institute of Geography of GDAS
Southern Marine Science and Engineering Guangdong Laboratory Guangzhou
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Guangzhou Institute of Geography of GDAS, Southern Marine Science and Engineering Guangdong Laboratory Guangzhou filed Critical Guangzhou Institute of Geography of GDAS
Priority to CN201910777738.1A priority Critical patent/CN110599001B/zh
Publication of CN110599001A publication Critical patent/CN110599001A/zh
Application granted granted Critical
Publication of CN110599001B publication Critical patent/CN110599001B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • E—FIXED CONSTRUCTIONS
    • E02—HYDRAULIC ENGINEERING; FOUNDATIONS; SOIL SHIFTING
    • E02B—HYDRAULIC ENGINEERING
    • E02B1/00—Equipment or apparatus for, or methods of, general hydraulic engineering, e.g. protection of constructions against ice-strains
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/15—Correlation function computation including computation of convolution operations
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
    • G06Q10/00—Administration; Management
    • G06Q10/06—Resources, workflows, human or project management; Enterprise or organisation planning; Enterprise or organisation modelling
    • G06Q10/063—Operations research, analysis or management
    • G06Q10/0635—Risk analysis of enterprise or organisation activities
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
    • G06Q50/00—Information and communication technology [ICT] specially adapted for implementation of business processes of specific business sectors, e.g. utilities or tourism
    • G06Q50/10—Services
    • G06Q50/26—Government or public services

Landscapes

  • Engineering & Computer Science (AREA)
  • Business, Economics & Management (AREA)
  • Physics & Mathematics (AREA)
  • Human Resources & Organizations (AREA)
  • General Physics & Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Strategic Management (AREA)
  • Tourism & Hospitality (AREA)
  • Economics (AREA)
  • Computational Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • General Business, Economics & Management (AREA)
  • Marketing (AREA)
  • Entrepreneurship & Innovation (AREA)
  • Development Economics (AREA)
  • Educational Administration (AREA)
  • Mathematical Physics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Computing Systems (AREA)
  • Civil Engineering (AREA)
  • Primary Health Care (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Health & Medical Sciences (AREA)
  • Mechanical Engineering (AREA)
  • Algebra (AREA)
  • Structural Engineering (AREA)
  • General Health & Medical Sciences (AREA)
  • Game Theory and Decision Science (AREA)
  • Operations Research (AREA)
  • Quality & Reliability (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明公开了一种泥石流流域系统的泥石流危险性评估方法及装置，其方法包括：获取待评估泥石流流域的地形子系统状态变量、固体物源子系统状态变量以及水动力子系统状态变量；利用所述地形子系统状态变量、所述固体物源子系统状态变量以及所述水动力子系统状态变量，构建泥石流流域系统状态变量；利用所构建的泥石流流域系统状态变量，构建泥石流流域系统信息熵模型；利用所构建的泥石流流域系统信息熵模型对所述待评估泥石流流域进行危险性评估。

Description

一种泥石流流域系统的泥石流危险性评估方法及装置
技术领域
本发明涉及灾害危险性评估技术领域，特别涉及一种泥石流流域系统的泥石流危险性评估方法及装置。
背景技术
泥石流是指在降水、溃坝或冰雪融化形成的地面流水作用下，在沟谷或山坡上产生的一种挟带大量泥砂、石块等固体物质的特殊洪流，其暴发突然、来势凶猛、历时短暂，具有强大的破坏力和冲击力，对山区城镇、村庄、交通、电力、通讯、水利、矿山、农业、生态等造成严重的威胁和危险。众所周知，泥石流灾害是系统问题，泥石流系统是一个开放的动态系统，系统能量分布状态和相互作用结果直接决定着流域是否为泥石流沟以及泥石流危险性大小等等。但是目前国内外对泥石流系统能量分布状态等问题未见相关研究和探讨，因此，亟需用系统科学的原理与方法去探讨泥石流灾害问题。
发明内容
根据本发明实施例提供的方案解决了现有潜在泥石流沟判识、泥石流危险性大小等防灾减灾工作中存在的关键判定评估问题。
根据本发明实施例提供的一种泥石流流域系统的泥石流危险性评估方法，包括：
获取待评估泥石流流域的地形子系统状态变量、固体物源子系统状态变量以及水动力子系统状态变量；
利用所述地形子系统状态变量、所述固体物源子系统状态变量以及所述水动力子系统状态变量，构建泥石流流域系统状态变量；
利用所构建的泥石流流域系统状态变量，构建泥石流流域系统信息熵模型；
利用所构建的泥石流流域系统信息熵模型对所述待评估泥石流流域进行危险性评估。
优选地，所述获取待评估泥石流流域的地形子系统状态变量包括：
获取待评估泥石流流域的全流域面积数据、最低高程数据和最高高程数据；
利用所获取的全流域面积数据、最低高程数据和最高高程数据，构建泥石流流域地形子系统状态变量；
其中，所述泥石流流域地形子系统状态变量的公式为：
Figure BDA0002175611950000021
其中，P1(x,t)满足：P1(x,t)≥0；
Figure BDA0002175611950000022
其中，所述P1(x,t)是指所述泥石流流域地形子系统状态变量；所述f(x,t)是指所述面积-高程曲线函数；所述
Figure BDA0002175611950000023
是指面积-高程曲线函数积分值。
优选地，所述获取待评估泥石流流域的固体物源子系统状态变量包括：
获取待评估泥石流流域的全流域面积数据、固体物源面积数据、最低高程数据和最高高程数据；
利用所获取的全流域面积数据、固体物源面积数据、最低高程数据和最高高程数据，构建泥石流流域固体物源子系统状态变量；
其中，所述泥石流流域固体物源子系统状态变量的公式为：
Figure BDA0002175611950000024
其中，P2(x,t)满足：P2(x,t)≥0；
Figure BDA0002175611950000025
其中，所述P2(x,t)是指所述泥石流流域固体物源子系统状态变量，所述g(x,t)是指所述固体物源面积-高程曲线函数；所述
Figure BDA0002175611950000031
是指固体物源面积-高程曲线函数积分值。
优选地，所述获取待评估泥石流流域的水动力子系统状态变量包括：
获取待评估泥石流流域的最低高程数据、最高高程数据、N个等高线数据以及每个等高线处的水动力数据；
利用所获取的最低高程数据、最高高程数据、N个等高线数据以及每个等高线处的水动力数据，构建泥石流流域水动力子系统状态变量；
其中，所述泥石流流域水动力子系统状态变量的公式为：
Figure BDA0002175611950000032
其中，P3(x,t)满足：P3(x,t)≥0；
Figure BDA0002175611950000033
其中，所述P3(x,t)是指泥石流流域水动力子系统状态变量，所述h(x,t)是指所述水动力-高程曲线函数；所述
Figure BDA0002175611950000034
是指水动力-高程曲线函数积分值。
优选地，所述利用所述地形子系统状态变量、所述固体物源子系统状态变量以及所述水动力子系统状态变量，构建泥石流流域系统状态变量包括：
构建泥石流流域系统的初始状态函数；
分别获取所述地形子系统状态变量P1(x1,t)对于x1的反函数、所述固体物源子系统状态变量P2(x2,t)对于x2的反函数以及所述水动力子系统状态变量P3(x3,t)对于x3的反函数；
根据所获取的所述地形子系统状态变量P1(x1,t)对于x1的反函数、所述固体物源子系统状态变量P2(x2,t)对于x2的反函数、所述水动力子系统状态变量P3(x3,t)对于x3的反函数以及所述泥石流流域系统的初始状态函数，构建泥石流流域系统状态变量；
其中，所述泥石流流域系统的初始状态函数的公式为：
Figure BDA0002175611950000041
其中，P(x1,x2,x3)≥0；∫∫∫P(x1,x2,x3)dx1dx2dx3＝1；
其中，所述a是指固体物源面积—高程曲线积分区间的第一端点值；所述b是指固体物源面积—高程曲线积分区间的第二端点值，且a<b；所述c1,c2,c3和c4分别是指x1,x2和x3的线性拟合系数。
优选地，所述泥石流流域系统状态变量的公式为：
Figure BDA0002175611950000042
其中，P(x1,x2,x3,t)满足：①P(x1,x2,x3,t)≥0；②∫∫∫P(x1,x2,x3,t)dx1dx2dx3＝1；
其中，所述an,an-1,…,a0是指拟合多项式为n时的n,n-1,…,0次项系数，所述bm,bm-1,…,b0是指拟合多项式为m时的m,m-1,…,0次项系数，所述f1是指地形子系统状态变量P1(x,t)关于自变量x1的反函数；所述f2是指固体物源子系统状态变量P2(x,t)关于自变量x2的反函数；所述f3是指水动力子系统状态变量P3(x,t)分关于自变量x3的反函数；所述t是指时间。
优选地，所述泥石流流域系统信息熵模型的公式为：
H(t)＝∫∫∫P(x1,x2,x3,t)lnP(x1,x2,x3,t)dx1dx2dx3；
其中，所述H(t)是指泥石流流域系统信息熵模型；所述x1＝f1(P1(x1,t),α)，指地形子系统状态变量对于变量x1的反函数；所述x2＝f2(P2(x2,t),an,an-1,…,a0)，指固体物源子系统对于变量x2的反函数；所述x3＝f3(P3(x3,t),bm,bm-1,…,b0)，指水动力子系统对于变量x3的反函数。
优选地，所述利用所构建的泥石流流域系统信息熵模型对所述待评估泥石流流域进行危险性评估包括：
利用所述泥石流流域信息熵模型对所述待评估泥石流流域的稳定程度和能量聚集程度进行定量评估；
根据所述稳定程度和能量聚集程度的评估结果，对所述待评估泥石流流域的危险性进行评估。
根据本发明实施例提供的一种泥石流流域系统的泥石流危险性评估装置，包括：
获取模块，用于获取待评估泥石流流域的地形子系统状态变量、固体物源子系统状态变量以及水动力子系统状态变量；
构建模块，用于利用所述地形子系统状态变量、所述固体物源子系统状态变量以及所述水动力子系统状态变量，构建泥石流流域系统状态变量，以及利用所构建的泥石流流域系统状态变量，构建泥石流流域系统信息熵模型；
危险性评估模块，用于利用所构建的泥石流流域系统信息熵模型对所述待评估泥石流流域进行危险性评估。
优选地，所述危险性评估模块具体用于利用所述泥石流流域信息熵模型对所述待评估泥石流流域的稳定程度和能量聚集程度进行定量评估，以及根据所述稳定程度和能量聚集程度的评估结果，对所述待评估泥石流流域的危险性进行评估；
其中，所述泥石流流域系统信息熵模型的公式为：
H(t)＝∫∫∫P(x1,x2,x3,t)lnP(x1,x2,x3,t)dx1dx2dx3；
其中，所述H(t)是指泥石流流域系统信息熵模型；所述x1＝f1(P1(x1,t),α)，指地形子系统状态变量对于变量x1的反函数；所述x2＝f2(P2(x2,t),an,an-1,…,a0)，指固体物源子系统对于变量x2的反函数；所述x3＝f3(P3(x3,t),bm,bm-1,…,b0)，指水动力子系统对于变量x3的反函数。
根据本发明实施例提供的方案，结合系统科学与泥石流科学，推导了泥石流流域系统及其三大子系统的状态方程的数学表达式，构建了泥石流流域系统及其信息熵模型，并提出了其物理意义。本方法思路明晰，建立的状态方程物理意义明确，可以实现泥石流流域系统状态的实时判别，可广泛应用到潜在泥石流沟判识、泥石流危险度评价、泥石流监测预警等领域。
附图说明
此处所说明的附图用来提供对本发明的进一步理解，构成本发明的一部分，本发明的示意性实施例及其说明用于理解本发明，并不构成对本发明的不当限定。在附图中：
图1是本发明实施例提供的一种泥石流流域系统的泥石流危险性评估方法流程图；
图2是本发明实施例提供的一种泥石流流域系统的泥石流危险性评估装置示意图；
图3是本发明实施例提供的泥石流流域系统状态方程及信息熵模型总体框架示意图；
图4是本发明实施例提供的泥石流流域系统组成示意图；
图5是本发明实施例提供的地形子系统的面积—高程曲线示意图；
图6是本发明实施例提供的地形子系统的泥石流流域示意图；
图7是本发明实施例提供的基于固体物源子系统的固体物源面积—高程曲线示意图；
图8是本发明实施例提供的固体物源子系统的泥石流流域示意图；
图9是本发明实施例提供的基于水动力子系统的固体物源面积—高程曲线示意图；
图10是本发明实施例提供的水动力子系统的泥石流流域示意图。
具体实施方式
以下结合附图对本发明的优选实施例进行详细说明，应当理解，以下所说明的优选实施例仅用于说明和解释本发明，并不用于限定本发明。
本发明的目的是针对现有潜在泥石流沟判识、泥石流危险性大小等防灾减灾工作中存在的关键判定评估方法问题，以系统能量为主线，基于信息熵理论和泥石流发生学原理，建立泥石流流域系统状态方程和信息熵模型，实现泥石流流域系统状态的实时判别，可广泛应用到潜在泥石流沟判识、泥石流危险度评价、泥石流监测预警等领域。本发明的技术方案是：首先，结合系统科学与泥石流科学，从小流域泥石流形成的三大要件(地形、固体物源和水动力)出发，定义泥石流流域系统；其次，基于系统能量分布状态，构建地形子系统、固体物源子系统和水动力子系统状态变量；在此基础上，基于概率论相关理论，建立泥石流流域系统状态方程和信息熵模型，整个实施方案技术路线如图3所示。
根据贝塔郎菲对系统(systems)的定义“相互作用的多元素的复合体”，本发明将泥石流流域系统(debris-flow basinsystems)定义为：相互作用的与泥石流相关的多个子系统(元素)的复合体。具体地，泥石流流域系统一般是指由地质系统、地层岩性系统、地形系统、水动力系统、固体物源系统、植被系统等相互作用的子系统构成的系统，该系统受地震、人类活动等影响。对于一个小流域系统，主要考虑泥石流形成的三大要件，则泥石流流域系由地形子系统、固体物源子系统和水动力子系统构成，而地质系统、地层岩性系统、植被系统等与大气系统、地壳系统、人文景观系统等被看作是环境背景，泥石流流域系统组成如图4所示。
图1是本发明实施例提供的一种泥石流流域系统的泥石流危险性评估方法流程图，如图1所示，包括：
步骤S1：获取待评估泥石流流域的地形子系统状态变量、固体物源子系统状态变量以及水动力子系统状态变量；
步骤S2：利用所述地形子系统状态变量、所述固体物源子系统状态变量以及所述水动力子系统状态变量，构建泥石流流域系统状态变量；
步骤S3：利用所构建的泥石流流域系统状态变量，构建泥石流流域系统信息熵模型；
步骤S4：利用所构建的泥石流流域系统信息熵模型对所述待评估泥石流流域进行危险性评估。
其中，所述获取待评估泥石流流域的地形子系统状态变量包括：获取待评估泥石流流域的全流域面积数据、最低高程数据和最高高程数据；利用所获取的全流域面积数据、最低高程数据和最高高程数据，构建泥石流流域地形子系统状态变量；其中，所述泥石流流域地形子系统状态变量的公式为：
Figure BDA0002175611950000081
其中，P1(x,t)满足：P1(x,t)≥0；
Figure BDA0002175611950000082
所述P1(x,t)是指所述泥石流流域地形子系统状态变量；所述f(x,t)是指所述面积-高程曲线函数；所述
Figure BDA0002175611950000083
是指面积-高程曲线函数积分值。
其中，所述利用所获取的全流域面积数据、最低高程数据和最高高程数据，构建泥石流流域地形子系统状态变量包括：从所述最低高程数据和所述最高高程数据中选取出N个等高线数据；利用所述全流域面积数据、所述最低高程数据、所述最高高程数据以及所述N个等高线数据，计算出每个等高线所对应的流域面积比重和流域高程比重；利用所述N个等高线的流域面积比重和流域高程比重，构建泥石流流域地形子系统状态变量；其中，N>1，且N为正整数。
具体地说，所述利用所述全流域面积数据、所述最低高程数据、所述最高高程数据以及所述N个等高线数据，计算出每个等高线所对应的流域面积比重和流域高程比重包括：根据所述最低高程数据和所述最高高程数据，计算所述最低高程数据与所述最高高程数据之间的最大流域高程差；根据所述N个等高线数据和所述最低高程数据，计算每个等高线数据的流域高程差；根据所述每个等高线数据的流域高程差和所述最大流域高程差，计算每个等高线所对应的流域高程比重。
具体地说，所述利用所述全流域面积数据、所述最低高程数据、所述最高高程数据以及所述N个等高线数据，计算出每个等高线所对应的流域面积比重和流域高程比重包括：根据所述全流域面积数据和所述N个等高线数据，计算出N+1个两两相邻等高线数据之间的流域面积数据；根据所述N个等高线数据和所述N+1个流域面积数据，计算出每个等高线数据以上的流域面积数据；根据所述每个等高线数据的以上的流域面积数据和所述全流域面积数据，计算每个等高线所对应的流域面积比重。
具体地说，所述利用所述N个等高线的流域面积比重和流域高程比重，构建泥石流流域地形子系统状态变量包括：根据所述N个等高线的流域面积比重和流域高程比重，构建泥石流流域地形子系统的面积-高程曲线函数；通过对所述泥石流流域地形子系统的面积-高程曲线函数进行积分处理，得到泥石流流域地形子系统的面积-高程曲线函数积分值；根据所述面积-高程曲线函数和所述面积-高程曲线函数积分值，构建泥石流流域地形子系统状态变量。
泥石流流域系统中地形子系统一般是指由沟坡坡度、地形坡向、集水区面积、沟谷形态等相互作用的多个子系统(要素)构成的系统。一般来讲，比较理想的地形子系统状态变量方程应该包含所有相互作用的子系统(要素)。考虑到地形在泥石流形成过程中主要贡献是为分布其上的松散固土物质提供一定势能，且已知通过流域面积大小能够反映流域集水情况，通过流域高程能够反映流域比降情况。为此，本发明在建立地形子系统状态变量方程时，主要考虑面积和高程这两个关键参数，本发明的具体方法与步骤如下：对于某一时刻t，分别以x和y为横坐标和纵坐标得到一系列点(x,y)，用曲线拟合各点，绘制面积—高程曲线，如图5所示，记为f(x,t)。
具体作法如下：如图6所示，流域最高高程为1000m，最低高程为100m，A1为100m-200m之间的流域面积；A2为200m-300m之间的流域面积；A3为300m-400m之间的流域面积；A4为400m-500m之间的流域面积；A5为500m-600m之间的流域面积；A6为600m-700m之间的流域面积；A7为700m-800m之间的流域面积；A8为800m-900m之间的流域面积；A9为900m-1000m之间的流域面积；因此，全流域面积为＝A1+A2+A3+A4+A5+A6+A7+A8+A9(km2)，其中
Figure BDA0002175611950000091
h表示流域等高线图上某条等高线与流域最低点的高差(m)，(即200m等高线与流域最低点的高差为：200-100＝100m；300m等高线与流域最低点的高差为：300-100＝200m；400m等高线与流域最低点的高差为300m；500m等高线与流域最低点的高差为400m；600m等高线与流域最低点的高差为500m；700m等高线与流域最低点的高差为600m；800m等高线与流域最低点的高差为700m；900m等高线与流域最低点的高差为800m)；H表示流域最高点与最低点的高差(m)，即1000-100＝900(m)，
Figure BDA0002175611950000101
a表示流域等高线图上相应等高线以上的面积(km2)，(即100m等高线的固体物源面积为
Figure BDA0002175611950000102
200m等高线的固体物源面积为
Figure BDA0002175611950000103
300m等高线的固体物源面积为
Figure BDA0002175611950000104
400m等高线的固体物源面积为
Figure BDA0002175611950000105
500m等高线的固体物源面积为
Figure BDA0002175611950000106
600m等高线的固体物源面积为
Figure BDA0002175611950000107
Figure BDA0002175611950000108
700m等高线的固体物源面积为
Figure BDA0002175611950000109
800m等高线的固体物源面积为
Figure BDA00021756119500001010
900m等高线的固体物源面积为
Figure BDA00021756119500001011
与x轴所围成的面积称为面积—高程曲线积分值，可以表示为如下形式：
Figure BDA00021756119500001012
其中，S表示面积—高程曲线积分值。
由图5可知，根据面积—高程曲线和面积—高程曲线积分值能够反映分布其上的松散固体物质具有势能分布状态的地形信息，也就是说，通过面积—高程曲线及其积分值可以反映地形子系统能量分布状态的地形信息。为此，下面将根据面积—高程曲线和面积—高程曲线积分值构造地形子系统能量分布状态的地形信息密度函数：
Figure BDA00021756119500001013
满足：①P1(x,t)≥0；②
Figure BDA00021756119500001014
(2)式满足密度函数性质。由于(2)式主要由面积—高程曲线和面积—高程曲线积分值构造，包含了反映势能能量分布状态的地形信息，因此，该式是称为地形子系统能量状态的地形信息密度函数，简称地形子系统密度函数，同时该式也能够表征地形子系统能量分布状态的地形信息，因此，该式称为地形子系统能量分布状态的地形信息变量，简称地形子系统状态变量。
其中，所述获取待评估泥石流流域的固体物源子系统状态变量包括：获取待评估泥石流流域的全流域面积数据、固体物源面积数据、最低高程数据和最高高程数据；利用所获取的全流域面积数据、固体物源面积数据、最低高程数据和最高高程数据，构建泥石流流域固体物源子系统状态变量；其中，所述泥石流流域固体物源子系统状态变量的公式为：
Figure BDA0002175611950000111
其中，P2(x,t)满足：P2(x,t)≥0；
Figure BDA0002175611950000112
所述P2(x,t)是指所述泥石流流域固体物源子系统状态变量，所述g(x,t)是指所述固体物源面积-高程曲线函数；所述
Figure BDA0002175611950000113
是指固体物源面积-高程曲线函数积分值。
其中，所述利用所获取的全流域面积数据、固体物源面积数据、最低高程数据和最高高程数据，构建泥石流流域固体物源子系统状态变量包括：从所述最低高程数据和所述最高高程数据中选取出N个等高线数据；利用所述全流域面积数据、固体物源面积数据、所述最低高程数据、所述最高高程数据以及所述N个等高线数据，计算出每个等高线所对应的固体物源面积比重和流域高程比重；利用所述N个等高线的固体物源面积比重和流域高程比重，构建泥石流流域固体物源子系统状态变量；其中，N>1，且N为正整数。
具体地说，所述利用所述全流域面积数据、固体物源面积数据、所述最低高程数据、所述最高高程数据以及所述N个等高线数据，计算出每个等高线所对应的固体物源面积比重和流域高程比重包括：根据所述最低高程数据和所述最高高程数据，计算所述最低高程数据与所述最高高程数据之间的最大流域高程差；根据所述N个等高线数据和所述最低高程数据，计算每个等高线数据的高程差；根据所述每个等高线数据的高程差和所述最大流域高程差，计算每个等高线所对应的流域高程比重。
具体地说，所述利用所述全流域面积数据、固体物源面积数据、所述最低高程数据、所述最高高程数据以及所述N个等高线数据，计算出每个等高线所对应的固体物源面积比重和流域高程比重包括：根据所述固体物源面积数据和所述N个等高线数据，计算出N+1个两两相邻等高线数据之间的固体物源面积数据；根据所述N个等高线数据和所述N+1个固体物源面积数据，计算出每个等高线数据以上的固体物源面积数据；根据所述每个等高线数据的以上的固体物源面积数据和所述全流域面积数据，计算每个等高线所对应的固体物源面积比重。
其中，所述利用所述N个等高线的固体物源面积比重和流域高程比重，构建泥石流流域固体物源子系统状态变量包括：根据所述N个等高线的固体物源面积比重和流域高程比重，构建泥石流流域固体物源子系统的面积-高程曲线函数；通过对所述泥石流流域固体物源子系统的面积-高程曲线函数进行积分处理，得到泥石流流域固体物源子系统的面积-高程曲线函数积分值；根据所述面积-高程曲线函数和所述面积-高程曲线函数积分值，构建泥石流流域固体物源子系统状态变量。
泥石流流域系统中的固体物源子系统一般是指由固体物源储量及分布等更小的相互作用的子系统(要素)构成的系统。泥石流固体物源类型丰富,有坍塌、崩塌、滑坡及人工破坏等多种成因类型,广泛分布于泥石流沟各区段(形成区、流通区和堆积区),其中以坍塌型和崩塌型最为普遍。目前泥石流固体物源储量的计算方法主要有：现场调查法、泥石流固体物源动储量法、可移动土体厚度法等。本发明通过固体物源面积和高程参数构造固体物源子系统状态变量方程，具体方法如下：
首先，获取研究对象(区域/单沟)固体物源面积的遥感影像，并应用ArcGIS软件对研究对象DEM数据进行遥感解译。然后，根据解译结果，获得固体物源面积数据和高程数据，并对所获得的数据进行处理，具体作法如下：如图8所示，流域最高高程为1000m，最低高程为100m，S1为100m-200m之间的固体物源面积，S1＝0；S2为200m-300m之间的固体物源面积，S2＝0；S3为300m-400m之间的固体物源面积，S3＝0；S4为400m-500m之间的固体物源面积；S5为500m-600m之间的固体物源面积；S6为600m-700m之间的固体物源面积；S7为700m-800m之间的固体物源面积；S8为800m-900m之间的固体物源面积；S9为900m-1000m之间的固体物源面积；因此，固体物源全面积为S＝S4+S5+S6+S7+S8+S9(km2)，A1为100m-200m之间的流域面积；A2为200m-300m之间的流域面积；A3为300m-400m之间的流域面积；A4为400m-500m之间的流域面积；A5为500m-600m之间的流域面积；A6为600m-700m之间的流域面积；A7为700m-800m之间的流域面积；A8为800m-900m之间的流域面积；A9为900m-1000m之间的流域面积；因此，全流域面积为A＝A1+A2+A3+A4+A5+A6+A7+A8+A9(km2)，对于某一时刻t，分别以x和y为横坐标和纵坐标得到一系列点(x,y)，其中
Figure BDA0002175611950000131
h表示流域等高线图上某条等高线与流域最低点的高差(m)，(即200m等高线与流域最低点的高差为：200-100＝100m；300m等高线与流域最低点的高差为：300-100＝200m；400m等高线与流域最低点的高差为300m；500m等高线与流域最低点的高差为400m；600m等高线与流域最低点的高差为500m；700m等高线与流域最低点的高差为600m；800m等高线与流域最低点的高差为700m；900m等高线与流域最低点的高差为800m)；H表示流域最高点与最低点的高差(m)，即1000-100＝900(m)，
Figure BDA0002175611950000132
a1表示流域等高线图上相应等高线以上的固体物源面积(km2)，(即100m等高线的固体物源面积为
Figure BDA0002175611950000133
200m等高线的固体物源面积为
Figure BDA0002175611950000134
300m等高线的固体物源面积为
Figure BDA0002175611950000135
400m等高线的固体物源面积为
Figure BDA0002175611950000136
500m等高线的固体物源面积为
Figure BDA0002175611950000137
Figure BDA0002175611950000138
600m等高线的固体物源面积为
Figure BDA0002175611950000139
Figure BDA0002175611950000141
700m等高线的固体物源面积为
Figure BDA0002175611950000142
800m等高线的固体物源面积为
Figure BDA0002175611950000143
900m等高线的固体物源面积为
Figure BDA0002175611950000144
用曲线拟合各点，这样的曲线称为固体物源面积—高程曲线，记为g(x,t)，如图7所示。
在x∈[a,b]且0≤a≤b≤1对g(x,t)进行积分，将得到的积分值称为固体物源面积—高程曲线积分值，记为M，可表示为：
Figure BDA0002175611950000145
其中，M表示固体物源面积—高程曲线积分值，a,b分别表示固体物源面积—高程曲线积分的端点。
由图7可知，固体物源面积—高程曲线积分值可以反映固体物源物质储量状态，同时对于具有相同的M值的流域，固体物源面积—高程曲线的形态有可能不一样，计算M/2对应的y值，记为K，通过K值的大小表示固体物源物质分布状态。其中K值有以下两种情况：当K≥0.5时，表示大部分固体物源物质分布于高程较大的区域(中上游)，即聚集区为中上游；当K<0.5时，表示大部分固体物源物质分布于高程较小的区域(中下游)，聚集区为中下游。
由此可见，通过固体物源面积—高程曲线积分值不但可以反映松散固体物质的储量状态，同时通过固体物源面积—高程曲线的形态可以反映松散固体物质分布状态，再结合与能量关系，可以表示流域内松散固体物质在空间任意位置具有能量状态，也就是说，固体物源面积—高程曲线和固体物源面积—高程曲线积分值是反映固体物源能量状态的固体物源储量和分布信息，因此，本文根据固体物源面积—高程曲线及其积分值，构造固体物源子系统能量状态的固体物源储量和分布信息密度函数：
Figure BDA0002175611950000146
满足：①P2(x,t)≥0；②
Figure BDA0002175611950000151
因此，(22)式满足密度函数性质。由于(22)式主要通过能够表征固体物源储量和分布信息的固体物源面积—高程曲线和固体物源面积—高程曲线积分值构造，结合与地形关系，就能够反映固体物源能量状态的固体物源储量和分布信息，因此，该式是固体物源子系统能量状态的固体物源储量和分布信息密度函数，简称固体物源子系统密度函数，又由于该式能够表征固体物源子系统能量分布状态的固体物源储量和分布信息，因此，该式又称为固体物源子系统能量状态的固体物源储量和分布信息变量，简称固体物源子系统状态变量。
其中，所述获取待评估泥石流流域的水动力子系统状态变量包括：获取待评估泥石流流域的最低高程数据、最高高程数据、N个等高线数据以及每个等高线处的水动力数据；利用所获取的最低高程数据、最高高程数据、N个等高线数据以及每个等高线处的水动力数据，构建泥石流流域水动力子系统状态变量；
其中，所述泥石流流域水动力子系统状态变量的公式为：
Figure BDA0002175611950000152
其中，P3(x,t)满足：P3(x,t)≥0；
Figure BDA0002175611950000153
所述P3(x,t)是指泥石流流域水动力子系统状态变量，所述h(x,t)是指所述水动力-高程曲线函数；所述
Figure BDA0002175611950000154
是指水动力-高程曲线函数积分值。
其中，所述利用所获取的最低高程数据、最高高程数据、N个等高线数据以及每个等高线处的水动力数据，构建泥石流流域水动力子系统状态变量包括：利用所获取的最低高程数据、最高高程数据、N个等高线数据以及每个等高线处的水动力数据，计算出每个等高线所对应的归一化后的水动力和流域高程比重；利用所述N个等高线的归一化后的水动力和流域高程比重，构建泥石流流域水动力子系统状态变量；其中，N>1，且N为正整数。
其中，所述利用所获取的最低高程数据、最高高程数据、N个等高线数据以及每个等高线处的水动力数据，计算出每个等高线所对应的归一化后的水动力和流域高程比重包括：根据所述最低高程数据和所述最高高程数据，计算所述最低高程数据与所述最高高程数据之间的最大流域高程差；根据所述N个等高线数据和所述最低高程数据，计算每个等高线数据的流域高程差；根据所述每个等高线数据的流域高程差和所述最大流域高程差，计算每个等高线所对应的流域高程比重。
其中，所述利用所获取的最低高程数据、最高高程数据、N个等高线数据以及每个等高线处的水动力数据，计算出每个等高线所对应的归一化后的水动力和流域高程比重包括：从所述N个等高线处的水动力数据中选取最大水动力数据和最小水动力数据；利用所述N个等高线处的水动力数据、所述最大水动力数据、所述最小水动力数据以及归一化无量纲公式，计算出每个等高线归一化后的水动力；其中，所述归一化无量纲公式为：
Figure BDA0002175611950000161
其中，所述xi是指i等高线归一化后的水动力；所述Wi是指i等高线处的水动力数据；所述Wmin是指最小水动力数据；所述Wmax是指最大水动力数据。
其中，所述利用所述N个等高线的归一化后的水动力和流域高程比重，构建泥石流流域水动力子系统状态变量包括：根据所述N个等高线的归一化后的水动力和流域高程比重，构建泥石流流域水动力子系统的水动力-高程曲线函数；通过对所述泥石流流域水动力子系统的水动力-高程曲线函数进行积分处理，得到泥石流流域水动力子系统的水动力-高程曲线函数积分值；根据所述水动力-高程曲线函数和所述水动力-高程曲线函数积分值，构建泥石流流域水动力子系统状态变量。
泥石流的发生与水的关系极为密切，泥石流发生的水源主要来自大气降水，其次为地下水和冰雪融水，降水对松散固体物质的稳定性有很大的影响，雨水能使松散固体物质内部的含水量发生变化，影响松散固体物质的内摩擦角和内聚力以及孔隙水压力，增加松散土体的自重，促进松散固土物质产生移动，从而为泥石流的产生和发展创造有利条件。本发明主要考虑泥石流流域系统水动力来自降水，目前针对降水型泥石流水动力的研究有：特征雨量法、土体含水量法、地表径流法等。比较理想的水动力子系统应该包含地表径流场和渗流场。本发明将泥石流流域系统的水动力子系统定义为包含渗流和地表径流流域内水循环有关的系统。考虑到山地降水量随海拔增高而增多，但存在一个最大降水量高度，超过此高度，山地降水不再随高度递增，而最大降水高度因气候干湿而异，也就是说水动力与高程密切相关，同时已知通过土体含水量和地表径流深度可以反映流域水动力情况，因此，本发明在构造水动力子系统状态变量方程时主要考虑水动力(土体含水量+地表径流深度)和高程这几个参数，具体方法如下：
如图10所示，流域最高高程为1000m，最低高程为100m，S1为100m-200m之间的固体物源面积，W100是100等位线处的水动力，W200是200等位线处的水动力，W300是300等位线处的水动力，W400是400等位线处的水动力，W500是500等位线处的水动力，W600是600等位线处的水动力，W700是700等位线处的水动力，W800是800等位线处的水动力，W900是900等位线处的水动力；，对于任意时刻t，建立水动力与高程关系曲线，分别以x和y为横坐标和纵坐标得到一系列点(x,y)，其中纵坐标为高程比重
Figure BDA0002175611950000171
h表示流域等高线图上某条等高线与流域最低点的高差(m)，H表示流域最高点与最低点的高差(m)；(即200m等高线与流域最低点的高差为：200-100＝100m；300m等高线与流域最低点的高差为：300-100＝200m；400m等高线与流域最低点的高差为300m；500m等高线与流域最低点的高差为400m；600m等高线与流域最低点的高差为500m；700m等高线与流域最低点的高差为600m；800m等高线与流域最低点的高差为700m；900m等高线与流域最低点的高差为800m)；H表示流域最高点与最低点的高差(m)，即1000-100＝900(m)，横坐标为水动力x，表示该等高线相应位置的归一化后的水动力(按照公式归一化无量纲公式为：
Figure BDA0002175611950000172
处理后的无量纲值)，用曲线拟合这些点，这样构造的曲线称为水动力—高程曲线，记为h(x,t)，对于某一时刻t＝t0，Wmin和Wmax是从所述10个等高线处的水动力数据中选取的最小水动力数据和最大水动力数据；因此，100等高线归一化后的水动力
Figure BDA0002175611950000181
200等高线归一化后的水动力
Figure BDA0002175611950000182
300等高线归一化后的水动力
Figure BDA0002175611950000183
400等高线归一化后的水动力
Figure BDA0002175611950000184
500等高线归一化后的水动力
Figure BDA0002175611950000185
600等高线归一化后的水动力
Figure BDA0002175611950000186
700等高线归一化后的水动力
Figure BDA0002175611950000187
800等高线归一化后的水动力
Figure BDA0002175611950000188
900等高线归一化后的水动力
Figure BDA0002175611950000189
1000等高线归一化后的水动力
Figure BDA00021756119500001810
用曲线拟合这些点，这样构造的曲线称为水动力—高程曲线，记为h(x,t)，对于某一时刻t＝t0，如图9所示。
其中，本发明实施例通过DEM获取流域高程和等高线数据，以及利用流域水文模型(如SHE/SWAT模型)计算流域的(水动力包括地表径流深和土体平均含水量)；本发明实施例将计算的水动力和高程数据导入到EXCEL中进行统计分析，按照等高距为100m计算流域高程比重和水动力。本发明的水动力包含地表径流和土体含水量两部分。所以水动力为两者之和。计算出每条等高线处的水动力之后，按数学中最常见的归一化公式
Figure BDA00021756119500001811
进行归一化处理，得到每条等高线的水动力，记为x。
在x∈[0,1]对h(x,t)进行积分，积分值称为水动力—高程曲线积分值，记为W，可表示为：
Figure BDA00021756119500001812
W取值有两种情况：①W→0：这时流域内几乎没有自由水，径流量→0，x→0，此时流域水动力不足。②W>0：W值越大说明水动力越充足，反之，W值越小说明水动力越不充足。由此可见，通过水动力—高程曲线积分值可以反映流域内水动力是否充足，同时通过水动力—高程曲线的形态可以反映流域内水动力的分布状态，再结合水动力与能量关系，就可以反映出流域内水动力能量状态，也就是说，水动力—高程曲线和水动力—高程曲线积分值是反映水动力的能量状态的时空分布信息，因此，下面根据水动力—高程曲线和水动力—高程曲线积分值，构造水动力子系统能量状态的水动力时空分布信息密度函数：
Figure BDA0002175611950000191
满足：①P3(x,t)≥0；②
Figure BDA0002175611950000192
因此，(32)式满足密度函数的性质。由于(32)式主要通过能够表征水动力能量状态的水动力时空分布信息的土体平均含水量和地表径流深度构造而成，结合与地形关系，就能够反映流域内水动力具有能量状态的水动力时空分布信息，因此，该式是水动力子系统能量状态的水动力时空分布信息密度函数，简称水动力子系统密度函数，又由于该式能够表征水动力子系统能量状态的水动力时空分布信息，因此，该式又称为水动力子系统能量状态的水动力时空分布信息变量，简称水动力子系统状态变量。
其中，所述利用所述地形子系统状态变量、所述固体物源子系统状态变量以及所述水动力子系统状态变量，构建泥石流流域系统状态变量包括：构建泥石流流域系统的初始状态函数；分别获取所述地形子系统状态变量P1(x1,t)对于x1的反函数、所述固体物源子系统状态变量P2(x2,t)对于x2的反函数以及所述水动力子系统状态变量P3(x3,t)对于x3的反函数；根据所获取的所述地形子系统状态变量P1(x1,t)对于x1的反函数、所述固体物源子系统状态变量P2(x2,t)对于x2的反函数、所述水动力子系统状态变量P3(x3,t)对于x3的反函数以及所述泥石流流域系统的初始状态函数，构建泥石流流域系统状态变量；
其中，所述泥石流流域系统的初始状态函数的公式为：
Figure BDA0002175611950000193
其中，P(x1,x2,x3)≥0；∫∫∫P(x1,x2,x3)dx1dx2dx3＝1；所述a是指固体物源面积—高程曲线积分区间的第一端点值；所述b是指固体物源面积—高程曲线积分区间的第二端点值，且a<b；所述c1,c2,c3和c4分别是指x1,x2和x3的线性拟合系数。
其中，所述泥石流流域系统状态变量的公式为：
Figure BDA0002175611950000201
其中，P(x1,x2,x3,t)满足：①P(x1,x2,x3,t)≥0；②∫∫∫P(x1,x2,x3,t)dx1dx2dx3＝1；所述an,an-1,…,a0是指拟合多项式为n时的n,n-1,…,0次项系数，所述bm,bm-1,…,b0是指拟合多项式为m时的m,m-1,…,0次项系数，所述f1是指地形子系统状态变量P1(x,t)关于自变量x1的反函数；所述f2是指固体物源子系统状态变量P2(x,t)关于自变量x2的反函数；所述f3是指水动力子系统状态变量P3(x,t)分关于自变量x3的反函数；所述t是指时间。
其中，所述泥石流流域系统信息熵模型的公式为：
H(t)＝∫∫∫P(x1,x2,x3,t)lnP(x1,x2,x3,t)dx1dx2dx3；其中，所述H(t)是指泥石流流域系统信息熵模型；所述x1＝f1(P1(x1,t),α)，指地形子系统状态变量对于变量x1的反函数；所述x2＝f2(P2(x2,t),an,an-1,…,a0)，指固体物源子系统对于变量x2的反函数；所述x3＝f3(P3(x3,t),bm,bm-1,…,b0)，指水动力子系统对于变量x3的反函数。
其中，所述利用所构建的泥石流流域系统信息熵模型对所述待评估泥石流流域进行危险性评估包括：利用所述泥石流流域信息熵模型对所述待评估泥石流流域的稳定程度和能量聚集程度进行定量评估；根据所述稳定程度和能量聚集程度的评估结果，对所述待评估泥石流流域的危险性进行评估。
图2是本发明实施例提供的一种泥石流流域系统的泥石流危险性评估装置示意图，如图2所示，包括：获取模块201，用于获取待评估泥石流流域的地形子系统状态变量、固体物源子系统状态变量以及水动力子系统状态变量；构建模块202，用于利用所述地形子系统状态变量、所述固体物源子系统状态变量以及所述水动力子系统状态变量，构建泥石流流域系统状态变量，以及利用所构建的泥石流流域系统状态变量，构建泥石流流域系统信息熵模型；危险性评估模块203，用于利用所构建的泥石流流域系统信息熵模型对所述待评估泥石流流域进行危险性评估。
其中，所述危险性评估模块203具体用于利用所述泥石流流域信息熵模型对所述待评估泥石流流域的稳定程度和能量聚集程度进行定量评估，以及根据所述稳定程度和能量聚集程度的评估结果，对所述待评估泥石流流域的危险性进行评估；其中，所述泥石流流域系统信息熵模型的公式为：
H(t)＝∫∫∫P(x1,x2,x3,t)lnP(x1,x2,x3,t)dx1dx2dx3；
其中，所述H(t)是指泥石流流域系统信息熵模型；所述x1＝f1(P1(x1,t),α)，指地形子系统状态变量对于变量x1的反函数；所述x2＝f2(P2(x2,t),an,an-1,…,a0)，指固体物源子系统对于变量x2的反函数；所述x3＝f3(P3(x3,t),bm,bm-1,…,b0)，指水动力子系统对于变量x3的反函数。
泥石流流域系统包括两部分：
(1)泥石流流域系统状态方程
泥石流流域系统能量分布状态的信息密度函数必定与泥石流三要素密切相关，下面将根据三个子系统构造泥石流流域系统能量分布状态的信息密度函数，具体做法如下：
已知f(x)、g(x)和h(x)对应的x分别为面积比重、固体物质面积比重和水动力，为了便于区分，分别用x1,x2,x3表示，代入f(x)、g(x)和h(x)中，将满足f(x1)＝g(x2)＝h(x3)的值构成集合，记为y，应用MATLAB的Stepwise函数作交互式逐步回归分析，得到y与x1,x2,x3多重线性关系式，用I(x1,x2,x3)表示，那么I(x1,x2,x3)数学形式可以表示为：
I(x1,x2,x3)＝c1x1+c2x2+c3x3+c4 (41)
其中，c1、c2、c3和c4为拟合函数I(x1,x2,x3)的系数，且均为常数。根据(14)式，构造函数P(x1,x2,x3)：
Figure BDA0002175611950000221
使得其满足：①P(x1,x2,x3)≥0；②∫∫∫P(x1,x2,x3)dx1dx2dx3＝1。
因此，(42)式满足密度函数性质。已知地形子系统状态变量为P1(x,t)、固体物源子系统状态变量为P2(x,t)和水动力子系统状态变量为P3(x,t)，如果将x1,x2,x3分别代入地形子系统状态变量P1(x,t)、固体物源子系统状态变量P2(x,t)和水动力子系统状态变量P3(x,t)中，那么有：
P1(x1,t)＝(1+α)(1-x1)α (43)
根据(43)式可得：
Figure BDA0002175611950000222
根据(44)式可得：
Figure BDA0002175611950000223
同理，当用多项式anxn+an-1xn-1+…+a0拟合g(x,t)时，可以得到固体物源子系统状态变量P2(x2,t)对于x2的反函数，记
x2＝f2(P2(x2,t),an,an-1,…,a0) (46)
当用多项式bmxm+bm-1xm-1+…+b0拟合h(x,t)时，可以得到水动力子系统状态变量P3(x3,t)对于x3的反函数，记
x3＝f3(P3(x3,t),bm,bm-1,…,b0) (47)
将(45)～(47)式代入(42)式，可得：
Figure BDA0002175611950000231
其中，a和b分别表示固体物源面积—高程曲线积分区间的端点值，an,an-1,…,a0表示拟合多项式为n时的n,n-1,…,0次项系数，bm,bm-1,…,b0表示拟合多项式为m时的m,m-1,…,0次项系数，f1,f2和f3分别表示地形子系统状态变量P1(x,t)、固体物源子系统状态变量P2(x,t)和水动力子系统状态变量P3(x,t)分别关于自变量x1,x2和x3的反函数，c1,c2,c3和c4分别表示拟合系数，t表示时间。
已知(48)式满足密度函数性质，并且(48)式是由三个子系统密度函数(状态变量)构成，充分包含了地形信息、固体物源信息以及水动力信息，也就是说，(48)式可以表征泥石流流域系统能量状态的信息，因此，(48)式既被称为泥石流流域系统能量状态信息密度函数，也被称为泥石流流域系统状态方程。
(2)泥石流流域信息熵模型
已知泥石流流域系统状态方程(48)式，且已知该状态方程是连续函数，结合连续熵定义，泥石流流域系统信息熵可以表示为：
H(t)＝∫∫∫P(x1,x2,x3,t)lnP(x1,x2,x3,t)dx1dx2dx3 (49)
其中，所述H(t)是指泥石流流域系统信息熵模型；所述x1＝f1(P1(x1,t),α)，指地形子系统状态变量对于变量x1的反函数；所述x2＝f2(P2(x2,t),an,an-1,…,a0)，指固体物源物源子系统对于变量x2的反函数；所述x3＝f3(P3(x3,t),bm,bm-1,…,b0)，指水动力子系统队医变量x3的反函数。
如果已知参数n,m,a,b,α,an,an-1,…,a0,bm,bm-1,…,b0，那么只需运用MATLAB中积分函数，对式(49)进行积分就可以求得泥石流流域系统信息熵H(t)的数学表达式。
下面以具体实施例对本发明实施例进行详细说明
案例区锅圈岩沟位于某市北部，距某市区约10km，是某河一级支流深溪沟左岸的一条支沟。该沟位于某国家级自然保护区内，地处某山断裂带的中南段，属于某地地震极震区(地震烈度为XI度)，流域面积为0.15km2，主沟长约580m，平均坡降270‰，流域最高海拔高程1222m，最低海拔高程943m，相对高差279m。经过调查，在地震之前，锅圈岩沟未曾发生泥石流；地震使沟内的岩土体松动，致使沟谷山体出现较大范围滑坡，形成大量的松散堆积体，为泥石流的活动提供了丰富的物质来源；同时，该地暴雨较频繁，雨量相对集中，为泥石流的产生提供了充足的水动力条件。正是这些因素的综合作用，使得锅圈岩沟在震后的每年都会暴发泥石流，如表1所示。
表1：地震前后锅圈岩沟泥石流发生频次表
Figure BDA0002175611950000241
按照本发明泥石流流域系统信息熵的方法，泥石流流域系统信息熵的计算过程如下：
锅圈岩沟地形子系统的拟合函数为P1＝(1-x)^1.27；固体物源子系统的拟合函数为P2＝-14.141x3+11.143x2-3.319x+0.9034；水动力子系统的拟合函数为P3＝-4.2628x3+7.9297x2-4.3037x+0.8169，利用matlab软件依次求得其反函数为：
x1＝1-x^(100/127)；
x2＝(((0.03536*x-0.01924)^2)^(1/2)-0.03536*x+0.01924)^(1/3)-0.009243/(((0.03536*x-0.01924)^2)^(1/2)-0.03536*x+0.01924)^(1/3)+0.2627；
x3＝(((0.03536*x-0.01924)^2)^(1/2)-0.03536*x+0.01924)^(1/3)-0.009243/(((0.03536*x-0.01924)^2)^(1/2)-0.03536*x+0.01924)^(1/3)+0.2627.
利用matlab软件Stepwise函数做回归分析，得到线性拟合函数I(x1,x2,x3)＝-4.2628x1+7.9212x2-4.3037x3+0.8168，可以得到c1、c2、c3和c4分别为-0.2628，7.9212，-4.3037，0.8186。a和b分别取0和1，将上述反函数和各参数取值代入到泥石流系统状态方程，简化得到锅圈岩沟流域系统状态方程表达式为：
P＝1.45*(((0.0354*x-0.0192)^2)^(1/2)-0.0354*x+0.0192)^(1/3)-0.0134/(((0.0354*x-0.0192)^2)^(1/2)-0.0354*x+0.0192)^(1/3)+0.105*x^(100/127)+0.603
将上式代入泥石流系统信息熵模型，并对其进行trapz梯形积分，可以得到泥石流系统信息熵值H＝0.2550.
本实施例所用的泥石流系统信息熵值与泥石流的危险性判别表如表2所示。
表2：泥石流流域系统信息熵与泥石流危险性判别表
Figure BDA0002175611950000251
可以看出，锅圈岩沟流域泥石流危险性为较高，泥石流的发育特点为：流域具备泥石流发生的地形、物源和水动力条件，能发生较大规模较大频率泥石流，潜在危害大。
根据本发明实施例提供的方案，利用所述泥石流流域信息熵模型实现泥石流流域系统状态的实时判别，可广泛应用到潜在泥石流沟判识、泥石流危险度评价、泥石流监测预警等领域。
尽管上文对本发明进行了详细说明，但是本发明不限于此，本技术领域技术人员可以根据本发明的原理进行各种修改。因此，凡按照本发明原理所作的修改，都应当理解为落入本发明的保护范围。

Claims (4)

1.一种泥石流流域系统的泥石流危险性评估方法，其特征在于，包括：
获取待评估泥石流流域的地形子系统状态变量、固体物源子系统状态变量以及水动力子系统状态变量；
利用所述地形子系统状态变量、所述固体物源子系统状态变量以及所述水动力子系统状态变量，构建泥石流流域系统状态变量，其包括：
构建泥石流流域系统的初始状态函数；
分别获取所述地形子系统状态变量P1(x1,t)对于自变量x1的反函数、所述固体物源子系统状态变量P2(x2,t)对于自变量x2的反函数以及所述水动力子系统状态变量P3(x3,t)对于自变量x3的反函数；
根据所获取的所述地形子系统状态变量P1(x1,t)对于自变量x1的反函数、所述固体物源子系统状态变量P2(x2,t)对于自变量x2的反函数、所述水动力子系统状态变量P3(x3,t)对于自变量x3的反函数以及所述泥石流流域系统的初始状态函数，构建泥石流流域系统状态变量；
其中，所述泥石流流域系统的初始状态函数的公式为：
Figure FDA0002882472820000011
其中，P(x1,x2,x3)≥0；∫∫∫P(x1,x2,x3)dx1dx2dx3＝1；
其中，所述a是指固体物源面积—高程曲线积分区间的第一端点值；所述b是指固体物源面积—高程曲线积分区间的第二端点值，且a<b；所述c1,c2,c3和c4分别是指x1,x2和x3的线性拟合系数；所述x1＝f1(P1(x1,t),α)，指地形子系统状态变量对于自变量x1的反函数；所述x2＝f2(P2(x2,t),an,an-1,…,a0)，指固体物源子系统对于自变量x2的反函数；所述x3＝f3(P3(x3,t),bm,bm-1,…,b0)，指水动力子系统对于自变量x3的反函数；所述t是指时间；I(x1,x2,x3)为3个反函数x1＝f1(P1(x1,t),α)、x2＝f2(P2(x2,t),an,an-1,…,a0)、x3＝f3(P3(x3,t),bm,bm-1,…,b0)的线性拟合函数；α为面积-高程积分曲线的拟合系数；an,an-1,…,a0是指固体物源面积-高程积分曲线拟合多项式为n时的n,n-1,…,0次项系数；所述bm,bm-1,…,b0是指水动力-高程积分曲线拟合多项式为m时的m,m-1,…,0次项系数；x1为地形子系统状态变量P1(x1,t)的自变量；x2为固体物源子系统状态变量P2(x2,t)的自变量；x3为水动力子系统状态变量P3(x3,t)的自变量；
利用所构建的泥石流流域系统状态变量，构建泥石流流域系统信息熵模型；
利用所构建的泥石流流域系统信息熵模型对所述待评估泥石流流域进行危险性评估；
其中，所述获取待评估泥石流流域的地形子系统状态变量包括：
获取待评估泥石流流域的全流域面积数据、最低高程数据和最高高程数据；
利用所获取的全流域面积数据、最低高程数据和最高高程数据，构建泥石流流域地形子系统状态变量；
其中，所述泥石流流域地形子系统状态变量的公式为：
Figure FDA0002882472820000021
其中，P1(x,t)满足：P1(x,t)≥0；
Figure FDA0002882472820000022
其中，所述P1(x,t)是指所述泥石流流域地形子系统状态变量；所述f(x,t)是指所述面积-高程曲线函数；所述
Figure FDA0002882472820000023
是指面积-高程曲线函数积分值；
其中，所述获取待评估泥石流流域的固体物源子系统状态变量包括：
获取待评估泥石流流域的全流域面积数据、固体物源面积数据、最低高程数据和最高高程数据；
利用所获取的全流域面积数据、固体物源面积数据、最低高程数据和最高高程数据，构建泥石流流域固体物源子系统状态变量；
其中，所述泥石流流域固体物源子系统状态变量的公式为：
Figure FDA0002882472820000031
其中，P2(x,t)满足：P2(x,t)≥0；
Figure FDA0002882472820000032
其中，所述P2(x,t)是指所述泥石流流域固体物源子系统状态变量，所述g(x,t)是指所述固体物源面积-高程曲线函数；所述
Figure FDA0002882472820000033
是指固体物源面积-高程曲线函数积分值；
其中，所述获取待评估泥石流流域的水动力子系统状态变量包括：
获取待评估泥石流流域的最低高程数据、最高高程数据、N个等高线数据以及每个等高线处的水动力数据；
利用所获取的最低高程数据、最高高程数据、N个等高线数据以及每个等高线处的水动力数据，构建泥石流流域水动力子系统状态变量；
其中，所述泥石流流域水动力子系统状态变量的公式为：
Figure FDA0002882472820000034
其中，P3(x,t)满足：P3(x,t)≥0；
Figure FDA0002882472820000035
其中，所述P3(x,t)是指泥石流流域水动力子系统状态变量，所述h(x,t)是指所述水动力-高程曲线函数；所述
Figure FDA0002882472820000036
是指水动力-高程曲线函数积分值；
其中，所述泥石流流域系统状态变量的公式为：
Figure FDA0002882472820000041
其中，P(x1,x2,x3,t)满足：①P(x1,x2,x3,t)≥0；②∫∫∫P(x1,x2,x3,t)dx1dx2dx3＝1；
其中，所述f1是指地形子系统状态变量P1(x1,t)关于自变量x1的反函数；所述f2是指固体物源子系统状态变量P2(x2,t)关于自变量x2的反函数；所述f3是指水动力子系统状态变量P3(x3,t)分关于自变量x3的反函数；
其中，所述泥石流流域系统信息熵模型的公式为：
H(t)＝∫∫∫P(x1,x2,x3,t)lnP(x1,x2,x3,t)dx1dx2dx3；
其中，所述H(t)是指泥石流流域系统信息熵模型。
2.根据权利要求1所述的方法，其特征在于，所述利用所构建的泥石流流域系统信息熵模型对所述待评估泥石流流域进行危险性评估包括：
利用所述泥石流流域信息熵模型对所述待评估泥石流流域的稳定程度和能量聚集程度进行定量评估；
根据所述稳定程度和能量聚集程度的评估结果，对所述待评估泥石流流域的危险性进行评估。
3.一种泥石流流域系统的泥石流危险性评估装置，其特征在于，包括：
获取模块，用于获取待评估泥石流流域的地形子系统状态变量、固体物源子系统状态变量以及水动力子系统状态变量；
构建模块，用于利用所述地形子系统状态变量、所述固体物源子系统状态变量以及所述水动力子系统状态变量，构建泥石流流域系统状态变量，以及利用所构建的泥石流流域系统状态变量，构建泥石流流域系统信息熵模型；
危险性评估模块，用于利用所构建的泥石流流域系统信息熵模型对所述待评估泥石流流域进行危险性评估；
其中，构建模块具体包括：
构建泥石流流域系统的初始状态函数；
分别获取所述地形子系统状态变量P1(x1,t)对于自变量x1的反函数、所述固体物源子系统状态变量P2(x2,t)对于自变量x2的反函数以及所述水动力子系统状态变量P3(x3,t)对于自变量x3的反函数；
根据所获取的所述地形子系统状态变量P1(x1,t)对于自变量x1的反函数、所述固体物源子系统状态变量P2(x2,t)对于自变量x2的反函数、所述水动力子系统状态变量P3(x3,t)对于自变量x3的反函数以及所述泥石流流域系统的初始状态函数，构建泥石流流域系统状态变量；
其中，所述泥石流流域系统的初始状态函数的公式为：
Figure FDA0002882472820000051
其中，P(x1,x2,x3)≥0；∫∫∫P(x1,x2,x3)dx1dx2dx3＝1；
其中，所述a是指固体物源面积—高程曲线积分区间的第一端点值；所述b是指固体物源面积—高程曲线积分区间的第二端点值，且a<b；所述c1,c2,c3和c4分别是指x1,x2和x3的线性拟合系数；所述x1＝f1(P1(x1,t),α)，指地形子系统状态变量对于自变量x1的反函数；所述x2＝f2(P2(x2,t),an,an-1,＝,a0)，指固体物源子系统对于自变量x2的反函数；所述x3＝f3(P3(x3,t),bm,bm-1,…,b0)，指水动力子系统对于自变量x3的反函数；所述t是指时间；I(x1,x2,x3)为3个反函数x1＝f1(P1(x1,t),α)、x2＝f2(P2(x2,t),an,an-1,…,a00、x3＝f3(P3(x3,t0,bm,bm-1,…,b00的线性拟合函数；α为面积-高程积分曲线的拟合系数；an,an-1,…,a0是指固体物源面积-高程积分曲线拟合多项式为n时的n,n-1,…,0次项系数；所述bm,bm-1,…,b0是指水动力-高程积分曲线拟合多项式为m时的m,m-1,…,0次项系数；x1为地形子系统状态变量P1(x1,t)的自变量；x2为固体物源子系统状态变量P2(x2,t)的自变量；x3为水动力子系统状态变量P3(x3,t)的自变量；
其中，所述获取待评估泥石流流域的地形子系统状态变量包括：
获取待评估泥石流流域的全流域面积数据、最低高程数据和最高高程数据；
利用所获取的全流域面积数据、最低高程数据和最高高程数据，构建泥石流流域地形子系统状态变量；
其中，所述泥石流流域地形子系统状态变量的公式为：
Figure FDA0002882472820000061
其中，P1(x,t)满足：P1(x,t)≥0；
Figure FDA0002882472820000062
其中，所述P1(x,t)是指所述泥石流流域地形子系统状态变量；所述f(x,t)是指所述面积-高程曲线函数；所述
Figure FDA0002882472820000063
是指面积-高程曲线函数积分值；
其中，所述获取待评估泥石流流域的固体物源子系统状态变量包括：
获取待评估泥石流流域的全流域面积数据、固体物源面积数据、最低高程数据和最高高程数据；
利用所获取的全流域面积数据、固体物源面积数据、最低高程数据和最高高程数据，构建泥石流流域固体物源子系统状态变量；
其中，所述泥石流流域固体物源子系统状态变量的公式为：
Figure FDA0002882472820000064
其中，P2(x,t)满足：P2(x,t)≥0；
Figure FDA0002882472820000065
其中，所述P2(x,t)是指所述泥石流流域固体物源子系统状态变量，所述g(x,t)是指所述固体物源面积-高程曲线函数；所述
Figure FDA0002882472820000066
是指固体物源面积-高程曲线函数积分值；
其中，所述获取待评估泥石流流域的水动力子系统状态变量包括：
获取待评估泥石流流域的最低高程数据、最高高程数据、N个等高线数据以及每个等高线处的水动力数据；
利用所获取的最低高程数据、最高高程数据、N个等高线数据以及每个等高线处的水动力数据，构建泥石流流域水动力子系统状态变量；
其中，所述泥石流流域水动力子系统状态变量的公式为：
Figure FDA0002882472820000071
其中，P3(x,t)满足：P3(x,t)≥0；
Figure FDA0002882472820000072
其中，所述P3(x,t)是指泥石流流域水动力子系统状态变量，所述h(x,t)是指所述水动力-高程曲线函数；所述
Figure FDA0002882472820000073
是指水动力-高程曲线函数积分值；
其中，所述泥石流流域系统状态变量的公式为：
Figure FDA0002882472820000074
其中，P(x1,x2,x3,t)满足：①P(x1,x2,x3,t)≥0；②∫∫∫P(x1,x2,x3,t)dx1dx2dx3＝1；
其中，所述f1是指地形子系统状态变量P1(x1,t)关于自变量x1的反函数；所述f2是指固体物源子系统状态变量P2(x2,t)关于自变量x2的反函数；所述f3是指水动力子系统状态变量P3(x3,t)分关于自变量x3的反函数；
其中，所述泥石流流域系统信息熵模型的公式为：
H(t)＝∫∫∫P(x1,x2,x3,t)lnP(x1,x2,x3,t)dx1dx2dx3；
其中，所述H(t)是指泥石流流域系统信息熵模型。
4.根据权利要求3所述的装置，其特征在于，所述危险性评估模块具体用于利用所述泥石流流域信息熵模型对所述待评估泥石流流域的稳定程度和能量聚集程度进行定量评估，以及根据所述稳定程度和能量聚集程度的评估结果，对所述待评估泥石流流域的危险性进行评估。
CN201910777738.1A 2019-08-22 2019-08-22 一种泥石流流域系统的泥石流危险性评估方法及装置 Active CN110599001B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910777738.1A CN110599001B (zh) 2019-08-22 2019-08-22 一种泥石流流域系统的泥石流危险性评估方法及装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910777738.1A CN110599001B (zh) 2019-08-22 2019-08-22 一种泥石流流域系统的泥石流危险性评估方法及装置

Publications (2)

Publication Number Publication Date
CN110599001A CN110599001A (zh) 2019-12-20
CN110599001B true CN110599001B (zh) 2021-04-30

Family

ID=68855301

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910777738.1A Active CN110599001B (zh) 2019-08-22 2019-08-22 一种泥石流流域系统的泥石流危险性评估方法及装置

Country Status (1)

Country Link
CN (1) CN110599001B (zh)

Families Citing this family (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111311879B (zh) * 2020-01-16 2021-04-20 广州地理研究所 一种泥石流的预警方法及装置
CN111241690B (zh) * 2020-01-16 2022-08-12 广东省科学院广州地理研究所 一种泥石流沟的判识方法及装置
CN114511990B (zh) * 2022-01-14 2023-08-11 中国科学院、水利部成都山地灾害与环境研究所 泥石流概率测算方法、泥石流多要素协同监测预警方法
CN116153029A (zh) * 2023-03-29 2023-05-23 广东省科学院广州地理研究所 一种基于物源信息熵和降雨量的泥石流预警方法

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104331744A (zh) * 2014-10-17 2015-02-04 中国科学院、水利部成都山地灾害与环境研究所 一种泥石流危险度评价方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7627491B2 (en) * 2003-01-07 2009-12-01 Swiss Reinsurance Company Method for evaluating flood plain risks

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104331744A (zh) * 2014-10-17 2015-02-04 中国科学院、水利部成都山地灾害与环境研究所 一种泥石流危险度评价方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
"耦合地貌信息熵和流域单元物质响应率的泥石流危险性评价";王钧;《科学技术与工程》;20170630;第17卷(第16期);第21-27页 *

Also Published As

Publication number Publication date
CN110599001A (zh) 2019-12-20

Similar Documents

Publication Publication Date Title
CN110599001B (zh) 一种泥石流流域系统的泥石流危险性评估方法及装置
Gregoretti et al. Runoff of small rocky headwater catchments: Field observations and hydrological modeling
CN110610302B (zh) 一种基于固体物源子系统的泥石流危险性评估方法及装置
Marutani et al. Influence of storm‐related sediment storage on the sediment delivery from tributary catchments in the upper Waipaoa River, New Zealand
Schmidt et al. Hydrologic conditions and terrestrial laser scanning of post-firedebris flows in the San Gabriel Mountains, CA, USA
Ghoneim et al. Assessing flash flood hazard in an arid mountainous region
Huo et al. Hydrological analysis of Loess Plateau highland control schemes in Dongzhi Plateau
Thouret et al. Volcanic hazards at Mount Semeru, East Java (Indonesia), with emphasis on lahars
CN110599002B (zh) 一种基于水动力子系统的泥石流危险性评估方法及装置
Lavigne Rate of sediment yield following small‐scale volcanic eruptions: a quantitative assessment at the Merapi and Semeru stratovolcanoes, Java, Indonesia
Al-Ansari et al. Flow of River Tigris and its effect on the bed sediment within Baghdad, Iraq
Ehiorobo et al. Monitoring of gully erosion in an urban area using geoinformation technology
CN110580573B (zh) 一种基于地形子系统的泥石流危险性评估方法及装置
CN106968217A (zh) 基于降雨频率的泥石流淤积危险范围划分方法
CN111311879B (zh) 一种泥石流的预警方法及装置
Modak et al. Gully erosion vulnerability modelling, estimation of soil loss and assessment of gully morphology: a study from cratonic part of eastern India
Singh et al. Remote sensing and GIS-based watershed prioritization for land and water conservation planning and management
Khadri et al. Detailed morphometric analysis of Man River basin in Akola and Buldhana districts of Maharashtra, India using Cartosat-1 (DEM) data and GIS techniques
Listyani et al. Groundwater occurrence prediction using regressions on morphometric variables in upstream Progo Watershed, Yogyakarta
CN109615195B (zh) 一种山区河流水文地貌分级评价方法
Burshtynska et al. Influence of geological structures on the nature of riverbed displacements for the rivers of the Dnister basin upper part
CN111241690B (zh) 一种泥石流沟的判识方法及装置
Liu et al. Geo-hazard initiation and assessment in the Three Gorges Reservoir
Dewidar Analysis of morphometric parameters using remote-sensing data and GIS techniques in the Wadi El Gemal basin, Red Sea coast, Egypt
Chaware et al. Google Earth imageries for validating groundwater recharge zone in basaltic region

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
CB02 Change of applicant information
CB02 Change of applicant information

Address after: 510070 one of the compound No. 100, Xianlie Middle Road, Guangzhou City, Guangdong Province

Applicant after: Guangzhou Institute of geography, Guangdong Academy of Sciences

Address before: 510070 one of the compound No. 100, Xianlie Middle Road, Guangzhou City, Guangdong Province

Applicant before: GUANGZHOU INSTITUTE OF GEOGRAPHY

TA01 Transfer of patent application right
TA01 Transfer of patent application right

Effective date of registration: 20210330

Address after: No.100, Xianlie Middle Road, Yuexiu District, Guangzhou City, Guangdong Province 510070

Applicant after: Guangzhou Institute of geography, Guangdong Academy of Sciences

Applicant after: Guangdong Provincial Laboratory of marine science and engineering of South China (Guangzhou)

Address before: 510070 one of the compound No. 100, Xianlie Middle Road, Guangzhou City, Guangdong Province

Applicant before: Guangzhou Institute of geography, Guangdong Academy of Sciences

GR01 Patent grant
GR01 Patent grant