基于gpu并行计算的浅水波运动数值模拟-许栋.pdf

《基于gpu并行计算的浅水波运动数值模拟-许栋.pdf》由会员分享，可在线阅读，更多相关《基于gpu并行计算的浅水波运动数值模拟-许栋.pdf（8页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、基于GPU并行计算的浅水波运动数值模拟许 栋。， 徐 彬，David PAyet， 白玉川， 及春宁(天津大学水利工程仿真与安全国家重点实验室，天津300072)摘 要：利用有限体积法求解描述水流运动的二维浅水方程组，模拟洪水波运动传播过程，并通过GPU并行计算技术对程序进行加速，建立了浅水运动高效模拟方法。数值模拟结果表明，基于本文提出的GPU并行策略以及通用并行计算架构(CUDA)支持，能够实现相比CPU单核心最高112倍的加速比，为利用单机实现快速洪水预测以及防灾减灾决策提供有效支撑。此外，对基于GPU并行计算的浅水模拟计算精度进行了论证，并对并行性能优化进行了分析。利用所建模型模拟了溃

2、坝洪水在三维障碍物问的传播过程。关键词：洪水波；二维浅水方程组；GPU并行计算；CUDA中图分类号：0351；024221 文献标志码：A 文章编号：10074708(2016)010113081 引 言浅水波运动数值模拟具有广泛的工程应用价值及科学意义。如通过模拟溃坝水流运动、海湾涌潮过程以及海洋中海啸的生成、传播及海岸带入侵等13，预测洪水波的传播过程及淹没范围4，进而指导防灾减灾科学预案编制和应急决策。洪水运动和传播一般通过求解平面二维浅水方程(shallow water equations)来预测。考虑到洪水传播速度通常较快，为有效指导洪水发生后的快速响应决策，数值模拟本身的运行效率应

3、尽可能提高，尤其是对于大范围复杂区域的二维水流运动，如滞洪区或城市区域“1等。除在数学模型及算法方面改进外，并行计算也是提高模拟效率的有效手段。如江春波等口采用信息传递接口(MPI)通讯库，实现了二维浅水方程有限元并行求解，在8核计算机上获得了约67的加速比(相对于单核)。左一鸣等63在3台计算机组成的简单集群上对非结构网格二维水动力模型并行计算，获得了27倍加速比。杨明等73利用8核并行计算对黄河二维水沙数学模型的加速比约为64。除上述基于CPU的并行技术外，基于图形处理器GPU(GraphicalProcessing Unit)的并行计算近来发展迅速。NVI收稿日期：201410-10；修

4、改稿收到日期：20141227基金项目：国家自然科学基金(51579175，51379144)；国家自然科学基金创新研究群体(51321065)；天津市应用基础与前沿技术研究计划(12JCQNJCOS600，12JCQNJC02600)资助项目作者简介：许栋。(1980-)，男，副教授(Email：xudongtjueduon)DIA公司最新Kepler架构Tesla GPU已达到中型计算机集群的计算能力，如单个Tesla K40 GPU拥有2880个计算核心，单精度计算速度为429Tflops，双精度为143 Tflops，全局内存为12 GB，数据带宽可达288 GBs。与CPU并行相比，

5、GPU并行具有多方面优势。(1)GPU并行硬件价格低廉，且功耗远低于同等计算能力的CPU。(2)GPU可以实现计算的本地化，节省了远程服务器计算后网络传输数据所费时间，这对于洪水模拟非常重要。(3)GPU通过全局内存(global memory)实现数据共享，避免MPI带来的数据通讯交互性能瓶颈。(4)GPU计算结果可以直接调入显存，并利用显示器终端展示，为开发实时可视化水动力模拟系统提供了方便。GPU数值计算局限性主要在于并行计算模型的复杂性以及显存容量的限制。随着硬件不断发展，目前NVIDIA Tesla系列GPU显存容量一般在5 GB以上，已能够满足大部分二维水流数值模拟要求。另外NVI

6、DIA公司推出的CUDA(Compute Unified Device Architecture)，即通用并行计算架构，大大降低了GPU并行计算程序编写复杂性。它具体包括CUDA指令集架构(ISA)以及GPU内部的并行计算引擎。CUDA采用C语言作为编程语言，提供了大量高性能计算指令开发能力，使开发者能够在GPU的强大计算能力基础上建立一种效率更高的密集数据计算解决方案。近年来，基于CUDA架构的GPU并行计算在多个数值计算领域得到了广泛应用。GPU非常适合细粒度、大量线程并行执行的数值计算，对数据16=仉叮一蛆一娼竺一嚣一一一a=n=一M一报mn=学丝毒一力竺算型一t暑|一C=期月一钉嚣lf

7、；|；万方数据计 算 力 学 学 报 第33卷耦合性弱和逻辑简单的模拟加速效果显著，如分子动力学啪、基于LBM的流体动力学模拟9以及基于SPH方法的水动力模拟1 0|。针对二维浅水方程的GPU并行计算研究起步较晚，近年来相关尝试正陆续展开11-18。Asuncion等15对Roe格式求解浅水方程进行MPICUDA联合并行，在GTX480 GPU和4个CPU核心支持下获得了近50倍加速比。Smith等14提出了浅水模拟GPU并行的通用框架，并指出32位精度相比64位可能带来较大的局部误差。赵旭东等11在CUDA开发平台下建立了非结构化网格的二维水动力模型，利用GTX460显卡将计算效率对比集群机

8、提升了一个量级。本文利用GPU并行计算技术对结构网格下的二维浅水方程组数值求解进行加速，模拟洪水波运动传播过程。其中数值离散采用有限体积法，通量计算采用Local Lax Friedrich格式。基于Linux操作系统开发了CUDA架构支持的并行程序，对基于GPU并行计算的浅水模拟精度进行论证，并对并行计算性能优化展开研究。2控制方程及数值求解塑+型+掣一Is(，y (1)8t。8x av 。、1”。 埘驹+2 砌 h口2+g(h2一b2)zjsc口，一sa+s，一Lg。+o：；量：一兰袅别为沿z和Y方向的地形梯度，S丘和Sn分别为为方便并行程序设计，采用结构化网格。所求解变量口定义在控制体中

9、心，如图1所示，控制体n玎的状态由变量q玎描述，沿z方向通量E和沿Y方向通量G定义在网格边界上。i一1 jl2，汁I2 i+1图1有限体积法控制体示意图Fig1 Schematics of the control volume for finite volume method控制体内在时间t。t。+。内的通量增量可通过控制体积分计算：r+1 f F(口)dr+1 r船Eqndsdt Eq(x删删Y)卜l F(q) 一I I。 川2，)一J Jan 一一 J咋一12Eq(xfl2，Y，t)Jdydt+r”1 rHv2CEq(z嘞+1，2，)卜Jr J。pl2Gq(xi，Yi-12，t)dxdt

10、(2)通量构建是有限体积法求解的关键。洪水传播、尤其是溃坝和海啸等非恒定水流运动具有强间断特征，较常用的通量求解方法有基于Riemann求解器的Godunov格式1，3123和Roe格式m3等，这些格式计算流程较为复杂，不利于GPU并行计算。本文采用简单的Local Lax Friedrich格式143进行通量计算，定义控制体积分量Qij可打Jn。J：g d，以z方向为例，其计算格式为Fj暑。一F(Qi)“+F(Q；+。)“2一(口；+。：2)Q；譬一Q? (3)式中ai+i2为波速，利用ai+l2一“灭瓦耳瓦丽来计算。相比Roe等格式，Local Lax Friedrich格式在GPU并行时

11、具有以下优势。(1)计算流程简单，且不存分支(if语句)，便于线程束(wraps)中线程并发执行。(2)涉及中间变量少，而这些中间变量一般存储在与线程对应的寄存器(register)上，每个计算块(block)寄存器数量非常有限，该格式节省了寄存器资源，为大量线程并发执行提供前提条件。时间格式方面，考虑采用高阶格式(如Rungekutta等)需要存储大量历史数据，为了节省有限的GPU全局内存(global memory)，采用一阶显示欧拉时间格式：万方数据第1期 许 栋，等：基于GPU并行计算的浅水波运动数值模拟 115Q铲一Q等=(AtAy)Eq(x件加，y，)一EEq(xi-11z，Y，)

12、+(AtlAx)Gq(x，Y川2f)一G口(zi，为一l2，) (4)计算过程不考虑淹没过程中干湿判断，固壁边界初始时刻即全部设定。固壁边界采用反射边界条件，以z方向负方向固壁边界为例，h叫一Il，U叫一一Ui+l，j口i，j。Vi+I。3 浅水方程GPU并行策略及CUDA实现GPU支持大量线程并发执行，但与CPU相比，GPU没有足够大的缓存，并且存储体系构成与管理方式相对简单。为适应GPU的这种特点，在浅水方程数值计算中采用一维数据存储结构，具体策略：(1)将二维计算域数据初始化，设置初始速度(“，口)及水深(h)。(2)将二维计算域数据在host(主机)端按地址映射到一维数组。(3)利用C

13、UDA函数cudaMemcpy将一维数组复制到device(设备)端，存储在GPU的global memory(全局内存)中。(4)在GPU上调用kernel function(核函数)进行数值计算，并在时间域向前推进，直至全部时间步计算完毕。(5)将计算结果数据由device(设备)端复制到host(主机)端。(6)将一维形式存储的数据重新映射到二维计算域，并输出最终计算结果。整体计算流程如图3所示。在GPU端，一个thread(线程)对应it个单元(竹1)，以的取值根据单元总数和线程总数确定。GPU和CPU各自拥有独立的物理内存，必须调用cudaMemcpy函数才能实现数据交互，而数据交互

14、速度受带宽限制，往往成为GPU程序加速的瓶颈。在本文并行策略中，从图2的计算流程可以看出，所有数值计算均在GPU上完成，因此无需在计算过程中在GPUCPU之间传输数据，这为获得高效并行提供了前提条件。针对所采用的数值计算方法，在GPU上运行的kernel function(核函数)主要有两个，功能分别为数值通量计算(cornputeflux)和时间步进(tstepping)。代码采用标准C语言编写，并利用gcc以及nvcc联合编译，在Linux Fedora x64操作系统上运行。硬件系统中主机采用至强6核心CPU(Xeon E5-2600)，主频为20 GHz；GPU采用NVIDIA Tes

15、la C2075，共有448个计算核心，约6GB DDR5全局内存，05Tflops双精度计算效率，103 Tflops单精度计算效率。本文对GPU的并行计算主要采用单精度。4数学模型验证利用经典一维溃坝水流进行数学模型验证，以解析解作为验证依据。考虑一维瞬间溃坝问题，计算域长度为i000 m，溃口位于z=500 m处，初始时刻上游水深为10 m，下游水深为2 m。数值模拟采用网格如一1 m，溃坝发生30 S后的水面线如图3(a)所示。为比较CPU和GPU之间、单精度和双精度之间的计算差别，数值模拟采用CPU单精度、CPU双精度和GPU单精度3种方法。从图2浅水方程的GPU并行策略示意图Fig

16、2 Schematics of the GPU parallelization strategy for shallow water equations万方数据计 算 力 学 学 报 第33卷xm(a)计算值与解析解的比较(a)Comparison ofthe numerical results and the analytical solutionJm(b)计算精度比较(b)Comparison ofcomputational accuracy图3溃坝水流数学模型验证Fig3 Verification of the numerical model for dam-break flow图3(a

17、)可以看出，3种方法模拟结果基本相同(三线重合)，模型计算波前传播速度(波前位置)与解析解符合较好，也与Liang等17的结果大致相同，但在水深方面与解析解存在较大偏差，波前处水深偏小约30，与本文所采用的Local Lax Friedrich格式简单、数值粘性偏大16有关。尽管如此，作为对洪水传播的粗略预测，其计算结果仍是可以接受的。另外，本文重点验证GPU并行加速效果，具有更高精度的数值计算格式(如Roe)可在此框架下进一步扩展。为了精细比较CPU单精度、CPU双精度和GPU单精度3种方法对溃坝水流的计算精度差异，将不同方法计算得到的水深值相减获得差值，如图3(b)所示。可以看出，CPU单

18、精度与GPU单精度计算结果差异非常小，最大为8810一m，约为当地水深的0001；而双精度和单精度之间差异较大，最大约为0003 m，约为当地水深的003。然而相比因数值计算格式引起的偏差(30)，这部分误差仍是可以忽略的。因此，本文采用效率更高的单精度数据进行洪水运动模拟，相比双精度数据，可以节省50的存储空间以及50以上的计算量。5城市洪水演进模拟51计算工况与GPU实现为验证所建模型对建筑物密集的城市区域洪水演进模拟的效果，测试GPU并行计算加速的可行性，本文针对一个虚拟构建的城市区域进行水流数值模拟。为简化起见，将建筑概化为交错放置在平整地面上的大量长方体障碍物，其长度和宽度均1000

19、曼500图4 城市区域洪水演进模拟计算域Fig4 Computational domain for the simulation ofurban flood propagation为100 m。计算域为长8192 m，宽1024 m的矩形区域，共设置建筑物30座，如图4所示。在z和Y两个方向均采用间距1 m的均匀网格。网格总数约为839万。边界条件设置为，四周均为固壁边界条件，水库位于左侧z32768 m的区域，初始水位为20 m；为避免计算过程中的干湿判断，右侧区域初始水位设为01 m，另外障碍物设置为固定边壁，不受淹没。计算过程中忽略底阻力影响，主要计算参数列入表1。在利用CUDA进行GP

20、U编程时，一块GPU对应一个计算网格(grid)，该grid分为多个计算块(block)，而一个block分为多个线程(thread)。在计算能力(computability，NVIDIA对不同型号GPU定义的可计算性)20以上的GPU上，block和thread均可以采用一维、二维或三维的组织形式。以本文采用的Tesla C2075为例，单个block所包含thread数量最大为1024，而一个grid上的block数量最大可达655363。Grid以及block的大小划分直接影响GPU并行效率，划分依据要综表1 主要计算参数Tab1 Major computational paramete

21、rs万方数据第1期 许 栋，等：基于GPU并行计算的浅水波运动数值模拟 117150010005000一_1l r 1目盈：0 32 64 96 128160192 224 256CUDA block ID一鬃jij=j ，瀚j(a)Block ID分配(a)BlockID1500、10005000jI 3 579一11 13 15 一IJEtl二1-。lj l ，一一J量一L一i o， CUDAthread IDl瀵蕤雾粼冀漱雾0 2000 4000 6000 8000X【b)Thread ID分配(b)Thread ID图5 GPU端thread和block在整个计算域分配举例Fig5 E

22、xample of the distribution of GPUs threads and blocks in the entire computational domain合考虑任务特性和GPU本身的硬件参数。本文将block和thread均按一维划分，其好处是，(1)thread和一维数组的对应关系简单，易于编程实现。(2)每个block上thread数目可设置为32的整数倍，而GPU恰好以32个线程组成一个线程束(wrap)，在流处理器(SM)上并发执行，因此能够最大程度地利用计算资源。以上述81921024的计算网格为例，thread和block的分配如图5所示，可以看出，此算例共设

23、置256个block，每个block含16个线程，因此每个线程计算的单元数为81921024(25616)=2048个。在GPU程序运行过程中，利用CUDA提供的GPU性能监视工具获得GPU运行状态。TeslaC2075 GPU共有CUDA核心448个，全局设备内存5376 MB，内存使用率为29，GPU利用率为99，PCI数据带宽利用率为24。为了分析GPU计算时间占用细节，利用NVIDIA CUDAProfiler工具监视程序运行过程。此算例共设置256个block，每个block含16个线程，通量计算核函数(computeflux)和时间步进核函数(tstepping)交替占用GPU，占

24、用率分别为504和496。而通量计算核函数的运算指令数远大于时间步进核函数，造成两者时间占用非常相近的原因是GPU的效率瓶颈为全局内存的读写(在设置的block和thread配置下)，而非运算本身。两个函数所操作的全局内存容量相等，均为6个变量，因此两者整体运行时间接近。52并行效率分析及优化为提高效率，计算全部在GPU端运行，因此数值模拟过程中无需调用cudaMemcpy频繁地进行CPUGPU数据交互。为实现这一方法，需要保证GPU内存空间足够存储整个计算域中所有变量数据，目前的GPU硬件水平基本能够满足这一要求。以本文所采用NVIDIA Tesla C2075 GPU为例，最大可容纳约13

25、44亿计算网格，相当于能够以1 m1 m的网格覆盖1344 km2的区域，足以满足绝大部分计算任务，具体数据列入表2。非结构网格以及其他更高精度的时间格式(如Runge-kutta)将增加所需的内存空间。NVIDIA GPU的共享内存(shared memory)具有远大于global memory的带宽，以TeslaC2075为例，其shared memory带宽为144 GBs，而global memory仅为20 GBs。本文未采用shared memory进行程序优化，主要基于以下考虑。(1)二维浅水方程每个单元计算所涉及节点数量仅5个(J，_，)，(J-一1，J)，(J+1，J)，(

26、I，J一1)和(I，+1)，涉及全局内存读写操作次数较少。(2)直接操作全局内存方便实现内存对齐(coalesced)来加速读取，即其存储地址连续，且开始地址是每个thread所存取大小的32倍。(3)要充分利用shared memory加速需要避免内存冲突(bank conflict)，而这在多维数值计算中往往难以实现。(4)大量线程的并行执行能够在一定程度上掩盖全局内存操作的延迟。表2 GPU全局内存空间容纳网格数量估算Tab2 Estimation of grid amount allocable on the global memory of GPU注：h，“和t，分别为水深和两个流速

27、分量，nF“和Fhv分别为沿z方向的三个通量，Gh，Ghu和Ghv分别为沿，方向的三个数值通量，Z6为底高程。最大网格划分按正方形区域估算。万方数据计 算 力 学 学 报 第33卷GPU上Grid和Block大小划分对并行效率有很大影响，是GPU并行优化关键因素之一。一般来说，GPU上并行执行的线程数量越大，计算效率越高。但线程数量受硬件限制，且与寄存器使用有关。以Tesla C2075为例，每个block上寄存器总量为32768个，若核函数定义中间变量(存储在寄存器中)32个，则单block上线程数不能超过32768321024个。为测试Grid和Block划分效果，对GPU内置变量thre

28、adDim和gridDim进行变化组合，变化范围分别为1，10242和-1，16384，定义加速比AC为相对于CPU单核(20 GHz)的计算速度，测试不同Grid和Block组合的加速比变化，如图6所示。虽然GPU单核心时钟频率超过CPU的50，为1147 GHz，但单线程计算效率仅为CPU的28，此时全局内存操作产生严重计算延迟。而随线程数量增加，线程并发能够在一定程度上掩盖内存延迟。值得注意的是，加速比并非随着线程总量增加而单调增大。为了更清楚地观察加速比变化，提取threadDim分别为256，512和1024这三种常用设置下，加速比随blockDim的变化，如图6(b)所示。可以看出

29、，blockDim在8512之间变化时，整体加速比基本不变，而随Grid dimension(a)整体加速比分布(a)Overall distribution ofacceleration ratioBlock dimension(b)典型线程数条件下的加速比曲线(b)Acceleration ratio cuwe for typical thread sizes图6 GPU并行加速比随线程及计算块数量的变化Fig6 Variation of the acceleration ratio of GPU parallelcomputation with the size of threads a

30、nd blocks着blockDim的进一步增大，加速比迅速增大，最大至112，此时，每个线程仅对应一个计算单元。53模拟结果及讨论针对城市洪水传播进行GPU并行数值模拟，结果如图7所示。模拟结果整体符合真实的物理情景。溃坝水流遇障碍物后在上游发生反射，反射波不断向水库方向传播。障碍物前方形成水位壅高，最大壅高约为5 ITI，且随距溃口距离增大而减小。水流传播到末端共需时约280 s，最大传播速度约为15 ms，出现在计算域的两侧。在本算例中，计算网格数量约为839万，硬件采用TeslaC2075，单步计算花费GPU时间0024 S，全部计算用时仅144 S(不计CPU存储用时)。随着地理信息

31、系统不断发展，虚拟城市(virtual city)E183技术、三维城市建模(3DCM)19技术近年来发展迅速，能够对城市进行多分辨率、多尺度的三维描述，用于模拟和表达城市地形地貌、建筑物布局等。这些技术能够为水流数值计算提供必要的边界条件及基础数据支持4。实际上针对1 m分辨率网xm(a)侧视图t=O(a)Side view att=0翡=!=3：要二工：王潮Hm：l 2 3 4 5 6 7 8 9 10(b)t=80 S时水深分布(b)Waterdepthdistribution atf=80 s翳=1=j：潮Hm：l 2 3 4 5 6 7 8 9 10(c)，-180 S时水深分布(c

32、)Waterdepth distribution atf_180 S露=i：3=王：互警融H，m：1 2 3 4 5 6 7 8 9 108006007-400200(d)t=280 s时水深分布(d)Waterdepth distribution at，=280 S玎=丑囡一4500 5000 5500 6000 6500 7000x111(e)f=280 S时流速矢量(e)Velocity vectors at t=280 s图7城市区域洪水演进过程数值模拟结果Fig7 Simulation results on the urban flood propagation万方数据第1期 许 栋

33、，等：基于GPU并行计算的浅水波运动数值模拟 119格的城市数字高程(DEM)数据目前已能够利用遥感等技术获取。在这样的背景下，利用GPU并行技术建立简易、可靠的洪水仿真模拟系统，进行城市洪水的精细、实时和大范围模拟预测是可行的。6 结 论基于平面二维浅水方程及有限体积法建立了描述浅水波传播的数值模拟模型，并采用基于CUDA的GPU并行计算技术对程序进行加速，实现了浅水运动高效模拟方法。主要结论如下。(I)基于GPU全局内存进行计算，在不使用共享内存的条件下，针对二维浅水方程能够获得高达112倍的计算加速比，高效并行的原因为先进硬件设备支持、简化的计算流程和算法以及优化GPU并行策略。(2)线

34、程数量的增加有利于掩盖GPU全局内存操作延迟，为获得高的加速比，block数量应在512以上，单个block线程数量应在256以上。(3)在本文算例中，双精度和单精度对水深的计算偏差可达003，但单精度计算结果对特定的研究对象(如洪水传播速度)仍是可以接受的。(4)通量计算采用的Local Lax Friedrich格式简单，有利GPU并行，但数值粘性较大，计算精度较低，需要在进一步研究中向更高精度格式扩展。总体来讲，本文的数值测试证明，利用单GPU完成大规模洪水实时模拟是可行的，能够有效指导洪水发生后的应急疏散和快速决策。参考文献(References)：E132333白玉川，许栋溃坝波在弯

35、曲河道中传播的二雏数值模拟J天津大学学报，2007，40(7)：771778(BAI Yu-chuan，XU DongTwodimensional numerical simulation of dambreak wave propagation inchannel bendsJJournal of Tianjin University，2007，40(7)：771-778(in Chinese)张大伟，权锦，何晓燕，等基于地理信息系统的一维溃坝洪水分析系统研发及应用J水利学报，2013，44(12)：14751481，1487(ZHANG Dawei，QUAN Jin，HE XiaoyanS

36、tudy on uniform flow inopen channelJJournal of Hydraulic Engineering，2013，44(12)：14751481，1487(inChinese)白玉川，许 栋，王玉琦，等二维渍坝波遇障碍物的水流泥沙数值模拟j水利学报，2005，36(5)：538543(BAI Yuchuan，XU Dong，WANG Yuqi2D numerical simulation for water flow and sedimenttransport of dam-break wave EJJournal of Hydraulic Engineeri

37、ng，2005，36(5)：538543(in Chinese)43张大伟，程晓陶，黄金池建筑物密集城区溃堤水流二维数值模拟J水利学报，2010，41(3)：272276(ZHANG Dawei，CHENG Xiaotao，HUANG Jin-chi2D numerical modeling of levee-hreach floods indense urban areasJJournal of Hydraulic Engineering，2010，41(3)：272-276(in Chinese)E53江春波，陈立秋，张庆海。浅水流动的并行计算EJ清华大学学报(自然科学版)，2002，34

38、(11)：15481551(JIANG Chun-bo，CHEN Li-qiu，ZHANGQing-haiParallel algorithm tO analyze shallow waterflowJJournal of Tsinghua University(Sci&nch)，2002，34(11)：1548-1551(in Chinese)6左一鸣，崔广柏二维水动力模型的并行计算研究J水科学进展，2008，19(6)：846850(ZUOYiming，CUI Guang-baiParallel computation for 2Dflow modelJAdvances in Water

39、Scinece，2008，19(6)：846-850(in Chinese)73杨明，余欣，姜恺，等水动力学数学模型并行计算技术研究及实现J泥沙研究，2007(3)：1-3(YANG Ming，YU Xin，JIANG KaiResearch andapplication of parallel computing technology onhydrodynamic modelJJourmal of SedimentResearch，2007(3)：1-3(in Chinese)83 费辉，张云泉，王 可，等基于GPU的分子动力学模拟并行化及实现J计算机科学，2011，38(9)：275278

40、，287(FEI Hui，ZHANG Yun-quan，WANGKeParallel algorithm and implementation for molecular dynamics simulation based on GPUJComputer Science，2011，38(9)：275278，287(in Chi-nese)93张眷泽，程永光，李勇昌二雏浅水波方程格子Bohzmann算法的GPU并行实现J水动力学研究与进展(A辑)，20ll，26(2)：194-200(ZHANG Chunze，CHENG Yongguang。LI YongchangLattice bohzman

41、n model for tWOdimensional shallow water flowsrealizing parallel computing on GPUJChineseJournal of Hydrodynamics，201 1，26(2)：194200(in Chinese)10王珏，邱流潮应用基于GPU的SPH方法模拟二维楔形体入水砰击问题J计算力学学报，2013，30(s)：174177(WANG Jue，QIU Liu-chaoApplication of the GPUbased SPH method to the simulationof water entry of 2

42、D wedge bodiesJChinese Journal of Computational Mechanics，2013，30(S)：174-万方数据120 计 算 力 学 学 报 第33卷177(in Chinese)11赵旭东，梁书秀，孙昭晨，等基于GPU并行算法的水动力数学模型建立及其效率分析J-I大连理工大学学报，2014，54(2)：204206(ZHAO Xu-dong，LIANG Shu-xiu，SUN ZhaochenFoundation andanalysis of computational efficiency for hydrodynamicmodel based

43、on GPU parallel algorithmJ-IJournalof Dalian University of Technology，2014，54(2)：204-206(in Chinese)123张大伟，程晓陶，黄金池，等基于Godunov格式的溃坝水流数学模型J水科学进展，2010，21(2)：167172(ZHANG Dawei，CHENG Xiaotao，HUANG Jin-chiNumerical model for dam-break flow based onGodunovmethodJAdvances in Water Scinece，2010，21(2)：167172

44、(in Chinese)13Roe P LApproximate Riemann solvers parametervectors，and difference schemesJJournal of Computational Physics，1981，43(2)：35737214Smith L S，Liang QTowards a generalised GPUCPU shallowflow modelling toolJComputers&Fluids，2013，88：334-34315De La Asuncion M，Castro M J，FernandezNieto ED，et a1E

45、fficient GPU implementation of a two wavesTVD-WAF method for the twodimensional one layershallow water system on structured meshesJComputers&Fluids，2013，80(i0)：441-45216Press W HNumerical Recipes 3“Edition：the Artof Scientific ComputingMCambridge UniversityPress，200717Liang D，Falconer R A，Lin BCompa

46、rison betweenTVDMacCormack and ADI-type solvers of the shallow water equationsJAdvances in Water Re-sources，2006，29(12)：1833-1845(in Chinese)18王星捷，李春花基于Unity 3D平台的三维虚拟城市研究与应用J计算机技术与发展，2013，23(4)：241244(WANG Xing-jie，LI Chun-huaResearchand application of 3D virtual city based on unity 3DJComputer Tec

47、hnology and Development，2013，23(4)：241244(in Chinese)193戴洪宝，徐长海，张 龙基于Skyline的数字城市3维可视化系统的研究J测绘与空间地理信息，2013，36(10)：3234(DAI Hongbao，XU Chang-hai，ZHANG LongResearch on 3D visualization systemof digital city based on skylineJGeomatics&Spatial Information Technology，2013，36(10)：3234(in Chinese)Numerical

48、 simulation of shallow water motion basedon parallel computation using GPUxu Dong，XU Bin，David PAyet， BAI Yuchuan， 儿Chunning(State Key Laboratory of Hydraulic Engineering Simulation and Safety，Tianjin University，Tianjin 300072，China)Abstract：The propagation of flood waves was simulated by solving the twodimensional shallow waterequations fo