对洪水淹没分析的若干思考
关键词：洪水淹没 DEM 格网模型
一、前言
我国是一个灾害十分频繁的国家，洪涝灾害一直严重威胁着人民生命财产安全和社会的稳定与。近几十年来，自然资源的开发利用不断扩大，城乡建设飞速发展，洪水出现的频率及其造成的损失也不断的增加。因此，快速、准确、科学地模拟、预测和显示洪水淹没范围，以便发挥防洪工程效益，并以非工程措施来减轻洪水危害，对防洪减灾具有重要意义。别是对于城市和蓄滞洪区，如果能够预先获知洪水的淹没范围和水深的分布情况，对于挽救人民的生命财产和减少洪涝灾害损失都具有十分重要的价值。获取洪涝灾害范围和该范围内相对准确的水深分布对洪涝灾害评估和经济损失评估也具有重要的意义。
二十世纪九十年代以来，利用GIS技术为手段进行洪水淹没研究一直是一个研究热点，但大多数的研究仍然比较粗略，没有充分利用GIS的空间分析功能，因而没有准确地反映地形的连通性和洪水的流向等。本文以数字高程模型（DEM）为基础，将地形的连通性和TIN和任意多边形格网模型技术相结合，尝试把三维地形能够较真实地反映地形地貌要素和二维GIS技术中矢量栅格一体化的空间分析功能这两大优势运用到洪水淹没范围的模拟研究中，采用平面模拟方法模拟淹没范围，能够更准确、科学地确定洪水淹没范围和水深分布，为洪水风险图制作、防洪指挥调度和洪涝灾害的损失评估提供准确的评判依据。
二、洪水淹没分析方法及其与洪水模拟演进的比较
洪水淹没是一个很复杂的过程，受多种因素的影响，其中洪水特性和受淹区的地形地貌是影响洪水淹没的主要因素。对于一个特定防洪区域而言，洪水淹没可能有两种形式，一种是漫堤式淹没，即堤防并没有溃决，而是由于河流中洪水水位过高，超过堤防的高程，洪水漫过堤顶进入淹没区；另一种是决堤式淹没，即堤防溃决，洪水从堤防决口处流入淹没区。无论是漫堤式淹没还是决堤式淹没，洪水的淹没都是一个动态的变化的过程。
针对目前防洪减灾的应用需求，对于洪水淹没分析的要求可以概化为两种情况，一是在某一洪水水位条件下，它最终会造成多大的淹没范围和怎样的水深分布，这种情况比较适合于堤防漫顶式的淹没情况。另外一种情况是在给定某一洪量条件下，它会造成多大的淹没范围和怎样的水深分布，这种情况比较适合于溃口式淹没。对于第一种情况，需要有维持给定水位的洪水源，这在实际洪水过程中是不可能发生的，处理的办法是可以根据洪水水位的变化过程，取一个合适的洪水水位值作为淹没水位进行分析。对于第二种情况，当溃口洪水发生时，溃口大小是在变化的，导致分流比也在变化。另外一般都会采取防洪抢险措施，溃口大小与分流比在抢险过程中也在变化，洪水淹没并不能自然地发生和完成，往往有人为防洪抢险因素的作用，如溃口的堵绝，蓄滞洪区的启用等。这种情况下要直接测量溃口处进入淹没区的流量是不大可能的，因为堤防溃决的位置不确定，决口的大小也在变化，测流设施要现场架设是非常困难也是非常危险的。所以实际应用时，考虑使用河道流量的分流比来进入淹没区的洪量。
归根到底，洪水淹没的机理是由于水源区和被淹没区有通道（如溃口、开闸放水等）和存在水位差，就会产生淹没过程，洪水淹没最终的结果应该是水位达到平衡状态，这个时候的淹没区就应该是最终的淹没区。基于水动力学模型的洪水演进模型可以将这一洪水淹没过程模拟出来，即在不同时间的洪水淹没的范围，这对于分析洪水的淹没过程是非常有用的。洪水演进模型虽然能够较准确地模拟洪水演进的过程，但由于洪水演进模型建模过程复杂，建模费用高，通用性不好，一个地区的模型不能应用到另外一个地区。特别是对于江河两侧大范围的地区模型的边界很难确定。所以上述两种概化的处理方法也是常用的。三、对二维与三维分析方法的比较
洪水的淹没分析人们往往想到三维的洪水淹没分析，但实际上更有实用价值的应该是二维的淹没分析，原因有下面几点：
（1）目前所说的三维不是真正的三维，严格地说应该是2.5维，即用二维的表现设备来表现三维，这种条件下的三维往往受视角和视点的影响，视觉效果并不如人意，不能满足实际的应用需求。虽然现在市场上出现了一些真三维的观察显示设备，但其在价格，对硬件性能的要求，实用性上尚不能达到普遍使用的目的。
（2）二维来表现我们三维的现实世界是我们概化问题一贯的处理方法，特别是对于场一类的事物用二维来表现具有直观、简单、明了的优点。
（3）目前洪水模拟演进、洪涝灾害评估等实用模型大都是基于二维的，所以采用二维的洪水淹没分析能够更好地与成熟的模型结合。
基于格网模型的洪水淹没分析就是一种二维的淹没分析方法。
四、基于格网模型的淹没分析思想基于DEM的洪水淹没分析可以解决上述两种洪水最终淹没范围和水深分布的问题，但由于DEM数据量大，对于较大范围的洪水淹没分析，在目前的机硬件技术水平上还不能较快地计算出结果，这对于减灾决策实施等方面，这种计算速度是不能忍受的。格网模型的思想很早就已经提出，并且在各个领域得到广泛的应用，如有限元计算的离散单元模型，目前所能见到的较先的洪水模拟演进模型（如陆吉康教授的水动力洪水演进模型）也是一种格网化的模型。基于空间展布式社会数据库的洪涝灾害损失评估模型也是基于格网模型的思想（见李纪人等人的相关文章）。由于格网本身对模型概化的优越性，同时考虑到与洪水演进和洪涝灾害损失评估模型更好地结合，所以采用基于格网的洪水淹没分析模型是比较好的选择。
由DEM可以较方便地生成TIN模型。生成的TIN模型，其三角网格的大小分布情况反映了高程的变化情况，即在高程变化小的区域其三角网格大，在高程变化大的区域其三角网格小，这样的三角格网在洪水淹没分析方面具有以下优点：
（1）洪水淹没的特性与三角格网的这种淹没特性是一致的，即在平坦的地区淹没面积大，在陡峭的区域淹没面积小，所以采用这种格网更能模拟洪水的淹没特性。
（2）洪水的淹没边界和江河边界等都是非常不规则的，采用三角形格网模型比规则的四边形格网模型等更能够模拟这种不规则的边界。
（3）三角形格网大小疏密变化不一致，既能满足模型物理意义上的需求，也能节省计算机的存储空间，提高计算速度。
五、基于三角格网模型的淹没分析方法针对一个特定地区的洪水淹没分析，为了减少数据量和便于分析，一般根据洪水风险，预先圈定一个最大的可能淹没范围，并且将沿江两岸分成左右两半分别进行处理分析，靠江边的边界处理为淹没区的进水边界。这样处理对于减灾来说是合理的，一般在区域，沿江两岸堤防建设的洪证率是不一样的，有重点地保护一些地区和放弃一些地区，所以需要将两岸分开处理。
目前国家局能够提供七大江河周边地区1：1万的DEM数据，在实际应用中需要根据特定的区域的微地形修正该DEM数据，以保证地形数据的准确，根据实测微地形（如堤防、工程等）数据修正DEM的GRID栅格高程值。将一系列实测数据进行自动快速修正GRID的程序已经在ARCINFO DESKTOP 8.1 平台上开发完成，可以直接进行交互式使用。将修正后的DEM数据用上面提到的洪水最大可能淹没范围进行剪裁，得到的区域就是所需要进行淹没分析研究的范围。
将DEM转换为TIN模型，提取三角格网，并对每个三角格网赋高程值，高程值按三个顶点从GRID上取得的高程值的取平均求得。该三角格网就是要进行洪水淹没分析的格网模型，如下图1所示。
图1 三角形格网模型
5．1给定洪水水位（H）下的淹没分析选定洪水源入口，设定洪水水位、选出洪水水位以下的三角单元，从洪水入口单元开始进行三角格网连通性分析，能够连通的所有单元即组成淹没范围，得到连通的三角单元，对连通的每个单元计算水深W，即得到洪水淹没水深分布，如下图2所示。
图2 三角单元水深分布图
单元水深的计算公式如下式（1）。
W=H-E （1）
式中：
W—单元水深；
H—水位；
E—单元高程。
5．2洪水淹没连通区域算法对洪水淹没区域连通性的考虑，在一些淹没分析中，仅考虑高程平铺的问题，即在任何地势低洼的区域都同时进水，实际上从洪水本身淹没的角度来说这是不准确的，洪水首先是从洪水源处开始向外扩散淹没，只有水位高程达到一定程度之后，洪水才能越过某一地势较高的区域到达另一个洼地。洪水淹没的连通性算法可以从投石问路法的原理来理解。
假定有一个探险家，他带着一个高程标准（水位高程）需要将这一高程以下所有能够相互连通的区域探寻出来，假定这片区域由不同大小的格网组成，格网是由边数一样多（当然也可以不一样多，但那样会使问题更复杂）的多边形组成，这里为讨论问题的方便我们假定为四边形（其他格网单元的多边形可作类似考虑），探险家前进的方向即为投石问路的石子，探险家背着一个袋子，袋子里装着前进方向的石子。开始，探险家只有一颗石子，某一个表明能够进入的边界单元的石子，能够从这一边界单元进入的条件是，他所带的高程标准表明这一单元的高程比高程标准低。探险家投出这颗石子从这一边界单元进入，进入该单元后（对该单元做标记，表明已经走过），又得到三颗石子，即三个可能前进的方向，需要对这三颗石子检验是否可以继续用于投石问路，首先检验石子指明方向的单元是否具有已走过的标志，如果有则丢弃之，如果没有则保留，继续下一步检验。继续检验的条件是石子指明前进方向的单元高程比所带的高程标准是高还是低，如果高则该石子不合格，丢弃之，是低则合格，放入袋子中，袋中石子个数自动增加。检验完后，判断袋子中的石头个数，如果不为零，则可以继续往下探寻，再从袋子中取出一颗石子（袋中石子个数减一），继续投石问路，直到袋子中没有石子为止。这样就能遍历整个区域，找出与入口单元相连的满足高程标准的连通区域。
从问题的收敛性上来看，这种算法是完全可以收敛的，因为探险家开始的本钱只有一颗石子，每前进一步，得到得石子个数可能为0，1，2，3（别的多边形数目可能不一样，一定包括零），但他一定得消耗一颗用于探路的石子，所以如此不断探寻下去，最后石子用完，连通区域也就找出来了。
5．3给定洪量（Q）条件下的淹没分析在进行灾前预评估分析时可以根据可能发生的情况给定一个洪量，或者取洪水频率对应的流量的百分数。在灾中评估分析时Q值可以根据流量过程曲线和溃口的分流比得到，有条件的地方，可以实测，不能实测的可以根据上下游水文站点的流量差，并考虑一定区间来水的补给误差计算得到。
在上述H分析方法的基础上，通过不断给定H条件下求出对应淹没区域的容积V与Q的比较，利用二分法等逼近算法，求出与Q最接近的V，V对应的淹没范围和水深分布即为淹没分析结果
一般：
（2）
简化计算式：
（3）
式中：
V—连通淹没区水体体积；
Ai—连通淹没区单元面积，由连通性分析求解得到；
Ei—连通淹没区单元高程，由连通性分析求解得到；
m—连通淹没区单元个数，由连通性分析求解得到。
定义函数：
（4）
显然该函数为单调递减函数，函数变化趋势如图3所示：
图3 F(H)函数变化趋势图
已知 
，H0为入口单元对应的高程，要求得一个H，使得F(H) 0。为利用二分逼近算法加速求解，在程序设计时考虑变步长方法进行加速收敛过程。需要预先求得一H1使
。H1的求解可以设定一较大的增量△H循环计算，直到，（
）。再利用二分法求算
在（H0，H1）范围内趋近于零的Hq。Hq对应的淹没范围和水深分布即为给定洪量Q条件下对应淹没范围和水深
 |
图4 Hq求解示意图
六、任意多边形格网模型的洪水淹没分析方法前面谈到利用TIN模型产生的三角单元格网来进行洪水淹没分析，这样的淹没分析方法是有一些缺点的，首先由DEM产生TIN模型时对于高程有一个概化过程，即在三角单元内认为高程是均匀的，在实际处理时由三个点的高程平均取得。
将DEM转化为多边形，处理时将具有相同高程并且相邻的单元合并为一个多边形，这样可以大大减少多边形的数量，同时又能保证DEM的高程精度完全不损失。这样得到的格网模型比较三角单元格网模型，单元数量要多得多，但单元的高程精度要比三角单元高，所以三角单元的格网模型可以用于较粗精度的分析，由DEM直接转化为多边形的格网模型可以用于较高精度的分析。
任意多边形格网模型的洪水淹没分析方法与三角单元格网模型相似，也可以采用投石问路算法，但相对于三角单元格网模型在算法上略作一些技巧上的处理，因为每一个单元相邻的单元数量是不确定的，在算法上将每个单元的相邻单元编号预先生成一个序列，在对每一个单元进行投石问路时，从预先生成的序列中提取出相邻单元的编号，完成投石问路的整个算法过程，每个单元的相邻单元数量虽然是不确定的，但是有限的，所以投石问路算法一定可以收敛。下图5是任意多边形格网模型洪水淹没分析的一个例子。
图5 任意多边形格网模型洪水淹没分析结果
七、遥感监测淹没范围水深分布分析遥感监测的手段对于洪水淹没范围的确定是非常有效的，对于水深的分布情况通常是很难确定的。
由DEM生成任意多边形网格模型，该模型保证了网格单元上的高程是均等的，将遥感监测洪水淹没范围与该多边形网格模型叠加，认为淹没边界线所在的单元水深为零，淹没边界线以内的单元水深即为边界单元高程减去所在单元的高程值（这种做法是在假定淹没边界单元上的高程是相等的，实际上可能不是这样，这时可以考虑求每一个淹没边界单元相对于该单元产生的水深，然后再用距离倒数平方和加权求得该点的水深）。下图6是这种方法的一个实例，洪水遥感监测的淹没范围通过圈定一个范围来模拟，粗线为模拟的洪水遥感监测的淹没范围，淹没范围内水深分布通过颜色梯度表现。
图6 遥感监测淹没范围水深分布分析结果
[1]刘仁义、刘南，基于GIS复杂地形洪水淹没区计算方法[J]，地报，第56卷第1期，2001年1月；
[2]葛小平、许有鹏等，GIS支持下的洪水淹没范围模拟[J]，水进展，第13卷第2期，2002年7月；