真實感水力侵蝕地形建模與繪制

鄭顧平，張 萌，邵緒強，伍朝輝

真實感水力侵蝕地形建模與繪制

鄭顧平1，張 萌1，邵緒強1，伍朝輝2

(1. 華北電力大學控制與計算機工程學院，河北 保定 071003；2. 交通運輸部科學研究院基礎研究創(chuàng)新中心，北京 100019)

水力侵蝕地表是現(xiàn)實生活及虛擬場景中常見的一種自然現(xiàn)象。為研究水力對地表的侵蝕現(xiàn)象，提出一種采用光滑粒子流體動力學(SPH)模擬流體對地形表面的濕潤以及侵蝕過程的真實感仿真模型。通過對自然界中真實濺蝕現(xiàn)象進行合理簡化描述，建立一種包含3個階段的地形水力濺蝕模型：首先，利用彈丸沖擊平板靶材變形模型模擬水流沖擊地表變形效果；其次，構建基于達西定律的水流滲透量計算方法以確保地表吸水與水分擴散的真實性；最后，提出將地表濕潤度因素對地表侵蝕的影響考慮在內的新型計算模型，實現(xiàn)與沖蝕相融合模擬積水繼續(xù)運動對地表沖蝕作用。實驗結果表明，該模型能夠真實地模擬出不同方向的水流對地表的侵蝕效果以及侵蝕過程中地表吸水與濕潤擴散現(xiàn)象。

水力侵蝕；光滑粒子流體動力學；濺蝕；達西定律；濕潤擴散

水力侵蝕地形現(xiàn)象時常發(fā)生在大自然中，大到泥石流、山洪，小到滴水穿石，都是流體和固體發(fā)生交互產生的結果。水力侵蝕地形真實感繪制在防治類似泥石流等自然災害、游戲與影視特效等多個方面有著廣泛的應用。隨著計算機性能的飛速提升，基于物理的侵蝕現(xiàn)象動畫模擬技術取得很大進展。但是如何在保證高效計算速度的前提下更為真實的模擬出侵蝕現(xiàn)象一直是這些年來的研究熱點。

水力侵蝕在地質學與物理學也被廣泛研究，多數(shù)研究者將重點集中于水力對地形土壤層結構與土壤流失量的影響或者是對生態(tài)環(huán)境的影響的評估與預測[1-3]。水力侵蝕現(xiàn)象的研究對象可以分為流體和固體。在計算機圖形學領域，基于物理的流體模擬方法可以分為歐拉網格法[4-5]和拉格朗日粒子法[6-8]。由于歐拉網格法是在規(guī)則的網格上進行數(shù)值求解，所以遇到較大的形變時計算效率很低并且在模擬效果上容易失真。相比而言，拉格朗日粒子法從研究流體的各個微團物理屬性(密度、壓強、速度等)著手，更容易模擬出如水花四濺等細節(jié)，其中光滑粒子流體動力學(smoothed particle hydrodynamics，SPH)方法[9]在流體模擬中應用最廣泛。

水力侵蝕過程中固體和流體交互伴隨著碰撞、濕潤擴散以及侵蝕現(xiàn)象的產生。在固液交互的研究中，柳有權等[10]利用GPU加速模擬出氣流在迷宮和城市中擴散的場景，但由于內存有限并不能進行大規(guī)模問題求解；NARAIN等[11]實現(xiàn)了固液之間的雙向交互，其僅將固體作為流體的運動邊界，并未考慮固液物理性質的變化；LENAERTS和PHILIP[12]提出了一種多孔介質和流體的耦合方法，實現(xiàn)了干旱的土壤被雨水浸潤成泥漿或沙結構被海浪侵蝕的現(xiàn)象；常元章等[13]將Navier-Stokes (N-S)方程添加彈性應力項實現(xiàn)了牛頓流體與粘彈性流體的統(tǒng)一模擬；SHAO等[14]通過將PCISPH和LSM方法進行耦合逼真的模擬出速度差異較大的流體與可變形固體之間的交互；TAMPUBOLON等[15]擴展了基于SPH方法的多相流體框架，包括可變形固體、顆粒材料和多種流體。

對于侵蝕現(xiàn)象的研究，?T’AVA等[16]提出了以水力為主要影響地表變化因素的侵蝕交互模擬方法，且將地表劃分成多層土壤層，在重力作用下模擬了滑移侵蝕的效果；KRI?TOF等[17]通過定義SPH粒子攜帶固體粒子的體積分數(shù)表示固液兩相流的狀態(tài)，并且提出供體-受體模型以求解地表土壤的溶解量，同時模擬溶解的物質由于水的流動而在不同的地方形成沉積的現(xiàn)象；MEI等[18]采用GPU加速計算腐蝕的方法模擬了大規(guī)模的地形侵蝕現(xiàn)象，在保證模擬逼真度的前提下提升了求解的效率；CORDONNIER等[19]從生態(tài)環(huán)境因素考慮，引入一種支持侵蝕和植被雙向交互的框架，能夠模擬出植被、砂石等對侵蝕現(xiàn)象的影響。

濺蝕是水力侵蝕的一種，與其他侵蝕不同的是濺蝕現(xiàn)象由于具有較大動能的流體與地表發(fā)生碰撞造成地表變化。當水滴從一定的高度落地，動能作用于地表土壤，導致土壤結構破壞，表層土壤孔隙減少或者堵塞，土壤團粒被分散、剝離，并使部分土粒沿拋物線的運動軌跡向四周飛濺。為了方便在計算機動畫中模擬此類現(xiàn)象，本文忽略侵蝕對土壤層結構的影響，以SPH方法為基礎，提出了包含3個階段的近似簡化模型：首先，水滴濺落到未濕潤的地表上，由于地表還未來得及吸收液體，細小土顆粒被水滴的動能沖撞開[20]；其次，隨著降水量的增大，地表持續(xù)吸收水分，固液混成的泥漿受水滴沖擊飛濺[21]；最后，表面粒子吸水到達飽和，產生積水，積水的流動對地表造成沖蝕的效果。與此同時利用一種流體在多孔固體顆粒滲透力學的方法來模擬地表吸水濕潤以及濕潤的擴散過程，并且在侵蝕的過程中考慮地表濕潤程度對地表侵蝕的影響。通過分析水流沖撞地面各個階段的成因，將濺蝕與沖蝕相融合來模擬不同階段地形濕潤、變形以及濕潤的擴散現(xiàn)象。

1 基于物理的固液建模

1.1 SPH方法流體模型

SPH方法是拉格朗日數(shù)值求解方法的一種，其將模擬的流體離散成粒子，每個粒子攜帶質量、位置、速度、密度等屬性，通過研究各個粒子團的運動來研究流體整體的運動。流體的運動狀態(tài)主要由N-S方程求解得到，N-S方程包括質量守恒方程與動量守恒方程，即

由式(1)和式(2)簡化可得

其中，m，分別為粒子的質量與密度；(–,)為徑向對稱的核函數(shù)；為光滑核半徑。屬性()的梯度算子可由核函數(shù)(–,)的梯度算子求得

同理可求得屬性()的拉普拉斯算子

1.2 多孔固體顆粒模型

由于地面的土壤縫隙間存在空氣，水滴浸潤時先填滿土壤中空氣所占的空間，達到一定飽和度后產生積水。為了模擬這種吸水與積水的過程，本文利用簡化多孔固體顆粒模型代替土壤顆粒，每個粒子屬性有質量、飽和度等，即

對于地形的模擬采用了高度圖方法，即一張二維灰度圖，每個像素的灰度值表示該點的高度，然后將其映射到三維空間中，即可得到地形。為了保證模擬效果并提高計算速度，需對地形表面進行多孔固體顆粒采樣。如圖1所示，將地形投影到平面的一個二維網格，在每個小正方形網格內均勻采樣5個多孔固體顆粒。

圖1 地表多孔固體顆粒采樣圖

流體剛開始接觸多孔固體顆粒時產生浸潤，同時顆粒中的水分向周圍擴散，這種現(xiàn)象從微觀上看是流體小分子在多孔固體顆粒之間傳輸?shù)倪^程。DEAR[22]從大量的實驗統(tǒng)計與總結中得出描述這種滲流力學現(xiàn)象的達西定律

其中，為滲透系數(shù)；為距離對滲水量的影響權重。多孔固體顆粒的毛細勢P與飽和度呈正相關，即

其中，c為大于0的毛細系數(shù)。對于任一多孔固體顆粒而言，當其與周圍顆粒毛細勢不相同時，高毛細勢和低毛細勢產生壓力差，迫使毛細勢高的多孔固體顆粒中的液體向周圍毛細勢低的多孔固體顆粒中滲透。因此利用這種模型可以宏觀的模擬出整個濕潤的擴散現(xiàn)象，水分子在多孔固體顆粒間擴散如圖2(b)所示。

當降水量達到一定程度，從地表脫落的沙土與流體混合形成泥漿。泥漿是一種典型的固液兩相流[26-27]，對于泥漿的密度mix、粘性系數(shù)mix等物理性質由泥漿的濃度所決定

其中，為泥漿濃度(與上述飽和度不同的是，飽和度是多孔粒子吸水量，而濃度是液體粒子攜帶固體量)；，分別為流體粒子的密度和粘性系數(shù)；，分別為多孔固體顆粒的密度和粘性系數(shù)。

2 侵蝕模型

與已有對侵蝕的研究相比，本文擴展了侵蝕的種類。在流體對地形侵蝕階段，基于KRI?TOF等[17]的沖蝕模型，提出了一種包含濺蝕過程的侵蝕模型。同時引入達西定律保證水分在地表擴散物理過程的真實性，并且提出的侵蝕模型計算方法能夠描述濕潤因素對地表變形影響。

2.1 濺蝕

地表變形階段可采用彈丸沖擊平面靶材變形模型進行模擬。如圖3所示，當一定速度的彈丸撞擊平面靶材時，由于沖擊力的作用平面靶材發(fā)生形變，形變的寬度a和深度d與速度之間關系的計算方法為[28]

其中，為彈丸的半徑；為平面靶材的屈服強度。對于濺蝕的某個區(qū)域地表變形與其周圍多孔固體顆粒的濕潤程度有關，地表越濕潤越容易發(fā)生形變，所以地表的屈服強度會隨著濕潤值的增大而減小，基于此性質提出了描述屈服強度與地表多孔固體顆粒飽和度之間相互關系的計算模型

濺起或者卷入水中的沙土與液體形成泥漿，其密度與粘滯力由式(10)和式(11)求得，同時對于泥漿粒子，若其與周圍粒子有濃度差，會產生濃度擴散現(xiàn)象，該現(xiàn)象可通過以下模型進行求解

其中，為擴散系數(shù)。由于泥漿粒子中的固體微顆粒所受的力主要是重力，所以不能忽略重力對沉淀物沉降作用的影響，文中采用供體-受體體系[17]來描述重力影響下的沉降現(xiàn)象。如圖4所示，對于每個漿體粒子而言，其與周圍粒子之間的供受體關系由當前泥漿中固體微顆粒沉降速度與粒子間距離夾角決定。若·≥0，即粒子的濃度向粒子方向沉降，反之則粒子向粒子方向沉降。重力影響下SPH形式的濃度擴散式為

圖3 彈丸沖擊平面靶材變形圖

泥漿與融入積水中的固體顆粒會在沉降過程中附著在地表形成堆積。對于堆積量的求解可由供體-受體體系分析得到，在此僅僅將地表處多孔固體顆粒作為受體，泥漿作為供體。堆積是由供體濃度向受體傳遞產生的，即計算·<0的所有項。堆積量變化率與濃度變化率之間的關系式為

地表多孔固體顆粒吸水達到飽和時產生積水，積水沿著地表流動的動能相對于從空中降落時減少了很多，這一階段的侵蝕主要是沖蝕。

2.2 沖蝕

沖蝕效果是由流體對地表的剪切力造成的。當流體與地形表面發(fā)生接觸，地表處多孔固體顆粒受到半徑為范圍內的流體粒子剪切力的作用，該力將固體顆粒卷入流體中形成泥漿。剪切力的計算方法為

其中，為剪切力的大小；K為剪切力常數(shù)；為剪切率，即

其中，K為侵蝕常數(shù)；τ為標準剪切力大小。故式(22)可以求得邊界處的多孔固體顆粒質量的減少

其中，L為多孔固體粒子之間的距離。

其中，A為小網格的面積。

3 實驗結果及分析

基于沖蝕和濺蝕的融合模型，本文實現(xiàn)了一系列的侵蝕現(xiàn)象的真實感繪制。實驗硬件配置為：Inter(R) Core(TM) i5-7400 3.00 GHz CPU，8.00 GB內存，NVIDIA GeForce GT 720顯卡。編程環(huán)境為：Windows 10，Visual Studio 2013，Unity3D 2018。流體浸潤及侵蝕地形過程如圖5所示，繪制的過程中，先利用Marching Cube算法提取出流體表面的三角形網格，然后利用Unity3D對提取的表面進行渲染。實驗中重要參數(shù)設置見表1。

圖5 浸潤及侵蝕算法流程圖

表1 重要參數(shù)取值表

圖6描述了水分浸潤材料以及逐漸被吸收的過程，實驗中不考慮材料變形，從第100幀開始，每隔100幀提取一張效果圖。從序列圖中可以看出，濕潤從中間向四周擴散，靠近中心的位置顏色深，邊緣較淺，顏色的深淺反映不同區(qū)域水分的飽和度。當水分完全被材料吸收并擴散開形成第1 200幀的效果，效果圖也證實本文計算方法中引入達西定律能夠保證濕潤及擴散的真實物理效果。

不同屈服強度下地表的變形效果如圖7所示。當σ=0時，地表干燥，沖擊形成的水坑較淺 (圖7(a))；當σ=1時，地表濕潤達到飽和，沖擊形成的水坑深(圖7(c))。本文計算方法效果如圖7(b)所示，坑的深淺介于前兩者之間，證明該方法可描述出濕潤對地表變形的影響。

圖8(a)和圖9(a)展示了水流剛開始接觸地表且濕潤地表現(xiàn)象，從圖8(b)與圖9(b)可以看出水流對地表的侵蝕與地表吸水及濕潤擴散的效果，圖8(c)和圖9(c)是地表完全干涸之后形成的水洼。由此可以看出本文中濺蝕模型具有較好的模擬效果。圖8和圖9中水流垂直沖撞地表時，其侵蝕現(xiàn)象主要是濺蝕，表2前2行數(shù)據(jù)也印證了本文濺蝕模型的可靠性。

圖6 水分擴散序列圖

圖7 不同屈服強度形成水坑深度圖

圖8 單股水流垂直方向侵蝕地表模擬

圖9 雙股水流垂直侵蝕地表模擬

表2 實驗水流速度及濺蝕與沖蝕造成地表變形量

單股水流和雙股水流從水平方向流出與地表交互的效果如圖10和圖11所示。圖10(a)和圖11(a)為水流濕潤流經的地表圖；從圖10(b)和圖11(b)可以直觀看出濕潤擴散與侵蝕現(xiàn)象；當?shù)乇硗耆珊院笮纬傻乃疁先鐖D10(c)和 圖11(c)所示，對比圖中的水溝可以發(fā)現(xiàn)，水流剛開始接觸地表時，由于具有一定的垂直初速度會有濺蝕的現(xiàn)象，水流繼續(xù)往前運動時，沖蝕效果起主要作用。圖11(c)2股流體發(fā)生了沖撞，所以侵蝕形成的水溝略短于圖10(c)，但由于 圖11(c)中沖撞區(qū)域水流較多，水溝較圖10(c)中相比較深。同時表2中后2行數(shù)據(jù)也證實了本文中模型較好的將濺蝕與沖蝕現(xiàn)象應用到侵蝕現(xiàn)象的模擬中。

圖10 單股水流水平方向侵蝕地表模擬

圖11 雙股水流水平侵蝕地表模擬

圖12展示了字母“JOG”濕潤并侵蝕地表，圖12(b)和(c)均能看到較為明顯的字母侵蝕痕跡。圖13(a)水流與書本發(fā)生碰撞順流而下，流經區(qū)域呈扇形，在書本的邊緣形成3條分支，圖13(c)可以看出，水分干涸后呈現(xiàn)的3條淺溝壑。在圖14(a)中，水流沖下與水管發(fā)生碰撞，由于水流初始位置偏向水管左邊，可以看出水管左側流水量多于右側，圖14(c)能夠清晰看出水管左側侵蝕痕跡比水管右側明顯。

圖12 字母“JOG”濕潤及侵蝕現(xiàn)象

圖13 水流與書本碰撞侵蝕地表模擬

圖14 水流與水管碰撞侵蝕地表模擬

4 結束語

本文通過簡化自然界中濺蝕的過程提出了一種濺蝕模型來模擬流體侵蝕地表的現(xiàn)象，同時利用多孔固體顆粒來模擬地表濕潤的擴散，并考慮地表濕潤度對地表變化的影響，通過將提出的濺蝕模型與沖蝕模型相融合實現(xiàn)了水流對地表真實感的侵蝕現(xiàn)象模擬。

在未來的工作中希望能夠從以下幾點繼續(xù)改進：在濕潤擴散過程中，不應該僅僅是在地表擴散，水流滲透到地表以下對不同土壤層產生影響，對于濕潤的地表還應該考慮溫度對水分蒸發(fā)的影響；此外對于固液兩相流本文并未將固體粒子單獨模擬，而是用濃度來表達固液混合狀態(tài)；最后在碰撞檢測及響應方面精度與速度的提高也是未來的研究熱點之一。

[1] 唐棟, 李典慶, 周創(chuàng)兵, 等. 考慮前期降雨過程的邊坡穩(wěn)定性分析[J]. 巖土力學, 2013, 34(11): 3239-3248.

[2] 吳宏偉, 陳守義, 龐宇威. 雨水入滲對非飽和土坡穩(wěn)定性影響的參數(shù)研究[J]. 巖土力學, 1999, 20(1): 2-15.

[3] 尹占娥, 許世遠, 殷杰, 等. 基于小尺度的城市暴雨內澇災害情景模擬與風險評估[J]. 地理學報, 2010, 65(5): 553-562.

[4] FOSTER N, METAXAS D. Realistic animation of liquids[J]. Graphical Models and Image Processing, 1996, 58(5): 471-483.

[5] STAM J. Stable fluids[EB/OL]. [2019-07-21]. http://graphics.cs.cmu.edu/courses/15-869-F08/lec/10/ns.pdf.

[6] REEVES W T. Particle systems a technique for modeling a class of fuzzy objects[J]. Computer Graphics, 1983, 17(3): 359-376.

[7] TAKAHASHI T, FUJII H, KUNIMATSU A, et al. Realistic animation of fluid with splash and foam[J]. Computer Graphics Forum, 2003, 22(3): 391-400.

[8] STAM J, FIUME E. Depicting fire and other gaseous phenomena using diffusion processes[C]//Proceedings of the 22nd Annual Conference on Computer Graphics and Interactive Techniques-SIGGRAPH’95. New York: ACM Press, 1995: 129-136.

[9] MüLLER M, CHARYPAR D, GROSS M. Particle-based fluid simulation for interactive applications[C]//Proceedings of the 2003 ACM SIGGRAPH, Eurographics Symposium on Computer Animation. New York: ACM Press, 2003: 154-159.

[10] 柳有權, 劉學慧, 吳恩華. 基于GPU帶有復雜邊界的三維實時流體模擬[J]. 軟件學報, 2006, 17(3): 568-576.

[11] NARAIN R, SEWALL J, CARLSON M, et al. Fast animation of turbulence using energy transport and procedural synthesis[J]. ACM Transactions on Graphics, 2008, 27(5): 1-8.

[12] LENAERTS T, PHILIP D. Mixing fluids and granular materials[J]. Computer Graphics Forum: Journal of the European Association for Computer Graphics, 2009, 28(02): 213-218.

[13] 常元章, 柳有權, 鮑凱, 等. 一種基于粒子的牛頓流體與粘彈性流體統(tǒng)一模擬方法[J]. 計算機學報, 2010, 33(7): 1286-1294.

[14] SHAO X, ZHOU Z, MAGNENAT-THALMANN N, et al. Stable and fast fluid-solid coupling for incompressible SPH[J]. Computer Graphics Forum, 2015, 34(1): 191-204.

[15] TAMPUBOLON A P, GAST T, KLáR G, et al. Multi- species simulation of porous sand and water mixtures[J]. ACM Transactions on Graphics, 2017, 36(4): 1-11.

[16] ?T’AVA O, BENE? B, BRISBIN M, et al. Interactive terrain modeling using hydraulic erosion[C]// Proceeding of the 2008 ACM SIGGRAPH/Eurographics Symposium on Computer Animation. Aire-la-Ville: Eurographica Association, 2008: 201-210.

[17] KRI?TOF P, BENE? B, K?IVáNEK J, et al. Hydraulic erosion using smoothed particle hydrodynamics[J]. Computer Graphics Forum, 2009, 28(2): 219-228.

[18] MEI X, DECAUDIN P, HU B G. Fast hydraulic erosion simulation and visualization on GPU[C]//The 15th Pacific Conference on Computer Graphics and Applications (PG’07). New York: IEEE Press, 2007: 47-56.

[19] CORDONNIER G, GALIN E, GAIN J, et al. Authoring landscapes by combining ecosystem and terrain erosion simulation[J]. ACM Transactions on Graphics, 2017, 36(4): 1-12.

[20] 鄭顧平, 邢玥, 張榮華. 一種基于GPU的彈坑實時繪制方法[J]. 圖學學報, 2016, 37(4): 451-456.

[21] 劉嘉瑞, 楊猛, 吳佳澤, 等. 基于SPH的雨滴打擊不規(guī)則邊界的模擬方法[J]. 圖學學報, 2018, 39(3): 411-418.

[22] BEAR J. Dynamics of fluids in porous media[J]. Engineering Geology, 1973, 7(2): 174-175.

[23] 史鑫. 布料濕潤過程的模擬[D]. 上海: 上海交通大學, 2015.

[24] FEI Y, MAIA H T, BATTY C, et al. A multi-scale model for simulating liquid-hair interactions[J]. ACM Transactions on Graphics, 2017, 36(4): 1-17.

[25] WANG R, XIAO S J. Motion simulation of wet cloth on contact surface considering air effects and humidity[C]//Proceedings of the 15th ACM SIGGRAPH Conference on Virtual-Reality Continuum and Its Applications in Industry-VRCAI ’16. New York: ACM Press, 2016: 89-95.

[26] 楊志亮, 王章野, 柯曉棣, 等. 多相流災害場景的真實感建模與繪制[J]. 計算機輔助設計與圖形學學報, 2008, 20(8): 1023-1032.

[27] 楊志亮. 泥石流現(xiàn)象的真實感建模與繪制[D]. 杭州: 浙江大學, 2008.

[28] 張遠西. 混合水射流沖擊強化不同幾何特征表面的研究[D]. 鄭州: 鄭州大學, 2018.

Realistic modeling and rendering of hydraulic erosion terrain

ZHENG Gu-ping1, ZHANG Meng1, SHAO Xu-qiang1, WU Zhao-hui2

(1. School of Control and Computer Engineering, North China Electric Power University, Baoding Hebei 071003, China; 2. Basic Research Innovation Center, Academy of Transportation Sciences, Beijing 100019, China)

Hydraulic erosion topography is a common natural phenomenon in both real life and virtual scenes. In order to study the surface erosion by hydraulic power, a realistic simulation model is proposed by applying the smoothed particle hydrodynamics (SPH) to simulate the wetting and erosion process of the surface affected by fluid. Based on the reasonable simplified description of the real splash erosion in nature, a three stages topographic hydraulic splash erosion model is proposed. Firstly, the effect of water flow on the surface is simulated by the deformation model of projectile impacting flat target. Secondly, the calculation method of water permeability based on Darcy’s law is used to ensure the authenticity of surface water absorption and water diffusion. Finally, a new calculation model is proposed, which takes into account the influence of surface moisture on surface erosion. The experimental results show that the model can simulate the erosion effects of water flow in different directions on surface and the phenomena of surface water absorption and moisture diffusion during the erosion process.

hydraulic erosion; smoothed particle hydrodynamics; splash erosion; Darcy’s law; wet diffusion

TP 391

10.11996/JG.j.2095-302X.2020020169

A

2095-302X(2020)02-0169-09

2019-07-26；

2019-09-13

國家自然科學基金項目(61502168)；中央高?；究蒲袠I(yè)務費專項(2018MS068)；北京市自然科學基金項目(4182018)

鄭顧平(1960–)，男，河北石家莊人，教授，博士。主要研究方向為計算機仿真、分布式系統(tǒng)和人工智能。E-mail：zhengguping@126.com

邵緒強(1982–)，男，山東泰安人，副教授，博士，碩士生導師。主要研究方向為計算機圖形學、虛擬現(xiàn)實等。E-mail：shaoxuqiang@163.com