長江口疏浚土利用的數(shù)理優(yōu)化

周 瑜，鄭偉安

（華東師范大學國際航運物流研究院，上海 200062）

1 前言

自19世紀40年代《南京條約》簽訂，以上海為代表的第一批通商口岸被迫開放以來，上海作為我國對外貿(mào)易往來的重要交通樞紐為眾人所熟知。特別是19世紀70年代，伴隨著黃浦江和蘇州河兩岸近代工業(yè)集聚區(qū)的興起，上海逐步奠定了其東方航運中心的地位。新中國成立后，借助于上海市及長江三角洲的緊密聯(lián)系，通過臨港工業(yè)的互動、國際航運中心建設的持續(xù)推動，上海港保持著上升的發(fā)展勢頭，維持著其強有力的競爭優(yōu)勢。

然而，在發(fā)展過程中，上海港及臨港產(chǎn)業(yè)正遭遇發(fā)展空間不足的瓶頸。

首先，在有限的地理維度中，上海和長江口深水岸線及土地已經(jīng)達到了近乎飽和的狀態(tài)。具體而言，目前上海港主力港區(qū)分布于外高橋、洋山、羅涇，而臨港產(chǎn)業(yè)主要集中在江口南岸、杭州灣金山、臨港新城及長江口長興島南岸。因而，重新尋找到完整的深水岸線及近岸土地已然成為制約發(fā)展的難題。

其次，上海港集裝箱碼頭和散雜貨碼頭的能力已趨超負荷。以外高橋港區(qū)為例，2011年外高橋港區(qū)集裝箱吞吐量已是設計能力的1.39倍，而洋山港區(qū)吞吐量也已達到了設計能力的1.41倍。有研究推算，上海港需要維持每年2×106～3×106TEU增量才能達到長江流域對外貿(mào)易發(fā)展的需求，而臨港臨海可用工業(yè)岸線土地資源已無法滿足。

再者，上海港作為建設中的國際航運中心，功能尚未完善，尤其水水中轉能力較弱。目前部分集裝箱的水水中轉功能主要由洋山港兼顧，大宗散貨的中轉均有羅涇港區(qū)和外高橋港區(qū)兼顧，但受各港區(qū)自身的條件限制，中轉能力較弱（特別是與長江流域各港口），無法使上海港成為完善的綜合運輸樞紐。

最后，上海周邊大型港口不斷崛起，上海國際航運中心的龍頭地位受到嚴峻挑戰(zhàn)。例如與鄰近的舟山港相比較，2011年，上海港包括內(nèi)河在內(nèi)的全港貨物吞吐量達到7.27×108t，而寧波—舟山港為6.94×108t，當年兩港增速分別為11.4%和12%，從增速看，勢均力敵。但2012年，寧波—舟山港增速快速領先上海港，達到7.2%，遠超上海港的1.1%，同時寧波—舟山港貨物吞吐量達到7.44×108t，首次突破7×108t，超過上海港800余萬噸，從而將上海港從雄踞多年的貨物吞吐量全球第一的位置上挑落下馬。從集裝箱增量看，2012年，上海港比2011年增加7.9×105TEU，而寧波—舟山增量為9.8×105TEU，此項上也是首次超過上海港。綜上所述，上海港面臨的問題集中在上海土地、岸線、航道等資源后備嚴重不足，制約了上海港及臨港產(chǎn)業(yè)發(fā)展空間及布局。在此港口資源發(fā)展壁壘的情況下，橫沙東灘的開發(fā)可謂是應運而生。

第一，橫沙東灘近期可提供的土地為1.7×105畝（1畝≈666.67m2大一區(qū)范圍），全部成陸后總共成陸的面積約5.8×105畝，可形成的岸線總長在90 km以上，并且可利用直接面向外海的有利條件，采用大型挖入式港池等手段，增加深水岸線，解決了上海港所無法滿足的土地和深水岸線資源的問題。

第二，橫沙東灘自身擁有非常好的優(yōu)勢，作為長江出海口的橋頭堡，可滿足江海聯(lián)運的緊迫要求。一方面，橫沙東灘通江達海，具有良好的地理位置。它位于長江出海口，扼守我國海岸線與長江黃金水道的T字形交點，若與長興島（海洋裝備島）用短距離隧道或橋梁連通后，即可經(jīng)滬崇蘇陸上通道直抵上海浦東和蘇北。此外，橫沙東灘與洋山深水港水域距離約100多千米，與外高橋港區(qū)（南港南岸）水域距離約30千米。借助其地理優(yōu)勢，可逐漸形成上海國際航運中心的港口群，確立“橫沙港—洋山深水港—南港”三足鼎立的格局，達到功能互補合理分配資源的目的。

另一方面，橫沙東灘可挖掘潛力大。第一，橫沙東灘已批準促淤圈圍的面積為112 km2（1.7×105畝，大一區(qū)）。在其東側還約有270 km2（4.1×105畝）的灘涂存在，可作為中遠期促淤圈圍規(guī)劃。第二，橫沙東灘北側有50 km以上岸線資源、緊貼北港航道，其中約14 km為目前10m以深的深水岸線，其余為目前7m左右深水岸線；南側約48 km岸線緊鄰長江口北槽12.5m航道。隨著北港航道整治規(guī)劃實施，岸線水深將不斷增加，存在巨大的開發(fā)潛力，東部及東側亦可形成連接深水15～20m的大型港池或人工島。第三，航道資源豐富。橫沙東灘地處橫沙島，南貼長江口北槽深水航道，北靠北港航道，西接長江黃金水道，東臨東海。長江口深水航道10m水深已向上貫通到南京長江大橋，2010年3月長江口的航道已達到12.5m水深，并將逐步向上延伸到南京長江大橋。北港航道規(guī)劃為10m航道。橫沙通道水深維持10m。第四，有效可利用泥沙資源充足。長江口豐水豐沙，大通站多年的年平均輸沙量為4.08×108t，近年來下泄泥量雖有所減少，但也有上億噸的泥沙下泄，長江口口門地區(qū)的灘涂淤漲能力依然較強。第五，成陸速度快。在橫沙東灘南側，長江口深水航道的年疏浚維護量可達7×107～8×107m3；在橫沙東灘北側，北港航道已在籌備建設，屆時也將產(chǎn)生大量的疏浚土。這些疏浚土就近吹填上灘，即可加快橫沙東灘的成陸過程，也可減少航道疏浚土的二次回淤現(xiàn)象，實現(xiàn)資源綜合利用。

因此，開發(fā)橫沙東灘，建設海港物流中心，完善國際航運中心功能成為迫在眉睫的發(fā)展需要。從近期看，開發(fā)橫沙東灘，形成新港區(qū)及臨港產(chǎn)業(yè)區(qū)（長江橋頭堡海港物流中心及臨海產(chǎn)業(yè)），與外高橋、洋山形成新的航運中心格局；從長遠看，橫沙東灘、九段沙、南匯邊灘這些長江口門區(qū)域的大片灘涂可綜合起來，形成陳吉余院士提出的長江口亞三角洲，建設上海真正的港口群及海洋經(jīng)濟帶。

綜上，為了配合上海國際航運中心的建設，利用長江流域的疏浚土在橫沙島吹填土地，建立新的港口或將成為一條非常行之有效的道路。那么將面臨的數(shù)學問題就是，如何將長江口岸的淤泥（圖1中兩條深色部分）吹填到它們所夾的粗線所圍區(qū)域內(nèi)，使得總工程耗費最少。

圖1 橫沙東灘及長江口淤泥示意圖Fig.1 East Hengsha shoaland the schematic diagram of siltat Yangtze River estuary

2 數(shù)學模型

將這個問題轉化為一個數(shù)學模型，就是困擾數(shù)學家232年的世界難題——“Monge-Kantorovich問題”。

1781年，Monge G提出了最優(yōu)輸運問題（即Monge問題）[1]：假如要挖一個給定形狀的土坑，搬到另一處堆出一個給定形狀的建筑物，怎樣運土可以使所作的功達到最小[2]，如圖2所示。

圖2 Monge問題背景Fig.2 M onge problem background

這個問題在一維的情況下早已經(jīng)解決。如果要把x處的土搬到y(tǒng)=f（x）處，記q（x）為x處的深度，p（y）為y處的設計高度，那么最優(yōu)的方案是f（x）是滿足p（y）dy=q（x）dx的單調(diào)遞增函數(shù)。但是在平面上的問題就要遠遠復雜了。主要原因是在平面上沒有像直線上的線性關系。

下面用數(shù)學語言把Monge問題重述如下：一個可分度量空間U，V上的測度μ、v和一個單位費用函數(shù)C（x，y），同時記

問題轉化為：尋找一個變換T0∈Γ1使得

其中，T0為最優(yōu)映射，Ⅰ（T0）為最優(yōu)運輸費用。

通常單位費用函數(shù)C（x，y）是非負的下半連續(xù)函數(shù)，常見的有：

線性費用：C(x,y)=|x-y|；

平方費用：C(x,y)=|x-y|2；

Lp費用：C(x,y)=|x-y|p,p∈(1,+∞)；

C(x,y)=?(|x-y|)，?為某一凸函數(shù)。

其中|x-y|表示x、y之間的直線距離，Monge最初研究的是線性費用情況。

自從Monge問題問世以來，它就成為了數(shù)學中幾個領域的數(shù)學家們所共同關心的經(jīng)典主題。微分幾何學家Appell P在1887年發(fā)表了一篇非常精彩的文章建立了一些最優(yōu)映射在平面和R3中的幾何性質[3]。但由于Monge問題的不適定性，即使在一般歐式空間中，其最優(yōu)輸運映射的存在性問題也不是那么顯而易見的。

為了解決這個困難，在1924年Kantorovich提出了輸運問題的弱解形式[4]，Kantorovich當時并不知道Monge的工作。Kantorovich問題簡述為：找到U×V上的概率測度γ0∈Γ={在U×V上邊際分布為μ。v的概率測度全體}，使得：

他的本質思想就是去尋找方案而不是輸運弱解，也就是說，去尋找空間U×V中的概率測定，而不是兩個空間之間的映射。在兩者之間的聯(lián)系中，可以注意到，任意的輸運映射T都可以推出一個方案γ，這個方案重點關注T在x×y中的圖譜，而且很容易證明反之亦成立。由于上述原因，任何一個輸運都可以推出一個具有相同成本的方案這就導致了

因此，所有的Γ內(nèi)的極值點都是由輸運推導而來的，那么可以直接從Kantorovich格式中得出輸運映射的存在性。Kantorovich格式可以看成是原來Monge問題的一個弱解形式。

至此，人們將上述最優(yōu)輸運問題統(tǒng)稱為“Monge-Kantorovich輸運問題”，該問題又引發(fā)了學術界更多的關注，Kantorovich與他的學生Rubinstein 得到了著名的對偶公式[5，6]，將“Monge-Kantorovich問題”很好地引入了經(jīng)濟學、自動化控制、運輸學、流體力學、幾何學、形狀優(yōu)化、氣象學和金融數(shù)學等眾多領域。1991年Brenier Y用凸函數(shù)的梯度刻畫了最優(yōu)映射[7]，他的這篇文章建立了最優(yōu)輸運問題與偏微分方程、概率論與泛函分析間的美妙聯(lián)系。從此最優(yōu)輸運問題變得極其流行，應用非常廣泛。

目前，這個問題的一維情況已經(jīng)被廣泛研究[2，8]，但是二維情況直到三年前才被鄭偉安和他的博士生沈銀芳所解決[9]，該文的想法是把原來的問題化為下列的二階擬線性橢圓偏微分方程的邊值問題。

這里未知函數(shù)u是一個二維的概率分布函數(shù)，ux與uy分別是u對x與對y的偏導數(shù)。藉此可以簡單地算出最佳的平方距離E|X-Y|2，它的平方根給出Monge原問題的一個均方意義下的解。用此可以估算出原問題中欲求的“功”。

對這個數(shù)學問題感興趣的讀者可以參考陳木法院士的文章[10].

3 數(shù)學證明的簡述

根據(jù)Kantorovich的表述，只需要找到兩個隨機的向量X=（X1，X2）與Y=（Y1，Y2），滿足下列條件：

1）X具有概率密度函數(shù)q（x1，x2），Y具有概率密度函數(shù)p（y1，y2）；

2）E|X-Y|2是所有滿足上述條件的隨機向量對中最小的。

文獻[9]中的想法是，假定（X，Y）就是欲找的最佳的一對。可以構造新的隨機變量Z=（X1，Y2）。于是

有趣的是，此時X與Z只有一個分量不同，Z與Y也只有一個分量不同。如果知道Z的密度函數(shù)，則根據(jù)前述的一維解法，X與Z的關系就固定了，同理Z與Y的關系也固定了。于是問題簡化為求Z的密度函數(shù)。更進一步，如果Z的概率分布函數(shù)知道了，那它對體積的導數(shù)就是Z的密度函數(shù)，所以問題又進一步簡化為求Z的概率分布函數(shù)u如下。

實際上，E|X-Z|2與E|Z-Y|2都可以用一個含有u的積分式表示出來。所以用數(shù)學中常用的變分法就可以從（5）得到（4）。

圖3 區(qū)域網(wǎng)格劃分Fig.3 M esh p lotsof the region

4 數(shù)值模擬算法及結果

本文用文獻[9]中的結論作為求解的依據(jù)再加上線性規(guī)劃的思想來模擬計算工程量。在設計算法時，首先將要求解的區(qū)域分成若干個小正方形單元，如圖3所示，圖中共劃分了20×30個單元，每個小正方形單元的邊長對應實際情況中的2 km。劃分好網(wǎng)格后，將該實際問題轉化為可求解的數(shù)學模型，如圖4所示。

圖4 數(shù)學模型示意圖Fig.4 Schematic diagram of themathematicalmodel

設計算法如下：

1）用常用的線性規(guī)劃尋找每一大塊（4 km2的方格）搬運的最優(yōu)方案；

2）在每個對應的小單元中，用前面所述的“Monge-Kantorovich輸運問題”解法分析每塊淤泥填到既定目標區(qū)的均方意義下的功。

線性規(guī)劃問題：初始給定兩個矩陣(aij)N×N,(bij)N×N的最優(yōu)化問題可以化為

式中，表示aij的部分移動到了bst,的花費為。

其拉格朗日函數(shù)為：

由于問題的數(shù)據(jù)量較大，通過KKT條件直接求解矩陣的逆，計算量會很大。因此使用交替方向法[11]和并行運算，求解該極小化問題：

初始化T0,λ0,μ0，選取步長h

取k=1,2,3,…

計算{

取j=1,2,3,…

分塊計算{

更新：

而在每個對應的小單元中，計算文獻[9]中的給出的變分問題：

其中

等價于求解極小化問題：

其中

r1，r2是常數(shù)。

使用交替下降法等計算方法，求解該變分問題的最優(yōu)值。算法如下：

初始化p0,q0,r0;

取k=1,2,3,…

計算{

更新：k=k+1．

在假設每平方千米單位深度的水里有a公斤泥沙，每頓泥沙搬運一千米需要b元的情況下，用上述算法初步估計出工程總耗為22 000ab元。當然這是在沒有任何水文資料與實際工程制約情況下的簡單估算，實際問題難得多。

5 模型的擴展

由于長江口的淤泥是隨著時間不斷增長的，而且還受到臺風的影響，所以在實際應用中會有更多的難度。因為這個問題還牽涉到了隨時間的變換。由于挖泥的形狀和次序影響以后新淤泥的堆積，目前還無法把它化為一個精確的數(shù)學題。但是可以從下列偏微分方程根據(jù)歷史數(shù)據(jù)確定懸浮泥沙的濃度。

式中，S為懸浮沙的濃度；u、v、w為x、y、z方向的流體速度；Kh、Kv為水平與豎直方向的擴散系數(shù)；Ws為泥沙沉降速度。

把這個方程與前面提到的邊界問題結合起來，就可以把問題分時間階段化為前面的“Monge-Kantorovich輸運問題”，例如把每三個月的淤泥形狀用前面的方法處理一下，再把它們累計起來，用計算機得到近似優(yōu)化值。

6 結語

隨著上海國際航運中心的龍頭地位日益面臨嚴峻挑戰(zhàn)，為了配合上海國際航運中心的建設，利用長江口流域的疏浚土在橫沙島吹填土地，建立新的港區(qū)或將成為一條非常行之有效的道路。本文以長江口疏浚土吹填橫沙東灘為背景，將其轉化為數(shù)學模型：“Monge-Kantorovich輸運問題”，以此設計算法，進而得到最有運輸條件下的總工程耗費。由于工程本身的規(guī)模巨大，通過尋找最優(yōu)的輸運方案，能夠最大限度地減少施工費用支出，降低工程造價成本。為了更好地模擬實際操作過程中所可能遇到的問題，全文最后還討論了模型進一步的擴展提升，提出了可將其轉化為一個與時間變化有關的問題，以加強模型的適用性。

[1] Monge G.Memoire sur la theorie des déblais et des remblais[J].Pairs ：Histoire de I′Acad.des Sciencesde Pairs，1781.666-704．

[2] Villani C.Optimal Transport[M].New York：Springer-Verlag，2009.

[3] Appell P.Mémoire sur les déblais et les remblais des systémes continus ou discontinues[J].Mémoires présentés par divers Savantsà l′Académie des Sciences de l′Institut de France，1887，29：1-208．

[4] Kantorovich L V.On a problem of Monge[J].UspehiMat Nauk，1948，3：225-226．

[5] Kantorovich L V.On the transfer of mass[J].Dokl Akad Nauk，SSSR.1942，37：227-229．

[6] Rubinstein G S，Kantorovich LV.On the space of completely additive functions[J].Vestnic Leningrad Univ，Ser MatMekh i Astron.1958，13（7）：52-59．

[7] Brenier Y.Polar factorization and monotone rearrangement of vector-valued functions[J].Comm Pure Appl.Math，1991，44：375-417．

[8] Rachev ST，Ruschendorf L.Mass Transportation Problems Volume I：Theory ，Volume II：Applications[M].New York ：Springer-Verlag，1998.

[9]Shen Y F，ZhengW A.On Momge-Kantorovich problem in the plane[J].CR Acad SciPairs，2010，1：267-271．

[10] 陳木法.談談概率論與其他學科的若干交叉[J].數(shù)學進展，2005，34（6）：661-672．

[11] Boyd S，Parikh N，Chu E，et al.Distributed optim ization and statistical learning via the alternating directionmethod ofmultipliers[J].Machine Learning，2010，3：1-122．