基于Excel的直剪試驗數據處理方法

莫林利，趙秀紹，王 旭，鄭 偉

(1.華東交通大學 軟件學院，江西 南昌 330013;2.華東交通大學 土木建筑學院，江西 南昌 330013)

基于Excel的直剪試驗數據處理方法

莫林利1，趙秀紹2，王 旭2，鄭 偉2

(1.華東交通大學 軟件學院，江西 南昌 330013;2.華東交通大學 土木建筑學院，江西 南昌 330013)

回歸分析中的最小二乘法是處理直剪試驗數據誤差最小的方法，但計算工作量大且容易出錯。應用Excel軟件給出可用于直剪試驗計算的5種方法，可顯著減少計算工作量并降低出錯率。實際工程數據計算對比表明，除基于點線距的規劃求解法外，其余4種方法均為基于擬合誤差的計算方法，兩種原理的計算結果差距可忽略不計。對于無黏性土情況，除INTERCEPT函數與SLOPE函數聯合求解法無法實現強制黏聚力為0外，其余方法均可實現。實際應用表明Excel圖解法、TREND函數法及截距和斜率函數聯合求解法操作簡單，易于應用。

直剪試驗 Excel內部函數 最小二乘法 規劃求解

邊坡土的抗剪強度是鐵路路基邊坡穩定性定量評價和計算滑坡推力的重要參數，有時內摩擦角相差1°～2°，滑坡推力就可能成倍增加［1］。馬顯春等［2］采用正交分析法對滑坡穩定性的影響因子進行了分析，得出滑坡影響因素敏感性順序為內摩擦角＞黏聚力＞地下水位＞重度＞坡高。可見合理確定巖土材料的抗剪強度指標對于工程設計、施工等十分重要［3-5］。

直剪試驗的數據處理實際上是求出土的抗剪強度參數，即黏聚力 c和內摩擦角 φ。斜截式計算法是在試驗中測得N，M兩點的垂直壓應力σN，σM，其對應的抗剪強度為τN，τM。通過 N，M的直線方程可以獲得土的抗剪強度指標［6］。此法的不足之處是非常明顯的，計算時僅能使用2組數據?！锻凉ぴ囼灱夹g手冊》［7］建議的方法為繪制抗剪強度與垂直壓力關系曲線，直線的傾角為內摩擦角φ，縱坐標的截距為黏聚力c。但是在繪庫侖破壞直線時，一般是繪圖者根據經驗和直觀判斷，盡量讓落在直線兩邊的點數大致相等，然后求出截距和斜率。由于它是人為地根據經驗和直觀判斷得出的，c，φ值往往帶有一定的人為性和不確定性。采用最小二乘法對試驗數據進行擬合求解［8-11］，計算量大且容易出錯。本文探討如何巧用Excel的各項功能簡化計算過程，并考慮直剪試驗中c值為負或強制c=0(無黏性土)的情況，使數據處理達到全自動化。

1 最小二乘法求解原理

為了便于說明各種求解方法的原理與實現方法，以江西撫州八景鎮火車站路基工程一組土樣實測數據(見表1)為例進行分析。

表1 撫州八景鎮路基土樣測試數據

直剪試驗求解中，當繪制測點的破壞包線時，由于土樣不均勻，試驗讀數誤差等原因，往往繪制的直線不能通過所有的測點，總是存在一些誤差，根據這些誤差的平方和為最小，可以確定一條近似的破壞包線。設直剪試驗中每組n個土樣的垂直應力、抗剪強度分別為σi，τi，利用庫侖公式 τ=σtanφ+c進行擬合，其中c，φ被視為待定系數。以偏差+c)］2最小為條件來確定c，φ值。運用高等數學的多元函數求極值的方法即可得到當 M為最小時的 c，φ值。分別對c，φ求一階偏導并令其等于0得［8］

根據式(1)，式(2)解得

如果用上述數學法求解，計算工作量大且易出錯。Excel軟件中有豐富的回歸分析函數和規劃求解方法，只要挖掘其功能，就可以簡單有效地完成上述計算。

2 直剪試驗的Excel求解方法

2.1 Excel圖解法

Excel軟件有強大的繪圖工具，先把試驗數據填入數據表格中，然后插入圖表，并為數據添加趨勢線和顯示公式，就可以得到直剪試驗圖。具體實現方法為:選擇正應力和剪應力數據，選擇“插入”→“圖表”→“散點圖”，然后選擇圖形中的點后“添加趨勢線”，在“選項”欄中選擇“顯示公式”，可直接得到擬合公式 y= 0.430 6x+3.260 9，進而得 c=3.260 9 kPa，tanφ= 0.430 6，即φ=23.297 1°。

2.2 TREND函數求解法

2.2.1 TREND函數求解的基本原理

由于直剪試驗各測點存在線性關系，這樣就可以方便地利用 Excel中的線性插值函數 TREND。TREND函數用途:返回一條線性回歸擬合線的一組縱坐標值(y值)。即找到適合給定的數組known_y's和known_x's的直線，并返回指定數組new_x's值在直線上對應的y值。TREND函數的語法為TREND(known _y's，known_x's，new_x's，const)。其中known_y's為已知的τ=σtanφ+c=σk+c中的τi值集合，known_x's為已知關系τ=σk+c中可選的σi值的集合，new_x's為需要函數TREND返回對應新σ值的新τ值，const參數為邏輯值指明是否強制常數項c為0。當擬合發現c值為負時，由于黏聚力不可能為負值，因此可以強制c值為0。

TREND函數求解原理如圖1所示，很顯然，對于抗剪強度包線τ=σtanφ+c，當new_x's=σ=0時，其返回的函數值為τ0=c，即為黏聚力值。當new_x's= σ=1時，函數返回值 τ1=1×tanφ+c，根據圖1中的幾何關系可得tanφ=(τ1－τ0)/1=τ1－τ0。

2.2.2 TREND函數求解法的實現

將表1中的原始數據填入Excel表中，如圖2(屏幕截圖)中所示。

圖1 TREND函數求解原理

圖2 TREND函數求解抗剪強度指標示例

在C2和D2單元格中，分別輸入“=TREND(B2: B5，A2:A5，0)”，“=ATAN(TREND(B2:B5，A2:A5，1)－TREND(B2:B5，A2:A5，0))*180/PI()”，即可求得該土樣的黏聚力c和內摩擦角φ。

TREND(B2:B5，A2:A5，0)在σ=0時可得τ0=c; TREND(B2:B5，A2:A5，1)在 σ=1 kPa時可得 τ1=c +k;TREND(B2:B5，A2:A5，1)－TREND(B2:B5，A2: A5，0)為系數k。

當土樣為無黏性土或計算出的黏聚力c值為負值時，可以強制TREND()函數中const項參數為0。對表1的測試數據求解結果為:c=3.260 9 kPa，φ= 23.297 1°。

2.3 INTERCEPT函數與SLOPE函數聯合求解法

2.3.1 INTERCEPT函數與SLOPE函數聯合求解原理

在Excel軟件中，函數INTERCEPT和SLOPE利用現有的x值和y值計算直線與y軸的截距和斜率。其中INTERCEPT(known_y's，known_x's)函數可以利用現有的σ值和τ值計算直線與τ軸的截距，即抗剪強度指標黏聚力 c;SLOPE(known_y's，known_x's)函數可以利用現有的σ值和τ值計算直線的斜率 k，即抗剪包線斜率tanφ。known_y's表示因變的觀察值或數據集合;known_x's表示自變的觀察值或者數據集合。

2.3.2 INTERCEPT函數與 SLOPE函數聯合求解法的實現

設數據在 Excel表中的排列仍如圖2所示，則在C2和D2單元格中分別輸入“=INTERCEPT(B2:B5，A2:A5)”，“=ATAN(SLOPE(B2:B5，A2:A5))* 180/PI()”。求解出的結果為 c=3.260 9 kPa，φ= 23.297 1°。

2.4 基于點線距的規劃求解法

2.4.1 基本原理

如圖3所示，對抗剪強度庫侖公式τ=σtanφ+c，令tanφ=k，則任意測點 A(σi，τi)到直線的距離為di［12］，見式(5)。

圖3 點線距求解基本原理示意

根據最小二乘法的原理，回歸分析時就是令 d2i和達到最小，即

式(6)中有k和c兩個變量，實際為線性規劃求解問題，求解目標函數的最小值。Excel規劃求解一般步驟為建立目標函數f(x)min，確定約束條件將決策變量用數學公式表達出來，使需要優化的目標及評判標準量化為這組變量的函數。

2.4.2 直剪試驗的規劃求解實現

Excel規劃求解首先要設計好一個工作表(如圖4所示)，根據建立的數學模型，將決策變量、約束條件、目標函數等依次排列，并在工作表中將相關固定信息初始化，具體實現步驟如下。

圖4 基于點線距的規劃求解示例

1)打開Microsoft Excel軟件后，按圖4基于點線距的規劃求解所示的格式建立工作表，將有關數據輸入工作表中(其中C2和D2為空白，有待計算)，E2單元格輸入“=((C2*A2－B2+D2)/SQRT(C2^2+1)) ^2+((C2*A3－B3+D2)/SQRT(C2^2+1))^2+ ((C2*A4－B4+D2)/SQRT(C2^2+1))^2+((C2* A5－B5+D2)/SQRT(C2^2+1))^2”。

2)在“工具”欄下拉菜單中找到“規劃求解”，如果沒有此選項，可以選擇“加載宏”將“規劃求解”加載進來。單擊“規劃求解”，出現“規劃求解參數”對話框，設置目標單元格:$E$2，在“等于”選項中選“最小值”，在“可變單元格”中輸入$F$2:$G$2，在“約束”中單擊“添加”后，逐個輸入下式:“$C$2＞=0”，“$D$2＞=0”，最后按“確定”。當需要強制 c值為0時，則條件變為“$C$2＞=0”，“$D$2=0”。求解結果為c=3.254 6 kPa，k=0.430 6，在C3單元格中求得φ=23.298 8°，此時目標函數最小，f(x)min=5.618 9。

2.5 基于測試值與擬合值誤差的規劃求解

直剪試驗數據σi，τi利用庫侖公式擬合直剪試驗中每組n個土樣的試驗數據時，c，φ被視為待定系數。以偏差最小為條件從而確定c，k值。這仍為線性規劃求解問題，求解目標函數的最小值。其實現方法同2.4中方法類似，僅需要根據公式更改其目標函數為“=(B2－A2*C2－D2)^2+(B3－A3*C2－D2)^2+ (B4－A4*C2－D2)^2+(B5－A5*C2－D2)^2”。

運行規劃求解后，選擇“保存極限值”即可同時求得c值和k值，然后根據k值可計算出內摩擦角φ值。根據上述方法，表1數據的求解結果為c=3.261 0 kPa，φ=23.297 1°，此時目標函數最小，M=6.660 9。

3 各種直剪試驗求解方法的對比

為了對比5種求解方法的有效性，以表1中測試數據作為第1組數據，以福建省福清市江陰鎮海邊固結淤泥測試數據作為第2組數據，第2組數據如表2所示。

表2 福清江陰鎮淤泥土樣測試數據

通過在Excel中求解，5種方法求解結果如表3所示。根據表3對2組試驗的計算結果，分析如下:

1)常規數據計算結果對比。雖然各種方法均采用了最小二乘法，在不強制黏聚力為0的情況下，除基于點線距的Excel規劃求解法外，其余方法求解結果均相同，說明其余4種方法均是基于擬合誤差的求解方法。從數據對比可知5種方法計算結果非常接近，

均可作為直剪試驗的求解結果。

表3 Excel中5種直剪試驗求解方法計算結果對比

2)黏聚力為負值情況處理對比。表3中福州省福清市江陰鎮由于取樣撓動等原因，多組測試數據測試結果擬合后，黏聚力均為負值，這與理論并不相符。因此對于c值為負的情況或無黏性土的黏聚力為0的情況，此時需要強制黏聚力為0。對于強制c=0的情況，INTERCEPT與SLOPE函數聯合求解法無法完成強制c=0，其余方法均可強制c=0。通過對比可得，對于無黏性土或黏聚力出現負值情況，不宜采用INTERCEPT與SLOPE函數聯合求解法。

3)操作簡易程度對比。從5種方法在數據更換后Excel操作繁簡比較，規劃求解需要重新求解，特別是試驗組數更改后，目標函數應隨著更改，而且會生成中間報告，故對計算機操作不太熟悉人員有些難度。其它方法均操作簡便，且計算結果會自動完成更新。

4 結論

根據最小二乘法的基本原理及前人的研究成果，提出了基于Excel圖解法、INTERCEPT與SLOPE函數聯合求解法、TREND函數法、基于點線距的規劃求解法和基于擬合誤差的規劃求解法，重點分析了TREND函數與基于點線距規劃求解法2種方法的求解原理，并針對每種求解方法在Excel軟件中設計了相應公式或實現方法。通過幾種直剪試驗的求解方法對比分析，可以得出以下結論:

1)Excel中的 TREND函數求解法，截距函數INTERCEPT與斜率函數 SLOPE聯合求解法均為Excel內部函數法，大大地減小了編輯公式的工作量，是一種操作簡便計算方法，在無黏性土或計算黏聚力＜0的情況下，TREND函數求解法最為方便。

2)5種求解方法基本原理均采用基于最小二乘法的回歸分析，其擬合誤差保證其與實測值差值的平方和最小，其精度較采用兩點計算的斜截求解法和手工繪圖更能反映各測點對試驗結果的影響，更能代表土的實際抗剪強度指標。

3)在5種求解方法中，當計算的黏聚力 c值為正值時，各種方法求解結果相近，均可作為試驗計算結果。當由于某種原因而引起擬合c值為負值或土樣為無黏性土時，c=0才正確。Excel圖解法、TREND函數法、兩種規劃求解法均可實現強制黏聚力為0，而INTERCEPT與SLOPE函數聯合求解法無法實現強制截距為0，故不能采用。

［1］李向國.滑帶土抗剪強度參數的確定［J］.西部探礦工程，2010(2):37-38.

［2］馬顯春，王雷，趙法鎖.滑坡穩定影響因子敏感性分析及治理方案探討［J］.地質力學學報，2008，14(4):381-388.

［3］劉遠明，劉杰，夏才初.不同節理表面形貌下非貫通節理巖體強度特性直剪試驗研究［J］.巖土力學，2014，35(5): 1269-1273.

［4］衡帥，楊春和，曾義金，等.基于直剪試驗的頁巖強度各向異性研究［J］.巖石力學與工程學報，2014，33(5):874-883.

［5］陳建勝，陳從新，魯祖德，等.強風化角巖力學—變形特性的直剪試驗研究［J］.巖土力學，2010，31(9):2869-2874.

［6］王曉斌，張鐵男，董志國.直剪試驗中的斜截式計算方法［J］.水利科技與經濟，1999，5(2):82-83.

［7］南京水利科學研究院土工研究所.土工試驗技術手冊［M］.北京:人民交通出版社，2003.

［8］談云志，王世梅.直剪試驗數據處理方法的對比分析［J］.三峽大學學報(自然科學版)，2005，27(2):132-134.

［9］林立軍.二次曲線擬合法在無機結合料穩定材料擊實試驗數據處理中的應用［J］.福建建筑，2012，30(5):83-85.

［10］余東明，姚海林，吳少鋒.三軸試驗抗剪強度參數值回歸分析法的區別與修正［J］.巖土力學，2012，33(10):3037-3042.

［11］陳輝，顏勝才，孟軍濤.K30試驗數據的曲線擬合分析方法［J］.鐵道建筑，2012(11):67-70.

［12］阮波，張向京，彭意.Excel規劃求解三軸試驗抗剪強度指標［J］.鐵道科學與工程學報，2009，6(5):57-60.

(責任審編 李付軍)

TU411

:ADOI:10.3969/j.issn.1003-1995.2015.09.29

2014-11-03;

:2015-02-10

國家自然科學基金項目(51208198);江西省教育廳項目(JXJG-13-5-15);華東交通大學科研基金(12RJ02)

莫林利(1977— )，女，河南安陽人，講師，碩士。

1003-1995(2015)09-0102-04