引言

在全球气候变暖，极端天气事件重发、频发、突发的大背景下(王会军等，2010)，夏季突发性短时强降水(王佳丽等，2012；谌芸等，2012；孙军等，2012)引发的城市内涝等次生灾害严重威胁着城市交通网的安全运行和中心城区老旧危房居民的生命安全(Shi et al, 2009)。尤其是，地下停车场、下凹式立交桥和地铁的大量修建，快速增加了新的城市内涝风险点。例如，2004年7月10日北京中心城区因立交桥积水，主干道交通多处中断。2012年“7·21”特大自然灾害中，立交桥又一次成了北京的阿克琉斯之踵，全市共形成积水点426处，中心城区道路积水点63处，公路中断39945条次。城市内涝灾害发生的主要原因包括：不透水下垫面比例过大、排水管网标准低、城市空间立体开发和突发性极端强降水增多等(洛塔·弗克斯等，2011)。前三者是城市化发展的必然结果，很难在短时间内有所改变。因此，精细化的短时强降水预报与内涝风险预警就成为应对城市内涝灾害的重要手段和突破口，能够为城市安全运行提供科学的决策依据。

目前，城市内涝数值模拟的主要研究方法有：基于水力学、水文学和气象学等的数值模拟(解以扬等，2004；2005)，水文和气象学相结合的统计方法(马晓群等，2002)，内涝起因综合分析法(Liu et al，2002)，气象和社会经济学相结合的方法(扈海波等，2013；尤凤春等，2013)，基于AVHRR和MODIS图像的分析方法(Huang et al，2008)。随着城市的发展，城市雨水排水系统已经由原来线状结构逐步演变成网络式，因此排水系统的水文学、水利学特性逐步显现。因此，构建城市内涝数值模拟模型，并结合精细化的降水预报，实现城市内涝风险预警是可行的。

美国在城市降水径流模型及城市排水系统的数值计算模型的开发上取得显著成绩，最有代表性的是城市暴雨雨水管理模型(SWMM)，对城市排水系统有很强的模拟计算功能(董欣等，2008；Lewis，2010)。我国自20世纪80年代开始将数值模拟方法用于城市内涝研究，虽然起步较晚，但发展迅速。1998年以来，天津气象科学研究所与中国水利水电科学研究院合作，研制了天津市城区内涝仿真模型(解以扬等，2004；2005)。2000年以后，该模型在南京、深圳、西安等十几个城市得到推广应用。由于各城市在地理特征、城市规划、排水系统等方面差别很大，城市内涝数值模型具备典型的城市个例特征(Quan et al，2010)。

北京的下凹式立交桥和地铁都是全国最多的，而这两种地下空间利用方式恰恰是城市内涝最大的风险点。在天津市城区内涝仿真模型的基础上，重点处理了下凹式立交桥和地铁口(权瑞松等，2011)的地理信息和物理过程，构建了适用于北京的城市内涝数值模型。同时，利用精细化的降水监测和预报结果，与城市内涝模型衔接在不同降雨个例、不同的降雨和排水情景下，对北京的城市内涝积水进行了数值模拟。

1 资料与模型介绍 1.1 资料

(1) 北京六环内49个高质量自动气象站的降水观测数据；(2) 桥区积水深度监测数据，监测方式包括积水监测站、视频巡检和交警巡检等；(3) 不同重现期6 h降水序列(马京津等，2012)；(4) 1:10000的地理信息图层，主要包括高程、河道、排水工程、建筑和道路等。

1.2 北京城市内涝数值模型

北京城市内涝数值模型(Beijing Urban Waterlogging numerical model，BUW)，根据北京复杂地形和大城市特点，以地理信息系统为支撑将各类空间信息剖分为6458个网格及相应的通道(图 1)。模型依托北京市气象局的精细化降水观测和预报，围绕城市地表、河道沟渠、排水管网等城市主要水文水动力学物理过程，模拟积水深度变化情况。

1.2.1 地表水文过程基本原理

有限元法可用来求解数学中具有初边值的偏微分方程，在地表水文模拟中得到了较广泛的应用(沈冰，1986)。BUW模型应用有限体积法的思想，根据地形、土地使用、路网和水系等地理信息数据，将研究区划分成若干个无结构不规则的网格，每个网格的水文特性是均一和连续的。模型以城市地表和明渠河道水流运动为主要模拟对象。地表水文过程的基本控制方程以平面二维非恒定流方程(解以扬等，2004)为骨架。

二维非恒定流基本方程如下：

$ \frac{\partial \mathit{H}}{\partial \mathit{t}}+\frac{\partial \mathit{M}}{\partial \mathit{x}}+\frac{\partial \mathit{N}}{\partial \mathit{y}}=\mathit{q} $

(1)

连续方程：

$ \begin{align} & \frac{\partial \mathit{M}}{\partial \mathit{t}}+\frac{\partial \left(\mathit{uM} \right)}{\partial \mathit{x}}+\frac{\partial \left(\mathit{vM} \right)}{\partial \mathit{y}}+\mathit{gH}\frac{\partial \mathit{Z}}{\partial \mathit{x}}+ \\ & \ \ \ \ \ \ \ \ \ \ \ \ \ \ g\frac{{{\mathit{n}}^{2}}\mathit{u}\sqrt{{{u}^{2}}+{{v}^{2}}}}{{{\mathit{H}}^{1/3}}}=0 \\ \end{align} $

(2)

动量方程：

$ \begin{align} & \frac{\partial \mathit{N}}{\partial \mathit{t}}+\frac{\partial \left(\mathit{uN} \right)}{\partial \mathit{x}}+\frac{\partial \left(\mathit{vN} \right)}{\partial \mathit{y}}+\mathit{gH}\frac{\partial \mathit{Z}}{\partial \mathit{y}}+ \\ & \ \ \ \ \ \ \ \ \ \ \ \ \ \ g\frac{{{\mathit{n}}^{2}}\mathit{v}\sqrt{{{u}^{2}}+{{v}^{2}}}}{{{\mathit{H}}^{1/3}}}=0 \\ \end{align} $

(3)

式中，H为水深；Z为水位; q为源汇项，也就是径流量；M和N分别为x和y方向上的单宽流量，且M=Hu，N=Hv；u和v分别为流速在x，y方向上的分量；n为糙率；g为重力加速度。

1.2.2 地下水文过程基本原理

排水管网及配套的泵、闸、管网出口等组成了自成体系的城市地下水文系统，对城市内涝有很大的影响。在降雨过程中，地面积水在管道内汇集后沿各自管道系统汇合至出口处，再经过出口处的闸门、泵站或淹没出流管道排到河道中，形成“雨水—地面积水—管道汇水—管道排水—河道汇水”的模拟过程。由地表进入管网的单元网格的排水能力根据雨强、积水深度以及网格的最大排水能力来确定：当降雨强度和积水深度小于网格的最大排水能力时，降水和积水全部排入地下管网；当最大排水强度小于积水深度和降水强度时，超出部分则转化为地表径流量。

为节约计算资源，模型将管道的属性概化为：(1) 经过网格中心，(2) 经过网格周边通道的中点与相邻网格相连(图 2)。获取全部城市排水管网有极大的难度，考虑到排水管网主要分布在道路下面，部分管网按道路长度概化长度，按道路等级概化管径。排水管网由于管道底部高程和管径不同，存在明渠流动和有压流动两种形式，模型中必须考虑在同一管网内明渠流动和有压流动的交替计算(解以扬等，2005)。

一维明渠流动：

$ 连续方程:\frac{\partial \mathit{y}}{\partial \mathit{t}}+\mathit{U}\frac{\partial y}{\partial \mathit{l}}+\mathit{y}\frac{\partial \mathit{U}}{\partial \mathit{l}}=0 $

(3)

$ 动量方程:g\frac{\partial \mathit{y}}{\partial \mathit{l}}+\frac{\partial U}{\partial \mathit{t}}+\mathit{U}\frac{\partial \mathit{U}}{\partial \mathit{l}}=-g\mathit{A}{{\mathit{S}}_{\mathit{f}}} $

(4)

一维有压流动：

$ 连续方程:\frac{\partial \mathit{H}}{\partial \mathit{t}}+U\frac{\partial H}{\partial \mathit{l}}+\frac{{{a}^{2}}}{g}\frac{\partial \mathit{U}}{\partial \mathit{l}}=0 $

(5)

$ 动量方程:g\frac{\partial \mathit{H}}{\partial \mathit{t}}+\frac{\partial U}{\partial \mathit{t}}+U\frac{\partial \mathit{U}}{\partial \mathit{l}}=-g\mathit{A}{{\mathit{S}}_{\mathit{f}}} $

(6)

式中，y为明渠中的水深，H为满流中的测压管水头，a为水击波速，l为A作为计算断面的过水面积，U为断面的单宽流量，S_f为摩阻坡降。

1.2.3 地势低洼处的涌出水处理

下凹式立交桥和居民小区低洼处等较四周高度明显偏低，极易产生严重积水。除了大范围的地表客水汇流外，地下管道涌出水也是重要来源。河道水位过高或下游排水堵塞时，雨水管道满溢，本应沿管道排走的水流在地势低洼处喷涌而出，会导致极为严重的积水，整个致灾过程往往在十几分钟内完成，让人措手不及。

北京有数十座下凹式立交桥，桥区跨度比较大，且比周边地势偏低数米以上，在两侧形成人工堤坝，是天然的“聚水盆”，极易在短时间内形成面积大、深度深的积水区。因此，在建立下凹式立交桥所处网格时，需要将网格高程从桥面高度降低到桥底高度，并扩大有效客水面积，使地表汇流和管道涌出水过程的物理描述更符合真实情况。

在地势低洼的地方，当流入管道的水量超过管道容积时，水体沿水井涌出地面，首先计算单元管道内的水体积

$ {{V}^{T+2{\rm{d}}t}}={{V}^{T}}+2\rm{d}\mathit{t}\times (\sum\limits_{\rm{i=1}}^{\mathit{N}}{\mathit{Q}_{\mathit{pi}}^{\mathit{T}\rm{+d}\mathit{t}}\rm{+}{{\mathit{Q}}_{\rm{1}}}}) $

(7)

式中dt为计算时间间隔，V^T和V^T+2dt分别为T和T+2dt时刻的单元管道内的水体积，N为单元网格内管网出口数，Q_pi^T+dt为T+dt时刻管网各出口的排水量，Q₁为降水后由该网格单元地面进入管网的排水量。

设网格单元内管道总体积为V_m，当V^T+2dt≤V_m，网格单元内的下水管道处于正常的泻水状态，不会上涌；当V^T+2dt＞V_m，下水管道向网格单元涌水，上涌的水体积dV为：

$ \rm{d}\mathit{V}={{\mathit{V}}^{\mathit{T}\rm{+d}\mathit{t}}}\rm{-}{{\mathit{V}}_{\mathit{m}}} $

(8)

2 BUW模型的性能检验 2.1 “7·21”降雨和积水实况

2012年7月21—22日，北京市出现历史罕见强降雨过程，为1951年以来最强的一次全市性特大暴雨过程，此次暴雨过程具有历时短、雨势强、范围广、山区雨量大等特点。降雨主要集中在21日10时至22日06时，在近20 h内全市平均降雨量170 mm，最大降雨量541 mm。就城区而言，平均雨量达215 mm，最大降雨量328.0 mm(模式口)。由图 3可见，降雨主要分为两个阶段，17时前以暖区降水为主(R1)，之后锋面系统移入北京，降水以锋面降水为主(R2)(孙军等，2012)。14—15时，城区的平均雨强超过20 mm·h^-1。第二阶段雨强显著增大，19—20时超过40 mm·h^-1。城区最大雨强出现在丰台体育中心站，为85.9 mm·h^-1。根据《北京市气象灾害预警信号与防御指南(2013年5月)》中1 h降雨量阈值, 城区93.9%的站点最大雨强达到橙色以上预警级别，30.6%的站点达到红色预警级别(图 4)。

全市共形成积水点426处，中心城区道路积水点63处(图 5)。从总体看城区南部积水比北部严重，东部比西部严重。积水点主要分布在ABCD四个区域，其中CD区域积水最严重。按照北京市的规定，下凹式立交桥和易积水路段积水深度达27 cm时，立即采取封路措施和车辆疏导工作。图 5中红点的积水深度均超过30 cm，尤其是莲花桥、广渠门桥、双营桥和肖村桥的积水深度达到或超过2 m，D区中广渠门桥因积水过深造成1人遇难。五环外，房山是重灾区，有38人遇难，其中青龙湖镇5人，河北镇4人，京港澳高速附近多人，多数都是城市内涝灾害的遇难者。

2.2 BUW模拟的积水情况

以北京六环内47个自动气象站的降水观测数据为驱动条件，BUW从7月21日10时开始模拟到22日01时，共16 h。最大积水深度如图 6所示，四个等级根据不同车型的涉水深度划定。从空间分布上看，能够模拟出五环内的严重内涝，也能表现出南多北少，东多西少的分布特征。图 5中A、B、C三个区域的积水情况能够得到较好的模拟，D区尽管也模拟出不少的积水，但和实况相比明显偏弱。同时，模型对五环外，尤其是房山部分地区的内涝积水也有很好的模拟能力。房山青龙湖和河北镇附近多个网格模拟出超过1 m的积水，部分网格超过2 m，京港澳高速周围也模拟出严重的积水，这和灾情分布是相符的。

下凹桥区极易成为周围降雨的客水汇流区或产生地下涌出水，是城市内涝的高风险点。表 1对北京市的21处桥区进行了更为细致的检验分析，表 1中的桥区主要集中在五环内，少数位于五环外。大部分桥区的模拟和实测值比较吻合，如五路桥、复兴门桥、安华桥、丽泽桥、六里桥、安贞桥、木樨园桥、赵公口桥和菜户营桥等。这些桥区的积水一般在1 m以下，模型对这种量级的积水有不错的模拟能力。莲花桥实测积水超过2 m，模拟值1.77 m，肖村桥积水2 m，模拟值1.23 m，其他积水达到2 m的桥区，模拟值都没有超过1 m。这说明，模型对2 m左右的深积水模拟能力还有限，这可能是因为模型中对下凹式立交桥两侧的陡直壁和大范围客水考虑不足的原因。

根据4个涉水等级，从表 1中各选1个代表桥区，检验模型对积水时间、积水过程和退水时间的模拟能力(图 7)。莲花桥从降水第一阶段(R1) 开始出现积水，一直在增加，到第二阶段(R2) 的19—20时达到最大积水，这也是雨强最大的时段，之后积水开始退去。同样，复兴门桥、六里桥和安贞桥的积水过程也和降雨的节奏是吻合的，能够模拟出积水和退水的时间。值得注意的是，因缺乏动态的积水深度数据，仅依靠雨强对积水时间、积水过程和退水时间模拟能力的检验是初步的，需要进一步加强。

除对“7·21”这种“长时间大雨强型”城市内涝进行模拟外，还对“长时间小雨强型”(2013年8月11日)，“短时间大雨强型”(2011年6月23日)、“短时间小雨强型”(2013年6月17日)等三种类型的城市内涝进行了模拟。模型对四种类型城市内涝的空间分布、积水时间、积水过程和退水时间，尤其是最大积水深度都有不错的模拟能力。

3 不同重现期强降雨的城市积水模拟

采用北京市观象台逐分钟降雨资料，应用广义偏态分布方法得出360 min历时不同重现期降雨量(马京津等，2013)，选取其中2年(62.4 mm)、10年(98.82 mm)、50年(157.1 mm)和100年(208.2 mm)重现期，构建不同降雨情景。每个网格都给予相同的降雨序列，前6 h为将360 min历时降雨量平均到6 h，后5 h无降雨量。如2年重现期降雨情景为降雨持续6 h，每小时降雨10.4 mm，无降水的5 h用于地表和地下水文过程的调整和退水模拟。

如图 8，2年重现期情景，能够代表北京城区现行的雨水排水管网的设计标准。此情景下，北京的城市内涝主要发生在四环之内，南多北少，积水深度一般在25 cm以下，部分网格超过25 cm。积水网格一般没有连成片，以孤立的积水点为主。六环外左上角也有一些积水出现，主要受山区地形的影响。10年重现期情景下，四环内积水深度明显加深，出现大片超过25 cm的积水区域，东部和北部的四环和五环之间开始出现明显积水。城市内涝在空间上表现出南多北少、东多西少的分布特点，南部开始出现成片的积水区域。50年重现期情景下，五环内城市内涝继续加重，积水深度以25~80 cm为主，五环到六环之间也开始出现零星的积水点。值得注意的是，房山青龙湖附近开始出现超过50 cm的积水。100年重现期，即雨强34.7 mm·h^-1的降雨情景下，整个五环内都出现严重的城市内涝，部分地区积水超过80 cm，南部，尤其是西南部出现大片超过50 cm的积水区域。六环外左上角出现也有一片超过80 cm的积水区域，可能是爆发山洪的结果。房山青龙湖附近也出现一些超过80 cm的积水区域，这和“7·21”比较类似。

4 不同排水能力下的城市积水模拟

“7·21”特大暴雨引发严重城市内涝灾害的原因，除雨强大时间长外，还有排水系统设计不合理等因素，比如管道直径太小等。在BUW数值模型中将北京市的管网直径人为拓宽，定量评估排水管网在城市积水当中的作用。1.2倍管网试验(1.2Pipe)中，除个别网格外，城市内涝积水并没有明显变化，也就是说直径拓宽20%并不能明显改善城市排水。1.6倍管网试验(1.6Pipe)中，四环到五环之间的积水明显减弱，四环以内的大部分积水减弱。但是东北二环附近却有新的超过0.8 m的积水出现，这可能是因为管网中的水太多，发生地下涌出水的原因。当管网拓宽到2倍时(2.0Pipe)，大部分城市积水都已经消失，仅在二环到四环之间还有一些较浅的积水，东北二环也没有出现地下涌出水。当管网拓宽到2.4倍时(2.4Pipe)，六环内大部分积水均已消失，即使面对“7·21”等级的特大暴雨也可以发挥排水系统应有的作用。值得注意的西北六环外和青龙湖附近的积水并没有减弱，它们产生的原因是地形地貌，而不是排水不畅。

5 结论与讨论

北京城市内涝数值模型(BUW)根据北京复杂地形和大城市的特点，围绕城市地表、河道沟渠、排水管网等城市主要水文水动力学物理过程，模拟积水深度变化情况。以精细化的降水监测为驱动条件，BUW可以较好地模拟出“7·21”城市内涝积水的空间分布，对重点桥区的积水深度，积水过程的模拟也比较贴近实际，表明模型具有较好的模拟性能。

2年重现期情景下，北京四环内就会产生一定的积水，以孤立的积水点为主。10和50年重现期下，积水的深度和范围都有所加重。100年重现期，整个五环内都出现严重的城市积水，南部出现大片超过50 cm的积水区域。六环外左上角出现也有一片超过80 cm的积水区域，可能是爆发山洪的结果。房山青龙湖附近也出现一些超过80 cm的积水区域，和“7·21”比较类似。

排水管网直径拓宽20%并不能明显改善城市排水能力。拓宽60%时，四环到五环之间的积水明显减弱，四环以内的大部分积水减弱。拓宽100%时，仅在二环到四环之间还有一些较浅的积水，拓宽140%时，六环内大部分积水消失。

需要指出的是城市内涝数值模式严重地依赖下垫面的地理特征和排水系统信息，因此必须不断对各种基础信息，包括泵站、蓄水池和调度策略等进行更新和调试才能持续获得良好的模拟效果。BUW模型模拟的下凹桥区的积水深度往往偏小，可能是未能精细描述桥区结构或对客水模拟过小引起的。考虑到下凹桥区在城市内涝中的重要性，针对每个桥区单独建模也是有必要的。根据不同降水和排水情景的模拟结果，发现现有的排水系统设计标准明显偏低，同时在拓建排水管网时需要统筹考虑，以免在别的地方引发地下涌出水或堵塞。此外，随着北京城区规模的不断扩大，西部和北部的山洪积水也必须纳入模拟范围，在边界衔接山洪模型或在BUW中发展山洪边界模块都是比较可行的做法。