基于Kriging模型的工程土方量計(jì)算方法分析

洪成晶

(中國(guó)路橋工程有限責(zé)任公司，北京100010)

土方量的計(jì)算是工程建設(shè)施工中重要的一環(huán)，直接關(guān)系到工程造價(jià)，是有效控制工程成本的重要依據(jù)[1]。目前土方量的計(jì)算方法主要由方格網(wǎng)法、等高線法、斷面法以及DTM(數(shù)字地面模型)等，眾多學(xué)者對(duì)這些方法進(jìn)行了詳細(xì)地研究，如常青[2]對(duì)上述四種方法的精度和適用性進(jìn)行了研究分析，并用實(shí)例進(jìn)行了證明；陳愛梅[3]等基于南方CASS軟件平臺(tái)對(duì)上述方法進(jìn)行了對(duì)比分析，并提出了在使用過程中應(yīng)注意的事項(xiàng)。

本文基于方格網(wǎng)法，針對(duì)方格網(wǎng)法的不足，減少外業(yè)的工作量，提出一種計(jì)算土方量的新方法，并通過具體案例進(jìn)行分析驗(yàn)證，證明了本文方法的正確性和精確性。

1 方格網(wǎng)法概述

方格網(wǎng)法[4]是工程土方量計(jì)算的一種常用方法，其基本原理就是將待測(cè)區(qū)域劃分為多個(gè)正方形網(wǎng)格，根據(jù)每個(gè)網(wǎng)格頂點(diǎn)的原地面高程和設(shè)計(jì)高程采用平均值法計(jì)算每個(gè)方形區(qū)域的填、挖方量，最后進(jìn)行匯總，從而得到整個(gè)區(qū)域的填方和挖方的總量。基本的計(jì)算式如下：

(1)

式中:Hij為第i行第j列的方格網(wǎng)的高差；a、b為方格網(wǎng)的邊長(zhǎng)；n式待測(cè)區(qū)域的方格數(shù)量。

方格網(wǎng)法適用于大面積的土方量計(jì)算，尤其是地形變化較小，地勢(shì)較為平坦的區(qū)域。對(duì)于地形起伏較大的區(qū)域，由于方格網(wǎng)法中假定同一網(wǎng)格線上高程同，最終土方量的計(jì)算精度會(huì)有所降低。因此在計(jì)算地形起伏較大區(qū)域的土方量時(shí)，往往需要對(duì)網(wǎng)格點(diǎn)進(jìn)行加密，外業(yè)的工作量會(huì)隨之大大增加，不利于工作效率的提高。

為了彌補(bǔ)方格網(wǎng)法的不足，將Kriging模型引入其中。其基本思路為：

(1)在待測(cè)區(qū)域均勻采集原地面的高程，獲取一系列樣本點(diǎn)。

(2)利用外業(yè)采集的樣本點(diǎn)建立Kriging模型，實(shí)現(xiàn)對(duì)原地面地形的模擬。

(3)根據(jù)模擬的地形條件，采用方格網(wǎng)法劃分網(wǎng)格并進(jìn)行土方量的計(jì)算。

上述方法的優(yōu)勢(shì)在于在采用方格網(wǎng)法計(jì)算土方量時(shí)，網(wǎng)格頂點(diǎn)處的原地面高程可以由模擬出的地形直接給出，當(dāng)網(wǎng)格需要加密時(shí)，不再需要進(jìn)行外業(yè)測(cè)量，大大減少了工作量。

2 Kriging模型

2.1 Kriging模型基本原理

Kriging模型是1951年由南非工程師Krige提出的一種給予隨機(jī)過程的無偏估計(jì)模型[5]，最早用于礦產(chǎn)儲(chǔ)量分布的估計(jì)。在Giunta[6]、Sacks[7]等學(xué)者的推動(dòng)下，Kriging模型目前已在氣象、水文地質(zhì)、地理信息系統(tǒng)、航空航天的多個(gè)領(lǐng)域獲得廣泛的應(yīng)用。

Kriging模型是一種插值模型，在已知樣本點(diǎn)x=[x1x2…xn]T和對(duì)應(yīng)的n個(gè)函數(shù)響應(yīng)值ys=[y1y2…yn]T的情況下，其插值結(jié)果由已知樣本點(diǎn)的函數(shù)值線性加權(quán)得到，即：

(2)

因此，只要能給出加權(quán)系數(shù)ω=[ω1ω2…ωn]T，便可以估計(jì)出任一點(diǎn)處的函數(shù)響應(yīng)。為此，Kriging模型引入統(tǒng)計(jì)學(xué)的假設(shè)，即將未知函數(shù)看作是某個(gè)高斯靜態(tài)隨機(jī)過程的具體體現(xiàn)。換言之，對(duì)于任意位置x，對(duì)應(yīng)的函數(shù)響應(yīng)y(x)被一個(gè)隨機(jī)函數(shù)Y(x)代替，y(x)只是Y(x)可能的結(jié)果之一，即：

Y(x)=β0+Z(x)

(3)

式中：β0是未知常數(shù)，代表Y(x)的數(shù)學(xué)期望值；Z(·)是均值為0，方差為σ2的靜態(tài)隨機(jī)過程。

在設(shè)計(jì)空間中，不同未知處的隨機(jī)變量之間的相關(guān)性可以用下式進(jìn)行描述。

Cov[Z(x),Z(x′)]=σ2R(x,x′)

(4)

其中：R(x,x′)為相關(guān)函數(shù)，表示不同位置處隨機(jī)變量之間的相關(guān)性。

為使Kriging模型預(yù)估值準(zhǔn)確，要求下式所示的均方根誤差：

(5)

最小，且滿足無偏估計(jì)的條件。

引入拉格朗日方法，將上述問題轉(zhuǎn)化為下式所示的求最小值問題：

minH(ω,λ)=σ2(1+ωTRω-2ωTrx)-λ(FTω-1)

s.t.FTω-1=0

(6)

式中：R為相關(guān)矩陣，由樣本點(diǎn)之間的相關(guān)函數(shù)計(jì)算得到，R=(R(xi,xj))i,j∈n×n；rx為相關(guān)向量，計(jì)算方法為rx=[R(x1,x),R(x2,x),…,R(xn,x)]；F為n維單位列向量，F(xiàn)=[1 1 … 1]∈n。

對(duì)式(6)進(jìn)行求解，可得：

λ=-2σ2(FTR-1F)-1(FTR-1rx-1)

(7)

ω=R-1(rx-F(FTR-1F)-1(FTR-1rx-1))

(8)

并將其帶入到式(2)中，可以得到：

(9)

其中:β0可根據(jù)最大似然估計(jì)得到，其值為：

β0=(FTR-1F)-1FTR-1ys

(10)

2.2 相關(guān)函數(shù)

相關(guān)矩陣R和相關(guān)向量rx的建立都與相關(guān)函數(shù)的選取有關(guān)。目前許多學(xué)者采用高斯相關(guān)函數(shù)，主要由以下幾類[8]：

(1)高斯函數(shù)：

Rk(xi,xj)=exp[-θk|xi-xj|2]

(11)

(2)各項(xiàng)同性高斯指數(shù)函數(shù)：

Rk(xi,xj)=exp[-θk|xi-xj|p] 1≤p≤2

(12)

(3)各向異性高斯指數(shù)函數(shù)：

Rk(xi,xj)=exp[-θk|xi-xj|pk] 1≤p≤2,k=1,2,3,…

(13)

式中:θ，p為相關(guān)函數(shù)參數(shù)。

3 土方量計(jì)算案例分析

3.1 案例選擇

為對(duì)本文中提出的方法進(jìn)行驗(yàn)證，假設(shè)有一方形場(chǎng)地，其尺寸為200m×200m，并且該場(chǎng)地內(nèi)的地形滿足式(14)所示的函數(shù)關(guān)系式。該方形場(chǎng)地的地形如圖1所示。

(14)

圖1 原地面地形

從圖1中可以看出，該場(chǎng)地內(nèi)地形起伏較大，在四個(gè)波峰處(坐標(biāo)分別為(-50,-50)、(-50,50)、(50,-50)、(50,50))原地面標(biāo)高最大，為32.4m，在(0,0)處原地面標(biāo)高最小，為0m。據(jù)此可以計(jì)算出該場(chǎng)地內(nèi)最大平均坡度約為15°，若采用傳統(tǒng)的方格網(wǎng)法進(jìn)行土方量的計(jì)算，則需要對(duì)網(wǎng)格進(jìn)行加密，工作量十分巨大。

3.2 均勻采樣與地形模擬

為保證采樣點(diǎn)能夠均勻分布于場(chǎng)地之中，本文中采用均勻設(shè)計(jì)法來進(jìn)行采樣點(diǎn)的選取。均勻設(shè)計(jì)法是方開泰于1978年提出的[8]。首先產(chǎn)生均勻設(shè)計(jì)表，然后根據(jù)均勻性的度量標(biāo)準(zhǔn)如中心偏差、卷積偏差等來產(chǎn)生使用表，即最終的樣本點(diǎn)。通過這一系列的過程，使樣本點(diǎn)在設(shè)計(jì)空間內(nèi)均勻分散。

圖2為100個(gè)樣本點(diǎn)的分布情況，從圖中可以看出，樣本點(diǎn)在場(chǎng)地內(nèi)的分布十分均勻，滿足均勻性的要求。

根據(jù)均勻設(shè)計(jì)獲得的樣本點(diǎn)之后，通過式(14)獲取樣本點(diǎn)處在該方形場(chǎng)地內(nèi)的標(biāo)高，然后采用Kriging模型對(duì)該場(chǎng)地的地形進(jìn)行模擬，模擬出的地形如圖3所示。在模擬過程中，各采樣點(diǎn)之間的相關(guān)函數(shù)采用各向同性高斯指數(shù)函數(shù)，取θ=0.1，p=1。

圖2 樣本點(diǎn)均勻分布情況

圖3 Kriging模型擬合出的地形

對(duì)比圖1和圖3可以發(fā)現(xiàn)，Kriging模型能夠較好地模擬出地形的起伏特點(diǎn)，模擬效果較好。

3.3 土方量計(jì)算結(jié)果分析

假設(shè)需要對(duì)該場(chǎng)地進(jìn)行整平，整平后標(biāo)高為15m。由于描述地形的函數(shù)已知，可以用積分方法算出填挖方的精確值。表1列出了本文中方法在不同采樣點(diǎn)數(shù)量下的填挖方計(jì)算值。

從表中計(jì)算結(jié)果可以看出，隨著樣本點(diǎn)數(shù)量的增加，Kriging模型對(duì)地形的模擬越精確，最終填方量和挖方量的計(jì)算越準(zhǔn)確，同時(shí)也標(biāo)明本文中提出的方法的可行性與精確性。當(dāng)樣本點(diǎn)數(shù)量為120時(shí)，即樣本點(diǎn)的平均間距為18.3m時(shí)，填方量和挖方量的計(jì)算誤差均在5 %以內(nèi)。

比較樣本點(diǎn)的平均間距為20m時(shí)的計(jì)算結(jié)果，可以發(fā)現(xiàn)根據(jù)本文方法計(jì)算結(jié)果比傳統(tǒng)方格網(wǎng)法的計(jì)算結(jié)果更加精確，與精確解相比，本文方法計(jì)算的填方量偏差為5.2 %，挖方量的計(jì)算偏差為2.7 %，而傳統(tǒng)方格網(wǎng)法計(jì)算的填方量偏差和挖方量偏差分別為16.8 %和14.8 %，本文中方法的計(jì)算結(jié)果遠(yuǎn)遠(yuǎn)優(yōu)于傳統(tǒng)方格網(wǎng)法。從表中的數(shù)據(jù)可以發(fā)現(xiàn)，在相同的計(jì)算精度下，本文中方法所需的外業(yè)工作量遠(yuǎn)遠(yuǎn)小于傳統(tǒng)方格網(wǎng)法。

表1 不同樣本數(shù)量下的計(jì)算結(jié)果對(duì)比

4 結(jié)論

本文基于方格網(wǎng)法，引入Kriging模型，提出了一種計(jì)算工程土方量的新方法，并用該方法對(duì)一起伏較大的地形進(jìn)行了土方量的計(jì)算，主要結(jié)論如下：

(1)隨著樣本數(shù)量的增加，土方量的計(jì)算精度也就越高，在本文的算例中，當(dāng)樣本數(shù)量大于等于120時(shí)，與精確解相比，計(jì)算誤差小于5 %。

(2)與傳統(tǒng)的方格網(wǎng)法相比，本文中提出的方法計(jì)算精度更高，在相同的計(jì)算精度下，所需的外業(yè)工作量更少。