黃河三角洲孤東近岸海床沖刷平衡深度預(yù)判

凡姚申，王廣州，王開榮，竇身堂，陳沈良，徐叢亮，陳俊卿

（1.黃河水利科學研究院水利部黃河下游河道與河口治理重點實驗室，河南 鄭州 450003； 2.華東師范大學河口海岸學國家重點實驗室，上海 200241； 3.黃河河口海岸科學研究所，山東 東營 257000）

受海平面上升和人類活動雙重脅迫，海岸侵蝕已成為全球性的嚴重問題［1-2］。 大河三角洲自然資源豐富、發(fā)展條件優(yōu)越、區(qū)域位置重要，其海岸侵蝕一直是學術(shù)界和管理者尤為關(guān)注的焦點，如密西西比河三角洲［3-4］、湄公河三角洲［5］和尼羅河三角洲［6］。 現(xiàn)代黃河三角洲屬于弱潮河控型河口三角洲，其由高輸沙能力的河流注入弱潮動力環(huán)境海域所形成，曾是世界上堆積速度最快、演變最劇烈的河口三角洲［7］。 近年來，在流域來水來沙減少和河口流路不斷出汊改道的新情勢下，除現(xiàn)行河口口門附近外，黃河三角洲多數(shù)區(qū)域都表現(xiàn)出不同程度的侵蝕［8-10］，其中孤東海岸在該區(qū)域侵蝕最為強烈。 受孤東大堤保護，該區(qū)域的岸線得以穩(wěn)定，但堤外海床從1996年開始一直處于強沖刷態(tài)勢，沖刷深度普遍為6～8 m［11］。 學者們針對該區(qū)域的海床沖刷演變過程與影響因素［12-13］、沖刷動力機制與防護策略［14-15］做了大量研究工作，但對其沖刷平衡深度及沖刷極限深度卻鮮有關(guān)注，也一直沒有定論。

突變是指自然界中物體狀態(tài)的突然且徹底轉(zhuǎn)換。突變理論最早出現(xiàn)于《Structural Stability and Morphogenesis》一書，書中提出了突變理論的相關(guān)概念，并提出結(jié)構(gòu)穩(wěn)定性、動態(tài)穩(wěn)定性和臨界集是解析常微分方程解最重要的3 個要素［16］，為日后突變理論的應(yīng)用與拓展打下了基礎(chǔ)。 近年來，隨著非線性科學研究的不斷深入，突變理論也獲得新的突破，在多個學科中都有涉及［17-22］。 Cubitt 等［23］運用突變理論解釋了大陸坡泥沙輸移，最早把突變理論用于泥沙問題研究。 王協(xié)康等［24］用尖點突變理論對泥沙起動的規(guī)律與機理進行探討。 拾兵等［25］基于尖點突變理論和泥沙運動力學機制，建立了波浪作用下的泥沙起動模式，并證實了該模式探討泥沙運動規(guī)律的有效性。 楊少鵬等［26］利用由泥沙突變理論建立起來的泥沙運動模式，推導(dǎo)了海底管線沖刷坑深度的預(yù)測公式。

海床沖刷涉及水深、波流和泥沙的共同作用，其中水流強度、泥沙要素以及水流與泥沙顆粒間的作用關(guān)系是最為重要的作用因素［27］。 波流剪切應(yīng)力與泥沙臨界起動應(yīng)力的大小關(guān)系決定了海床泥沙的起動、輸移，也是海床沖刷發(fā)生與否的重要因素［28］。 當波流剪切應(yīng)力大于泥沙臨界起動應(yīng)力時，就會引起泥沙的起動。 臨界起動應(yīng)力很大程度上與底質(zhì)泥沙粒徑、沙粒密度、水深等有關(guān)，對于粗顆粒泥沙，泥沙粒徑越細、水越淺，需要的臨界起動應(yīng)力越小，但對于細顆粒泥沙，顆粒間的黏結(jié)力使泥沙發(fā)生絮凝，粒徑越小反而需要的臨界起動應(yīng)力越大［29］。 海床沖刷往往先沖走較易起動的泥沙，留在床面的泥沙越來越難起動，水深也越來越深，泥沙起動需要的剪切力就越來越大，直至剪切力小于臨界起動力時，所有泥沙無法起動，沖刷不再發(fā)生，海床沖刷達到平衡深度［30］。 因此，本研究選取海床沖刷達到穩(wěn)定狀態(tài)后的沖刷深度作為狀態(tài)變量，以泥沙突變理論為基礎(chǔ)，探討泥沙運動的突變關(guān)系，預(yù)判黃河三角洲孤東近岸海床沖刷平衡深度，以期為新情勢下黃河三角洲海岸侵蝕防護與綜合治理提供科學依據(jù)。

1 研究區(qū)域與資料來源

孤東海岸位于現(xiàn)行黃河口西北，是黃河口尾閭河道改道清水溝流路初期，即行河于北股河時期［31］，入海泥沙堆積形成的淤積質(zhì)海岸（見圖1）。 該地區(qū)石油資源豐富，是勝利油田石油開采海上戰(zhàn)略的主陣地。這里油井密布，防波堤道路縱橫，海灘油井圍堤連片，油井生產(chǎn)安全是迫在眉睫的大事。 然而，1980年黃河入海流路穩(wěn)定在北股河以南的區(qū)域后，孤東海岸失去黃河泥沙的直接補給，逐漸成為侵蝕型海岸。 受孤東大堤的防護，該區(qū)域的海岸線處于固定狀態(tài)，不再向陸蝕退，但堤外海床沖刷強烈。 1987年孤東大堤建成時，堤外潮灘灘面高程為0.5 m 左右，2003年測得孤東海堤迎海面海床平均高程約為-4.5 m，可見沖刷速度之快［32］。 合理的海堤結(jié)構(gòu)、設(shè)計標高等可以保證其穩(wěn)定性可靠，但海堤穩(wěn)定性除與自身條件有關(guān)外，也受海床沖淤等影響。 在強烈的海床沖刷下，海床等深線向岸移動、岸灘坡度變陡、高灘喪失［33］，不僅增加波浪、風暴潮對大堤的作用強度，還威脅海堤根基穩(wěn)固。 孤東海岸強烈的海床沖刷不斷影響到孤東大堤安全，部分海岸工程建筑物（如海堤、丁壩、灘海路堤等）受到?jīng)_刷毀壞，孤東油田和海港安全受到嚴重威脅。 從圖1 中的兩張現(xiàn)場照片可以看出孤東大堤建筑物工程受損的嚴重性。

海床沉積物資料來源于2013年7月在黃河三角洲臨近海域用帽式采泥器采取的表層沉積物，采樣深度為5 ～10 cm，采樣點位置見圖1。 收集的孤東近岸水深剖面資料共27 個剖面，平均間隔500 m，由水利部黃河水利委員會所屬水文水資源勘測隊測量得到。水深測量基面為黃?；鶞拭?，利用SDH-13D 型測深儀測水深，定位采用GPS 信標機，精度在3 m 以內(nèi)。

圖1 黃河三角洲孤東海岸位置及近岸海域測深剖面分布

以上資料用于預(yù)判和驗證孤東海床沖刷平衡深度，而平衡深度的預(yù)判方法是基于尖點突變模型建立起來的，因此下文著重介紹近岸海床泥沙尖點突變模型及平衡深度判別式的建立方法。

2 近岸海床泥沙尖點突變模型的建立

2.1 尖點突變模型的基本方程

在自治動力學系統(tǒng)的基礎(chǔ)上可以建立突變理論中的尖點突變模型，針對其系統(tǒng)狀態(tài)變化現(xiàn)象，可以通過分析其勢函數(shù)變化來解釋［20-22］。 尖點突變模型的標準形式如下：

式中：x為狀態(tài)變量；u、v為控制變量。

尖點突變模型中臨界點由一階導(dǎo)數(shù)來確定，一階導(dǎo)數(shù)為0 是我們關(guān)心的，即

式（2）是三次方程，其根的個數(shù)由以下判別式?jīng)Q定：

當Δ＞0 時，式（2）有一個實根；當Δ＜0 時，其有 3個互異的實根；當Δ＝0 時，其有一個二重根（u、v均不為 0）或一個三重根（u、v均為 0）。

式（1）的二階導(dǎo)數(shù)為

在尖點突變模型中，突變流形是式（2）確定的曲面M。 式（2）和式（4）聯(lián)立消去x得到方程：

他在國際上首先發(fā)現(xiàn)波里馬油菜細胞質(zhì)雄性不育，為國際雜交油菜實用化鋪平了道路；他在國內(nèi)首次育成甘藍型油菜自交不親和系及其雜種，育成優(yōu)質(zhì)高產(chǎn)雜交油菜，為中國油菜生產(chǎn)作出了重要貢獻。

分叉集是式（5）所確定的曲線B。 畫出M和B的圖形見圖2（a），該圖表明分叉集B突變流形的褶皺可以很好地投影在u—v平面上。 由路徑Ⅰ可以看出突跳和滯后，跳躍的那一段表明某些x值是不可達到的。發(fā)散性可由參數(shù)平面上兩條很接近的路徑Ⅱ和Ⅲ在突變流形上到達的截然不同的位置可以看出。

圖2 泥沙運動尖點突變理論模型示意（修改自文獻［25］）

2.2 沖刷平衡深度計算模型

構(gòu)建沖刷平衡深度計算模型，首先需確定穩(wěn)定變量和狀態(tài)變量。 本文根據(jù)海床泥沙起動規(guī)律，定義海床達到穩(wěn)定狀態(tài)的沖刷平衡深度h為狀態(tài)變量，泥沙粒徑的無量綱參數(shù)d*和希爾茲數(shù)φ作為控制變量，見圖2。 將水流作用下的希爾茲數(shù)φ作為反映水流影響的控制因素，將泥沙粒徑的無量綱參數(shù)d*作為反映泥沙及其與水流邊界層關(guān)系的控制變量。 尖點突變模型是三維突變模型，將該三維突變模型轉(zhuǎn)化成3 個變量，φ和d*作為參數(shù)，h作為狀態(tài)變量。d*可以用泥沙平均粒徑（dm）和泥沙起動粒徑（d0）的比值表示，即d*＝mdm／d0，其中m為修正系數(shù)，這里用泥沙密實系數(shù)表示，是一個動態(tài)變量，m＝0.75-0.65／（2＋d60／d10），d60／d10在這里表示泥沙不均勻程度。φ和d*構(gòu)成了海床沖刷的交叉集，h、φ和d*構(gòu)成了突變流形，則得到波流聯(lián)合作用下的泥沙起動尖點突變模型。

參考文獻［25］建立的波浪作用下的泥沙起動模式，將尖點突變模型的標準方程坐標系Q（x，u，v）轉(zhuǎn)換為泥沙起動方程坐標系O（h，d*，φ），轉(zhuǎn)換方程組為

式中：a、b、c為點Q在O（h，d*，φ）坐標系中的位置坐標值；li、mi、ni（i＝1，2，3）為原坐標對應(yīng)新坐標的方向余弦。

如圖2（b）所示，h坐標軸與x坐標軸的夾角為 0，其方向余弦為1，即l1＝1。 注意到，在尖點突變中勢函數(shù)為

在O（h，d*，φ）坐標系中有勢函數(shù)方程：

坐標變換附加條件如下：

同為右手坐標系時需滿足：

因曲面經(jīng)過原點，故勢函數(shù)值為0。 將上述條件代入式（6），可得

其中β＝arctana，聯(lián)立f（0，0，0）＝ 0 及式（9） ～式（13）可得

當泥沙靜止時，沖刷深度為 0，故a＝0，β′為泥沙mdm／d0的最小代表值，本文取 0.6。 將式（11）代入式（6），聯(lián)立β＝arctana和β′＝b／cosβ可得

該式即為海床沖刷平衡深度的預(yù)測公式。

3 孤東近岸海床沖刷平衡深度

3.1 沖刷平衡深度預(yù)判

研究區(qū)海床底質(zhì)沉積物樣品實測粒度數(shù)據(jù)顯示，該區(qū)域海床沉積物中黏土含量占比為19.8%～27.9%，粉砂占比為62.1～77.0%，砂占比為2.2%～15.7%，沉積物類型以黏土質(zhì)粉砂為主，有少部分的砂質(zhì)粉砂和粉砂（見圖 3），泥沙中值粒徑在 6.31 ～18.40 μm 范圍內(nèi)。 由此可見該區(qū)域泥沙粒徑分布較廣，d*和φ不能完全確定。 根據(jù)泥沙實測數(shù)據(jù)可以推測：該區(qū)域φ在0.1～16 范圍內(nèi)；泥沙平均粒徑是起動粒徑的1.8 ～3.3倍，即dm／d0取值在1.8 ～3.3 范圍內(nèi)；根據(jù)黃河口多年實測數(shù)據(jù)可以確定m取值在0.69 ～0.73范圍內(nèi)，變動范圍不大，本文取平均值0.71。 根據(jù)這些取值，推算出d*的取值范圍為1.3～2.3。

圖3 孤東海床底質(zhì)沉積物類型

將以上取值代入式（15），確定d*與φ的范圍后，該式可以看作一元三次方程。 關(guān)于其解可根據(jù)下式判斷。

將d*與φ的取值代入式（15），發(fā)現(xiàn)所有的取值組合中都滿足Δ＞0，即式（15）只有一個實根，且為正值，此實根即沖刷平衡深度。

根據(jù)以上取值，計算出該區(qū)域沖刷平衡深度的范圍，如圖4 所示。 圖4 中橫軸表示不同的希爾茲數(shù)φ，不同的線型表示不同泥沙粒徑的無量綱參數(shù)d*，縱軸表示兩者搭配取值時計算出的平衡深度。 由圖4 可知，當d*一定時，沖刷平衡深度h隨著φ的增大而增大；同樣，當φ一定時，沖刷平衡深度h隨著d*的增大而增大。 從圖4 還可以看出，d*較小時平衡深度h隨d*變化的幅度大于d*較大時的。 如當d*＝1.3 時，平衡深度h隨φ的變化從6.2 m 增大到7.3 m，變化幅度為1.1 m；當d*＝2.3 時，平衡深度h隨d*的變化從7.0 m增大到7.6 m，變化幅度為 0.6 m，是前者的54.5%。φ較小時平衡深度h隨d*變化的幅度大于φ較大時的。 如當φ＝0.2 時，平衡深度h隨d*的變化從6.2 m 增大到 7.0 m，變化幅度為 0.8 m；當φ＝16時，平衡深度h隨d*的變化從7.3 m 增大到7.6 m，變化幅度為0.3 m，僅是前者的37.5%。 由此可見，d*和φ越小，平衡深度h變化幅度越大。

圖4 沖刷平衡深度計算結(jié)果

綜合以上分析可以看出，雖然泥沙粒徑不唯一造成d*和φ在一定范圍內(nèi)分布，但計算的該區(qū)域沖刷平衡深度分布卻比較集中。 計算值分布統(tǒng)計顯示，90%的計算值（圖4 中灰色陰影覆蓋的區(qū)域）集中在6.5 ～7.5 m范圍內(nèi)。

3.2 沖刷平衡深度檢驗

海岸剖面是近岸泥沙與水動力相互作用在地形上的表現(xiàn)形式，從海岸剖面演變的視角，可以發(fā)現(xiàn)海床沖淤的未來趨勢［34］。 分析孤東近岸的海岸剖面演變過程，比較其未來發(fā)展趨勢與計算得到的平衡深度預(yù)測結(jié)果，進而可以達到檢測沖刷平衡深度預(yù)測結(jié)果的目的。 從前人研究分析來看，孤東海岸可以分為三部分，北大堤、南大堤和向海凸出處，北大堤包含剖面J1 ～J14，向海凸出處包含剖面J15 ～J19，南大堤包含剖面J20～J27，見圖 1（b）。

從剖面水深變化分析來看，孤東海岸剖面水深不斷增大，這也證實了該區(qū)域海床不斷發(fā)生沖刷。 值得注意的是，2012年以后孤東海岸各剖面水深增大趨勢有所減緩，特別是2015年和2017年距離起點稍遠的區(qū)域，剖面水深幾乎沒有變化，這一現(xiàn)象在孤東北大堤J5 剖面、向海凸出處J16 剖面和南大堤J23 剖面尤為典型（見圖5）。 J5 剖面距岸2.2 km 以內(nèi)水深下降明顯，而距岸2.2 km 以外2012年以后水深穩(wěn)定在7.0 m左右，基本處于動態(tài)平衡狀態(tài)。 以距岸2.2 km 劃分J5剖面附近的海床，以內(nèi)海床依然存在沖刷，以外海床2012年以后水深穩(wěn)定在7.0 m 左右。 J16 剖面距岸3.3 km以內(nèi)水深下降明顯，而距岸3.3 km 以外2012年以后水深穩(wěn)定在7.2 m 左右，基本處于動態(tài)平衡狀態(tài)。以距岸3.3 km 劃分J16 剖面附近的海床，以內(nèi)海床依然存在沖刷，以外海床2012年以后水深穩(wěn)定在7.2 m左右。 孤東南大堤J23 剖面也有同樣的演變趨勢。J23 剖面距岸 3.6 km 以內(nèi)水深下降明顯，而距岸3.6 km以外2012年以后水深穩(wěn)定在6.8 m 左右，基本處于動態(tài)平衡狀態(tài)。 以距岸3.6 km 劃分J23 剖面附近的海床，以內(nèi)海床依然存在沖刷，以外海床2012年以后水深穩(wěn)定在 6.8 m 左右。 由此可見，7.0、7.2、6.8 m水深分別是J5、J16、J23 剖面的沖刷平衡深度。這3 條剖面大致位于孤東北大堤、向海凸出處和孤東南大堤海岸中部，所在位置基本可以涵蓋孤東海岸各種動力環(huán)境，它們的沖刷極限深度在孤東海岸有一定的代表性。

圖5 孤東近岸J5、J16 和J23 剖面水深套匯

為了進一步驗證孤東海岸剖面的沖淤變化特征與沖刷平衡深度，本文計算了剖面沖淤變化速率Vd，其計算公式如下：

式中：DLm為前一個年份（用ym表示）的水深；DLn為同一條剖面后一個年份（用yn表示）同一起點距的水深。 這里計算了 2002—2005年、2005—2008年、2008—2012年、2012—2014年和 2014—2017年 5 個時間段的 J5、J16 和J23 剖面沖淤變化速率，計算結(jié)果（見圖6）表明，2012—2014年各斷面深水區(qū)剖面沖淤變化速率明顯減小，2014—2017年水深大于7.0 m的區(qū)域剖面沖淤變化速率接近于0。

圖6 孤東近岸J5、J16 和J23 剖面年均沖淤變化速率

海岸剖面演變過程綜合分析表明，6.8 ～7.2 m 水深是該區(qū)域沖刷的極限深度，在預(yù)判公式計算的6.5 ～7.5 m 的海床沖刷平衡水深范圍內(nèi)。 J5、J16 和 J23 剖面位于海床沉積物取樣點附近（見圖1（a）），預(yù)判公式計算出來的海床沖刷平衡深度最能代表它們的實際沖刷極限深度。 由此可見，基于泥沙尖點突變理論推導(dǎo)沖刷平衡深度計算公式，能夠滿足波流和泥沙共同作用下海床沖刷平衡深度的預(yù)測，將該理論應(yīng)用到海岸演變預(yù)測中是可行的。

4 結(jié) 語

孤東近岸海床是黃河三角洲海岸侵蝕最強烈的區(qū)域，為更好地評估該區(qū)域的海床沖刷發(fā)展趨勢，基于尖點突變理論推導(dǎo)出海床沖刷平衡深度預(yù)判公式，通過預(yù)判公式計算得到孤東近海海床沖刷平衡深度范圍為6.5～7.5 m，與分析海岸剖面演變過程得到的沖刷極限深度6.8～7.2 m 一致，說明該預(yù)測公式能滿足波流和泥沙共同作用下海床沖刷平衡深度的預(yù)測。 孤東海岸堤前海床目前水深為5 m 左右，距離該預(yù)判的平衡深度還有2 m，該區(qū)域還將持續(xù)沖刷，加強海堤防護工程日常維護十分必要。