ECOMSED模型动边界处理技术研究

赵海亮，张艳飞

（1.河北省水运工程规划设计院，天津 300074；2.天津大学建筑工程学院，天津 300072）

ECOMSED模型动边界处理技术研究

赵海亮1，张艳飞2

（1.河北省水运工程规划设计院，天津 300074；2.天津大学建筑工程学院，天津 300072）

在ECOMSED数值模式的水动力模块基础上，引入了动边界处理技术，建立了数值试验，并采用单一斜坡地形和圆形浅滩斜坡地形对该技术进行了验证。数值试验结果显示，改进后的水动力计算模式可成功再现由涨落潮引起的潮间带漫滩现象，而且模拟得到的水面线较为光滑，且未产生零水深和负水深。

ECOMSED；动边界；数值模拟

Biography：ZHAO Hai-liang（1972-），male，engineer.

在天然海岸及河口地区，潮流涨落会形成广阔的潮间带，其中涨潮和落潮的水陆边界可相距数公里。在水动力数值模拟中，如果计算水陆边界不断移动，将会导致计算区域不断变化。在具有宽阔潮间带的近岸海域，动边界的合理处理不仅直接关系到流场计算的准确性，而且对准确模拟泥沙搬运和污染物输移扩散也至关重要，成为当前数值模拟中的热点问题。

在当前动边界技术（或称漫滩技术）的研究中，存在多种处理形式，其中有干湿法［1］、冻结法［2］、切削法［3］、线边界法［4］、开挖法［5］、窄缝法［6］、自适应网格法［7］等，各有优劣。干湿判断法是最为常见的漫滩处理方法，其理论简单，但较为粗糙；冻结法适用于宽浅、底坡平坦的露滩问题，对潮滩相间的海岸河口水域，则因水量的过分“冻结”而失真；线边界法适用于建筑物尺度远小于网格空间步长的条件；开挖法适用于滩地面积占整个计算域面积比例较小的情况，而对存在大范围浅滩的区域则易失真；自适应网格法是任意运动曲线坐标系下的可变动边界运动学边界条件，但该方法仅针对二维浅海潮波方程组，难以拓展至三维；窄缝法将动边界隐式处理，适合于岸边滩的处理，但对水域内的浅滩则不适用。实际上，干湿判断法虽为最常规的漫滩处理技术，但其容易实现，且稳定性良好，Oye［8］在对POM模式动边界的改进中，采用干湿判断法，并取得了良好的应用效果，本文继续采用干湿判断法的思路作为指导。

ECOMSED是一个模拟水动力、波浪和泥沙输运的三维数值模式［9］，其求解技术与POM类似，采用内外模式分裂法，网格划分采用平面正交曲线和垂向σ坐标系，可用以研究海岸河口区域的潮流、波浪、温盐、示踪物、粘性及非粘性泥沙的时空分布等现象，被广泛的应用于近岸海域的水流、泥沙数值模拟中。然而，由于ECOMSED原始代码中采用固定边界，无法模拟近海由于涨落潮而引起的漫滩流。综上所述，本文在ECOMSED水动力模块的基础上，进一步开发动边界处理技术，并利用数值试验、经典地形对所建模式进行测试，并论证其效果与可行性。

1 动边界技术的引入

ECOMSED模式中引入动边界时，首先必须区分陆地和滩地。为此引入绝对陆地边界的概念（absolute land boundary），在绝对陆地边界以外的单元定义为陆地，用FSM=0标识；在绝对陆地边界以内的单元则定义为或干或湿的单元，用FSM=1标识。为了便于判断，引入二维数组WETMASK，当FSM=0时，单元为陆地；当FSM=1且WETMASK=0时，单元为干点；当FSM=1且WETMASK=1时，单元为湿点。

其次判断单元的干湿状况。根据区域的水深变化范围和模拟精度要求，设定初始时刻滩地水深Hdry，近岸海域 Hdry通常取 10～20 cm。假设水位点总水深为 Di，j，水位为 ηi，j，Hi，j为平均海平面起算的水体深度，则 Di，j=Hi，j+ηi，j。具体判断方法如下：

（1）在二维外模式计算过程中，对每个时间步长计算域各单元的水位和断面平均流速进行扫描，并作以下判断。若单元水深 Di，j≤Hdry，则 WETMASKi，j=0；否则 WETMASKi，j=1；若湿单元周围（四边中点）都为干单元且Di，j≤1.01Hdry，则定义其为干单元，即WETMASKi，j=0；若（Di，j+Di-1，j）/2≤Hdry，则i，j=0；如果（Di，j+Di，j-1）/2≤Hdry，则i，j=0；若 WETMASKi-1，j=0 并且i，j＞0，或者 WETMASKi，j=0 且i，j＜0，则i，j=0；若 WETMASKi，j-1=0并且，j＞0，或者 WETMASKi，j=0 且i，j＜0，则i，j=0。

（2）三维内模式计算过程中，对三维流速进行判断：若 WETMASKi，j×WETMASKi-1，j=0，则 Ui，j，k=0；若WETMASKi，j×WETMASKi，j-1=0，则 Vi，j，k=0。

2 模型算例

2.1 单一斜坡地形

该算例地形为单一斜坡，沿岸长度为14 400 m，离岸长度为12 000 m。左端为水位开边界（纵向第2个网格中心），水深4.0 m，右端为闭边界，水深0.5 m。网格离散采用空间步长dx=dy=1 200 m，外模式时间步长为1.0 s，内模式为10.0 s，并持续模拟3.0 d。左端开边界处设置正弦潮波，振幅为2.0 m，周期为12.0 h。模拟中忽略柯氏力和表面风应力的影响，干湿网格的判断标准为0.05 m。

模拟结果（横向第6个网格）见图1。由图1可见，模拟所得的水面线光滑，未产生明显的数值震荡，并且涨落潮过程中干单元水深（最大值0.499 9 m）均小于干湿判断标准，且未出现零水深和负水深。

图1 单一斜坡地形下涨落潮过程干湿判断模拟Fig.1 Simulation of inundation for single-slope beach caused by flood-ebb process

2.2 圆形浅滩斜坡地形

本算例地形为圆形浅滩斜坡，横向12 000 m，纵向14 400 m，左端为水位开边界，水深zbl=6.0 m，右端为闭边界，水深zbr=0.5 m，其他参数设置与算例1相同。其详细地形如下

由图3可见，采用本文提出的干湿判断方法，对该地形下由涨落潮引起的潮间带进行模拟时，落潮过程水面线较光滑，涨潮过程在浅滩右侧产生了一定程度的振荡，但可以看出此振荡由浅滩地形所引起，并非干湿判断方法所致。由图4也可看出，本文的干湿判断方法能够较好地处理潮间带涨落潮过程。

图2 圆形浅滩斜坡地形示意Fig.2 Sketch of sloping circular shoal

图3 圆形浅滩斜坡地形下涨落潮过程干湿判断模拟Fig.3 Simulation of inundation for sloping circular shoal caused by flood-ebb process

图4 圆形浅滩斜坡地形下涨落潮流场Fig.4 Flow fields of sloping circular shoal caused by flood-ebb process

3 结论

（1）针对ECOMSED无法处理漫滩过程的缺点，提出了一种适合ECOMSED水动力模块的动边界处理技术。该方法物理概念明确、方法简单，能够保证计算域各单元的物质守恒，同时不会产生负水深，保证计算顺利进行。

（2）采用了单一斜坡地形、圆形浅滩斜坡地形，对所提出的动边界处理技术进行了验证，结果显示模拟得到的涨潮、落潮过程水面线均光滑，能够有效模拟潮间带的涨潮落潮过程。

［1］LEENDERTSE J J.A water quality simulation model for well-mixed estuaries and coastal seas［M］.Santa Monica：Rand，1970.

［2］程文辉，王船海.用正交曲线网格及“冻结”法计算河道流速场［J］.水利学报，1988（6）：18-25.

CHENG W H，WANG C H.The calculation of flow pattern in rivers by the orthogonal curvilinear coordinates and“condensation”technique［J］.Journal of Hydraulic Engineering，1988（6）：18-25.

［3］严以新，宋志尧.西江水道工程整治效果数值研究［C］//成功大学.第十七届海洋工程研讨会暨1995两岸港口及海岸开发讨论会论文集.台湾：成功大学，1995.

［4］孙英兰，陈时俊，王丽霞，等.胶州湾环流和污染扩散的数值模拟Ⅴ、胶州湾变空间步长的嵌套模型［J］.青岛海洋大学学报，1989，19（1）：1-18.

SUN Y L，CHEN S J，WANG L X，et al.Numerical modeling of the circulation and the pollution dispersion in Jiaozhou Bay，an inland model with grid-nets of different meshes［J］.Journal of Ocean University of Qingdao，1989，19（1）：1-18.

［5］周发毅，郭勇，虞和莹.大规模多岛屿海域潮流场的数值模拟［J］.海洋工程，1995，13（4）：61-69.

ZHOU F Y，GUO Y，YU H Y.The numerical simulation of the tidal flow in a large sea area with multi-islands［J］.The Ocean Engineering，1995，13（4）：61-69.

［6］陶建华.波浪在岸滩上的爬高和破碎的数学模型［J］.海洋学报，1984，6（5）：692-700.

TAO J H.Numerical model of the wave swash and breaking on the shore［J］.Acta Oceanologica Sinica，1984，6（5）：692-700.

［7］史峰岩，林文心，魏更生.运动侧边界海洋问题的自适应网格模拟方法［J］.海洋学报，1997，19（2）：1-9.

SHI F Y，LIN W X，WEI G S.The self-judgment method of the movable lateral boundary for ocean problems［J］.Acta Oceanologica Sinica，1997，19（2）：1-9.

［8］OYE L Y.A wetting and drying scheme for POM［J］.Ocean Modelling，2005，9（2）：133-150.

［9］HydroQual Inc.A premier for ECOMSED Ver.1.3 users manual［M］.New York：HydroQual Inc.，2002.

Development of inundation scheme for ECOMSED

ZHAO Hai-liang1，ZHANG Yan-fei2
（1.Hebei Port and Waterway Engineering Planning and Designing Institute，Tianjin300074，China；2.School of Civil Engineering，Tianjin University，Tianjin300072，China）

On the basis of the hydrodynamic module of ECOMSED program，an inundation scheme was developed.Two numerical experiments，including a single-slope beach and a sloping circular shoal，were applied to test the proposed inundation scheme.The results show that the inundation phenomenon caused by tidal flood and ebb process can be successfully reconstructed by the usage of this scheme，and the water surface profile obtained by simulation is smooth，which do not generate zero/negative water depth.

ECOMSED；moving boundary；numerical simulation

TV 148.2；O 242.1

A

1005-8443（2010）05-0395-04

赵海亮（1972-），男，河北省唐山人，工程师，主要从事港口、海岸及近海工程研究。