一种基于GIS的大范围泥石流数值模拟及数值处理方法技术

本发明专利技术属于数据处理技术领域，公开了一种基于GIS的大范围泥石流数值模拟及数值处理方法，通过现场观测获得地形数据，建立数值标高模型；采用有限差分法求解泥石流的连续方程和运动方程，利用地球信息系统GIS中的数字高程模型DEM，自动生成地形变量，利用数字高程模型中的栅格网作为有限差分的网格；利用野外试验或实际工程验证模型和进行参数调整；使用基于DEM栅格网的有限差分公式可以实现数值解。本发明专利技术在牛顿流体运动模型基础上，推导宾汉流体的连续方程和运动方程，尤其是宾汉流体启动条件，并将其深度积分以适用于有限差分网格模拟计算。并将宾汉流体模型与GIS耦合，便于利用大范围复杂地形的计算。

A GIS based numerical simulation and processing method for large-scale debris flow

【技术实现步骤摘要】
一种基于GIS的大范围泥石流数值模拟及数值处理方法
本专利技术属于数据处理
，尤其涉及一种基于GIS的大范围泥石流数值模拟及数值处理方法。
技术介绍
目前，最接近的现有技术：基于牛顿模型和Herschel-Bulkley模型为代表的复杂黏塑性模型进行求解的连续介质计算方法。其中牛顿模型没有考虑泥石流的黏性特征，无法真实地模拟泥石流的运动特点；复杂黏塑性模型通常有过多的变量，通常只能针对某一实例进行模拟，没有普遍适用性。本专利技术利用了有限差分的宾汉模型对泥石流进行模拟，计算消耗少，并能充分考虑泥石流的黏性特征。此外对于连续介质常用于小尺度范围的局限，本专利技术与GIS结合起来，巧妙利用DEM栅格网作为有限差分网格进行计算，可以处理大范围复杂地形上的泥石流运移。泥石流是山区常见的自然灾害，每年在世界各地都造成大量的人员伤亡和财产损失。泥石流的流动距离可能很长，泛滥的范围可能很大，而且在非常平缓的斜坡上也具有较高的速度和传播能力。泥石流灾害在防灾减灾工作中已成为一项亟待解决的问题。泥石流运动数值模拟方法可以分为连续介质、离散介质和混合介质三大类：(1)连续介质计算方法，又分为基于网格和基于粒子进行离散化，进而求解计算的方法。其中基于网格的计算方法主要有：有限差分法、有限元法和有限体积法。基于网格是指将计算域划分为网格或者体元，控制方程可以转化为含有节点变量信息的代数方程组。求微分方式的控制方程为有限差分法，求积分方式的控制方程为有限元法和有限体积法。对于不同类型的泥石流常采用不同的流变本构关系。以黏性应力为主的泥石流(通常拜格诺数小于40，且雷诺数小于500)剪胀和湍动特征都不明显，可以基于宾汉模型、Cross模型、Herschel-Bulkley模型以及Voellmy模型等黏塑性本构模型进行求解；而水石流(通常拜格诺数大于450，且相对流深小于10)主要以固体颗粒间的相互作用(碰撞、摩擦)为主，剪胀效应明显，可以采用剪胀体模型进行分析。(2)离散介质计算方法，主要是离散元法和格子玻尔兹曼法。采用离散元法模拟三维颗粒流动，颗粒间接触力考虑弹黏性和阻尼作用，流体按欧拉法模拟，流固耦合考虑相间力，得到泥石流流动过程中的应力分布和能量变化，提出了DEM模拟泥石流的潜在优势。基于格子玻尔兹曼法，结合泥石流流团模型的特点建立特殊格式的格子玻尔兹曼(LB)模型，可以从平衡和非平衡态分布函数两个层次刻画泥石流复杂的流变关系和运移特征，尤其是泥石流的堆积过程。(3)混合介质计算方法，是连续介质与离散元耦合的方法。包括混合的计算流体动力学－离散元法(CFD-DEM)、大涡模拟与离散元法的耦合(LES-DEM)，还包括有限差分法、有限体积法、光滑粒子流体动力学、格子玻尔兹曼法与离散元法的耦合。耦合算法涉及的变量多且算法复杂，当模型较大时，每一时步的离散元计算结果与其他算法间交换的数据量很大，对计算机资源的消耗很大。对于小尺度的泥石流问题，常采用配合模型试验的方法，常用离散元和连续介质耦合方法进行分析。由于这些模型应用了过多的流变参数，在实际应用中受到了限制。离散元的计算速度跟单元体网格划分有很多关系，过多的接触可使计算量变大，花费大量时间，在大范围的泥石流计算中不存在优势。综上所述，现有技术存在的问题是：(1)现有连续介质方法所采用的本构模型通常无法准确描述泥石流的运动过程或者难以描述所有类型泥石流的运动过程。(2)离散元的计算速度跟单元体网格划分有很多关系，过多的接触可使计算量变大，花费大量时间，在大范围的泥石流计算中不存在优势。(3)采用配合模型试验的方法，常用离散元和连续介质耦合方法应用了过多的流变参数，这些流变参数通常难以获得，而且耦合算法复杂，计算数据量很大，在实际应用中受到了限制。解决上述技术问题的难度：上述技术问题中，复杂参数的获取通常难以实现，针对不同形态的泥石流不能普遍使用，大范围计算造成的计算成本大幅增加，这些都是制约泥石流数值计算的难点。解决上述技术问题的意义：本试验提出的一种基于GIS的大范围泥石流数值模拟及数值处理方法，基于宾汉模型的有限差分计算可以实现对各种形态泥石流的模拟；所需变量仅为地形数据、泥石流密度、粘聚力和摩擦角和速度修正系数，不受过多经验参数准确性的影响；而且作为一种连续介质的方法算法成熟，在大范围泥石流数值模拟计算中可根据需要调整网格大小，计算消耗易于接受。综上，可以完美解决上述技术问题。
技术实现思路
针对现有技术存在的问题，本专利技术提供了一种基于GIS的大范围泥石流数值模拟及数值处理方法。本专利技术是这样实现的，一种基于GIS的大范围泥石流数值模拟及数值处理方法，所述基于GIS的大范围泥石流数值模拟及数值处理方法包括以下步骤：第一步，通过现场观测获得地形数据，建立数值标高模型；第二步，采用有限差分法求解泥石流的连续方程和运动方程，利用地球信息系统GIS中的数字高程模型DEM，自动生成地形变量，利用数字高程模型中的栅格网作为有限差分的网格；利用野外试验或实际工程验证模型和进行参数调整；第三步，使用基于DEM栅格网的有限差分公式可以实现数值解。进一步，所述第一步泥石流是均匀连续、不可压缩、非定常的牛顿流体，流动遵循下列方程，其连续性方程为：运动方程为：其中p为水压力，g为重力加速度，u，v和w分别为x，y，z方向的平均流速，ρd是泥石流的平均密度；ρd＝ρsνs+ρwνw，其中ρs和ρw分别是固体颗粒和水的密度，νs和νw分别是混合物中固体颗粒和水的体积浓度；根据浅水流运动原理，u,v>>w，方程简化为：p＝ρdg(η-ηb)；其中η是流体上表面的标高，ηb是流体底面的标高；流体底面的运动边界条件写为：流体上表面的边界条件写为：得到泥石流的连续方程为：其中h为泥石流流深，和分别为x和y方向的平均流速；对方程进行深度上的积分，得到：最终形式的深度积分运动方程为：其中α为速度修正系数，提供了有关垂直速度分布与均匀性偏差的信息，β为垂直正应力的比值，H为自由面的高度，ν＝μ/ρ为运动粘度系数，μ是粘度，θx和θy分别是地形沿x和y方向的倾斜角度，为动摩擦系数，t是流动时间。进一步，所述第二步在基于栅格的DEM分析中，每个像元都有左、右、上、下以及四个对角线八个流动方向。进一步，所述第三步使用基于DEM栅格网的有限差分公式可以实现数值解，控制方程的有限差分形式为：x方向的有限差分形式为：宾汉启动条件：宾汉流体在没有达到一定切应力(τ0)之前不会由任何剪切速率；描述宾汉流体需要两个参数：屈服应力τ0和粘度μ：流体当达到：τ>τ0；时开始流动，对于层流，根据粘性流体的牛顿法则，水平方向的切应力：...

【技术保护点】
1.一种基于GIS的大范围泥石流数值模拟及数值处理方法，其特征在于，所述基于GIS的大范围泥石流数值模拟及数值处理方法包括以下步骤：/n第一步，通过现场观测获得地形数据，建立数值标高模型；/n第二步，采用有限差分法求解泥石流的连续方程和运动方程，利用地球信息系统GIS中的数字高程模型DEM，自动生成地形变量，利用数字高程模型中的栅格网作为有限差分的网格；利用野外试验或实际工程验证模型和进行参数调整；/n第三步，使用基于DEM栅格网的有限差分公式可以实现数值解。/n

【技术特征摘要】
1.一种基于GIS的大范围泥石流数值模拟及数值处理方法，其特征在于，所述基于GIS的大范围泥石流数值模拟及数值处理方法包括以下步骤：
第一步，通过现场观测获得地形数据，建立数值标高模型；
第二步，采用有限差分法求解泥石流的连续方程和运动方程，利用地球信息系统GIS中的数字高程模型DEM，自动生成地形变量，利用数字高程模型中的栅格网作为有限差分的网格；利用野外试验或实际工程验证模型和进行参数调整；
第三步，使用基于DEM栅格网的有限差分公式可以实现数值解。


2.如权利要求1所述的基于GIS的大范围泥石流数值模拟及数值处理方法，其特征在于，所述第一步泥石流是均匀连续、不可压缩、非定常的流体，流动遵循下列方程，其连续性方程为：



运动方程为：









其中p为水压力，g为重力加速度，u，v和w分别为x，y，z方向的平均流速，ρd是泥石流的平均密度；ρd＝ρsνs+ρwνw，其中ρs和ρw分别是固体颗粒和水的密度，νs和νw分别是混合物中固体颗粒和水的体积浓度；
根据浅水流运动原理，u,v>>w，方程简化为：
p＝ρdg(η-ηb)；
其中η是流体上表面的标高，ηb是流体底面的标高；
流体底面的运动边界条件写为：



流体上表面的边界条件写为：



得到泥石流的连续方程为：



其中h为泥石流流深，和分别为x和y方向的平均流速；
对方程进行深度上的积分，得到：






最终形式的深度积分运动方程为：






其中α为速度修正系数，提供了有关垂直速度分布与均匀性偏差的信息，β为垂直正应力的比值，H为自由面的高度，ν＝μ/ρ为运动粘度系数，μ是粘度，θx和θy分别是地形沿x和y方向的倾斜角度，为动摩擦系数，t是流动时间，以及不同泥石流试验或实例中密度ρ。


3.如权利要求1所述的基于GIS的大范围泥石流数值模拟及数值处理方法，其特征在于，所述第二步在基于栅格的DEM分析...

【专利技术属性】
技术研发人员：史啸，王永岩，秦楠，
申请(专利权)人：青岛科技大学，
类型：发明
国别省市：山东;37