以分離元素法探討板岩邊坡變形機制

全文

(1)國立高雄大學土木與環境工程學系碩士論文. 以分離元素法探討板岩邊坡變形機制 A study on the deformational mechanism of slate slopes using discrete element method. 研究生：莊庭鳳撰指導教授：翁孟嘉博士. 中華民國 103 年 09 月.

(2) 誌謝首先誠摯的感謝翁孟嘉老師的指導，不時地指導並指點我正確的方向，使我在這些年中獲益匪淺，成功完成我的論文。這兩年多來，感謝實驗室的同伴與學長姊、學弟妹的幫忙與陪伴。感謝洪銘謙、蔡尚豪學長不厭其煩地指出我研究中的缺失，更在我迷惘時為我解惑。也感謝蔡承翰在我電腦出問題時，替我解決問題，感謝蘇威豪在我需要幫忙時，能熱心的幫助我，讓我順利度過這兩年。當然也不會忘記實驗室的吳政賢學弟、蘇育萱學妹，感謝你們辛苦的幫忙黏試體與幫忙做實驗，讓我省了很多事，真的非常感謝你們。最後，我要將這篇論文獻給我的家人，感謝他們一直以來的支持，讓我能專注於課業與研究。.

(3) 目錄第一章緒論 ..................................................................................................... 1 1.1 研究動機 ............................................................................................. 1 1.2 研究目的 .............................................................................................. 2 1.3 研究方法及範疇 ................................................................................. 2 第二章文獻回顧 ............................................................................................. 4 2.1 臺灣板岩 ............................................................................................. 4 2.2 邊坡變形型態 ..................................................................................... 6 2.3 PFC3D 理論介紹 ................................................................................. 13 2.4 數值模擬研究 ................................................................................... 28 第三章研究方法 ........................................................................................... 30 3.1 物理模型試驗 ................................................................................... 30 3.2 數值模擬分析 ................................................................................... 40 第四章物理試驗暨數值模擬驗證成果 ....................................................... 43 4.1 物理模型試驗 ................................................................................... 43 4.2 數值模擬物理試驗 ........................................................................... 51 4.3 小結 ................................................................................................... 60 第五章參數之影響 ....................................................................................... 61 5.1 摩擦係數 ........................................................................................... 61 I.

(4) 5.2 球勁度 ............................................................................................... 66 5.3 平行鍵結勁度 ................................................................................... 70 5.4 平行鍵結鍵結強度 ........................................................................... 72 第六章坡度及葉理角度之影響 ................................................................... 74 6.1 模型建置 ........................................................................................... 74 6.2 探討同坡度、不同葉理情形下之影響 ........................................... 78 6.3 探討同葉理、不同坡度情形下之影響 ........................................... 89 第七章 Burger’s model 模擬 ....................................................................... 100 7.1 探討依據 .......................................................................................... 100 7.2 Viscosity for Maxwell section 之影響 ............................................. 101 7.3 Stiffness for Maxwell section 之影響.............................................. 103 7.4 Viscosity for Kelvin section 之影響 ................................................ 105 7.5 Stiffness for Kelvin section 之影響 ................................................. 107 第八章新人崗現地模擬 ............................................................................. 109 8.1 研究區域之基本介紹 ..................................................................... 109 8.2 模型設定 ......................................................................................... 110 8.3 成果探討 .......................................................................................... 113 第九章結論與建議 ..................................................................................... 115 參考文獻 ....................................................................................................... 117 II.

(5) 圖目錄圖 1 臺灣大部分板岩分布區域 ....................................................................... 4 圖 2 現地調查點位 ........................................................................................... 7 圖 3 高傾角葉理角度邊坡 ............................................................................... 7 圖 4 低傾角葉理角度邊坡 ............................................................................... 8 圖 5 岩石潛變四大模態(Chigira, 1992) .......................................................... 9 圖 6 佳陽地區劈理位態分布圖(邱, 2000) .................................................... 10 圖 7 高屏溪現地調查成果圖(陳等, 2012) .................................................... 11 圖 8 大雪山林道邊坡所見之板岩重力變形(李, 2012) ................................ 12 圖 9 PFC3D 運算流程圖(PFC3D manual) ...................................................... 14 圖 10 顆粒與顆粒接觸之關係圖(PFC3D manual) ......................................... 15 圖 11 顆粒與牆接觸之關係圖(PFC3D manual) ............................................. 16 圖 12 顆粒與牆接觸之法線向量舉例說明圖(PFC3D manual) ..................... 17 圖 13 接觸鍵結模式行為 ............................................................................... 22 圖 14 平行鍵結模式示意圖(PFC3D manual) ................................................. 22 圖 15 PFC3D 之 Burger’s model 示意圖(PFC3D manual) .............................. 25 圖 16 研磨石 ................................................................................................... 31 圖 17 片狀試體 ............................................................................................... 31 圖 18 研磨石之直接剪力實驗結果 ............................................................... 31 III.

(6) 圖 19 研磨石之破壞包絡線 ........................................................................... 32 圖 20 片狀試體之直接剪力實驗結果 ........................................................... 32 圖 21 片狀試體之破壞包絡線 ....................................................................... 33 圖 22 樹脂之拉應力 ....................................................................................... 33 圖 23 物理試驗模型 ....................................................................................... 34 圖 24 物理試驗箱，(a)乾燥狀態，(b)泡水弱化狀態 ................................. 35 圖 25 各參數物理意義示意圖 ....................................................................... 36 圖 26 試體製備步驟 ....................................................................................... 38 圖 27 片狀試體依序排列情形 ....................................................................... 39 圖 28 數值模擬模型範例圖 ........................................................................... 41 圖 29 DT1B3-S30FO90 試驗過程，(a)試驗前，(b)經過 30 秒 .................. 43 圖 30 DT1B3-S30FO90 試驗過程，(c)1 分鐘後，(d)5 分鐘後 .................. 44 圖 31 DT1B3-FO90，(a)坡度 17 度，(b)坡度 30 度，(c)坡度 45 度 ........ 44 圖 32 DT1-S30FO90 試驗成果，(a)寬度 3 顆，(b)寬度 12 顆 .................. 45 圖 33 DT1B3-S30 試驗成果 .......................................................................... 45 圖 34 試體斷面圖 .......................................................................................... 46 圖 35 DB12-S30FO90 試驗成果，(a)弱面間距 1 層，(b)弱面間距 2 層 .. 46 圖 36 WT1B3-S30FC60 試驗過程，(a)試驗前，(b)經過 1 小時 ............... 47 圖 37 WT1B3-S30FC60 試驗過程，(c)經過 2 小時，(d)經過 40 小時 ..... 48 IV.

(7) 圖 38 WT1-S30FC60 試驗成果，(a)寬度 3 顆，(b)寬度 12 顆 .................. 48 圖 39 WT1B3-S17FC60 試驗成果，(a)側翻前，(b)側翻後 ....................... 49 圖 40 WT1B3-S45FC60 試驗成果，(a)側翻前，(b)側翻後 ....................... 49 圖 41 WT1B12-FC60 試驗成果(由左至右分別為坡度 17、30、45 度) .... 50 圖 42 試驗成果，(a)WT1B3-S30FC60，(b)WT1B3-S30FO90 ................. 50 圖 43 DT1B3-S17FC60 模擬成果 ................................................................. 52 圖 44 DT1B3-S30FC60 模擬成果 ................................................................. 52 圖 45 DT1B3-S45FC60 模擬成果 ................................................................. 52 圖 46 DT1B12-S17FC60 模擬成果 ............................................................... 53 圖 47 DT1B12-S30FC60 模擬成果 ............................................................... 53 圖 48 DT1B12-S45FC60 模擬成果 ............................................................... 53 圖 49 DT1B3-S17FO85 模擬成果 ................................................................. 55 圖 50 DT1B3-S30FO85 模擬成果 ................................................................. 55 圖 51 DT1B3-S45FO85 模擬成果 ................................................................. 55 圖 52 DT1B12-S17FO85 模擬成果 ............................................................... 56 圖 53 DT1B12-S30FO85 模擬成果 ............................................................... 56 圖 54 DT1B12-S45FO85 模擬成果 ............................................................... 56 圖 55 DT1B3-S17FO60 模擬成果 ................................................................. 57 圖 56 DT1B3-S30FO60 模擬成果 ................................................................. 57 V.

(8) 圖 57 DT1B3-S45FO60 模擬成果 ................................................................. 58 圖 58 DT1B12-S17FO60 模擬成果 ............................................................... 58 圖 59 DT1B12-S30FO60 模擬成果 ............................................................... 58 圖 60 DT1B12-S45FO60 模擬成果 ............................................................... 59 圖 61 WT1B3-S30FC60 之模擬成果 ............................................................ 59 圖 62 不同參數下影響行為之探討因子 ....................................................... 64 圖 63 逆向坡模型改變摩擦係數之 h/H、l/L 比值圖 .................................. 64 圖 64 逆向坡模型改變摩擦係數之變形角度 β 比較圖............................... 65 圖 65 逆向坡模型改變球勁度之 h/H、l/L 比值圖 ...................................... 68 圖 66 逆向坡模型改變球勁度之變形角度 β 比較圖................................... 69 圖 67 探討因子示意圖，(a)泡水弱化情形，(b)乾式情形 ......................... 77 圖 68 坡度 17 度各葉理角度之 t/T、h/H 比值 ............................................ 79 圖 69 坡度 30 度各葉理角度之 t/T、h/H 比值 ............................................ 80 圖 70 坡度 45 度各葉理角度之 t/T、h/H 比值 ............................................ 80 圖 71 坡度 60 度各葉理角度之 t/T、h/H 比值 ............................................ 81 圖 72 各坡度各葉理角度之 t/T 比較圖 ........................................................ 82 圖 73 各坡度各葉理角度之 h/H 比較圖 ....................................................... 83 圖 74 各坡度各葉理角度之變形角度 α 比較圖........................................... 83 圖 75 坡度 17 度各葉理角度之 h/H、l/L 比值 ............................................ 85 VI.

(9) 圖 76 坡度 30 度各葉理角度之 h/H、l/L 比值 ............................................ 86 圖 77 坡度 45 度各葉理角度之 h/H、l/L 比值 ............................................ 86 圖 78 坡度 60 度各葉理角度之 h/H、l/L 比值 ............................................ 87 圖 79 各坡度各葉理角度之 h/H 比較圖 ....................................................... 88 圖 80 各坡度各葉理角度之 l/L 比較圖 ........................................................ 88 圖 81 各坡度各葉理角度之變形角度 β 比較圖........................................... 89 圖 82 葉理 60 度各坡度之 t/T、h/H 比值 .................................................... 90 圖 83 葉理 70 度各坡度之 t/T、h/H 比值 .................................................... 91 圖 84 葉理 80 度各坡度之 t/T、h/H 比值 .................................................... 91 圖 85 各葉理角度各坡度之 t/T 比較圖 ........................................................ 92 圖 86 各葉理角度各坡度之 h/H 比較圖 ....................................................... 93 圖 87 各葉理角度各坡度之變形角度 α 比較圖........................................... 94 圖 88 葉理逆向 60 度各坡度之 h/H、l/L 比值 ............................................ 95 圖 89 葉理逆向 70 度各坡度之 h/H、l/L 比值 ............................................ 96 圖 90 葉理逆向 80 度各坡度之 h/H、l/L 比值 ............................................ 96 圖 91 葉理逆向 85 度各坡度之 h/H、l/L 比值 ............................................ 97 圖 92 各葉理角度各坡度之 h/H 比較圖 ....................................................... 98 圖 93 各葉理角度各坡度之 l/L 比較圖 ........................................................ 98 圖 94 各葉理角度各坡度之變形角度 β 比較圖........................................... 99 VII.

(10) 圖 95 探討依據示意圖 ................................................................................. 100 圖 96 Viscosity for Maxwell section 之 h/H、l/L 比值 ............................... 102 圖 97 Viscosity for Maxwell section 之變形角度比較圖............................ 102 圖 98 Stiffness for Maxwell section 參數更改成果之 h/H、l/L 比值........ 104 圖 99 Stiffness for Maxwell section 參數更改成果之變形角度比較圖 .... 104 圖 100Viscosity for Kelvin section 參數更改成果之 h/H、l/L 比值 ......... 106 圖 101Viscosity for Kelvin section 參數更改成果之變形角度比較圖 ...... 106 圖 102 Stiffness for Kelvin section 參數更改成果之 h/H、l/L 比值 ......... 108 圖 103 Stiffness for Kelvin section 參數更改成果之變形角度比較圖...... 108 圖 104 現地模擬調查點位示意圖 ............................................................... 109 圖 105 調查區域現地照片(Lo and Feng, 2014) .......................................... 110 圖 106 新人崗之現地模擬模型 (Lo and Feng, 2014) ................................ 111 圖 107 本研究之新人崗現地模擬模型示意圖 ........................................... 111 圖 108 本研究之新人崗現地模擬模型示意圖.......................................... 114. VIII.

(11) 表目錄表 1 臺灣雪山山脈帶北部及中央脊樑山脈帶地層分類及對比表 .............. 5 表 2 試驗規劃表 ............................................................................................. 39 表 3 基本參數單位表 ..................................................................................... 40 表 4 顆粒參數設定 ......................................................................................... 42 表 5 順向坡模型改變摩擦係數 ..................................................................... 62 表 6 逆向坡模型改變摩擦係數 ..................................................................... 63 表 7 順向坡模型改變球的勁度 ..................................................................... 66 表 8 逆向坡模型改變球的勁度 ..................................................................... 67 表 9 順向坡模型改變平行鍵結勁度 ............................................................. 70 表 10 逆向坡模型改變平行鍵結勁度 ........................................................... 71 表 11 順向坡模型改變平行鍵結鍵結強度 ................................................... 72 表 12 逆向坡模型改變平行鍵結鍵結強度 ................................................... 73 表 13 乾式情況模擬成果 ............................................................................... 75 表 14 泡水弱化情況模擬成果 ....................................................................... 76 表 15 Viscosity for Maxwell section ............................................................. 101 表 16 Stiffness for Maxwell section 之成果................................................. 103 表 17 Viscosity for Kelvin section 之成果 ................................................... 105 表 18 Stiffness for Kelvin section 之成果 .................................................... 107 IX.

(12) 表 19 顆粒參數設定 ..................................................................................... 112. X.

(13) 以分離元素法探討板岩邊坡變形機制指導教授：翁孟嘉博士國立高雄大學土木與環境工程所學生：莊庭鳳國立高雄大學土木與環境工程所摘要有鑑於現地觀察到之板岩邊坡重力變形型態，再加上近年來高強度、長延時之降雨頻率漸增，導致許多坡地災害的發生，突顯出板岩邊坡崩塌發生潛勢之事前判斷與評估的重要性。為了模擬板岩邊坡之重力變形行為，本研究將現地所觀察到的邊坡特性簡化至室內物理模型試驗，並探討三種影響因素下之運動行為，分別為邊坡角度、葉理角度與岩性弱化。再利用物理模型試驗中各力學性質與所得到之重力變形行為，代入以分離元素法為理論基礎之數值分析軟體 PFC3D 進行模擬。綜合以上試驗與模擬結果，結論如下：(1)由物理實驗得知，試體在變形過程中主要分為兩大類，翻覆破壞與褶皺破壞。翻覆破壞主要分布於高葉理角度之逆向坡；褶皺破壞主要發生於低葉理角度之順向坡，這種破壞形式不只因邊坡受到重力的影響，更因遭受水的弱化，於坡趾處產生折彎變形。 (2) 根據數值模擬所得出的模擬成果與物理實驗成果進行比對，確立此種模擬方法之可行性，並整理出邊坡在不同角度、葉理角度、弱化下之變形型態，供日後滑動分析參考。關鍵字：分離元素法、重力潛移、板岩. XI.

(14) A study on the deformational mechanism of slate slopes using discrete element method Advisor: Dr. Meng-Chia Weng Institute of Civil and Environmental Engineering National University of Kaohsiung Student: Ting-Feng Chuang Institute of Civil and Environmental Engineering National University of Kaohsiung ABSTRACT According to the gravitational deformation of slate slope in the field, this study explored the influencing factors,including the slope angle, the foliation angle, and the material deterioration on the deformation characteristics. To simulate these phenomena, a series of physical model experiments with simplified environmental conditions were performed in the laboratory. Then, the discrete element software, PFC3D was used for simulations to further our understanding of deformation behavior in slate slopes under the long-term influence of gravity and material deterioration. According to the physical model tests and numerical simulations the results possesses the following characteristics:(1) Based on the model tests, two kinds of gravitational deformation of slate slope are identified, the first one is overturning destruction, it appears in the obsequent slope with high angle foliation. The other one is fold extrusion near the slope toe, it exhibits in the consequent slope with low angle foliations. It is not only influenced by gravitation but also by wet deterioration.(2) The results of physical model test provide the basic case for the verification of discrete element simulation. The simulated deformed patterns agree well with the actual deformation under different slope angles. Therefore, the model of discrete element method and associated parameters are verified, and they could be used for the subsequent analysis in the future.. Keywords: discrete element method, gravitational deformation, slate. XII.

(15) 第一章緒論 1.1 研究動機 2009 年莫拉克風災重創臺灣，風災過後，高屏溪上游出現為數不少之深層崩塌或具發生潛勢之邊坡，造成多處坡地災害。此外，廬山地區的地滑情形亦趨嚴重。Chigira (2011)於莫拉克颱風後曾提及，臺灣近中央山脈之板岩區邊坡常有大規模潛移(creep)之特徵，造成深層滑動面之形成，當極端豪雨事件發生時，水易滲入岩坡弱面內，導致深層滑動之發生。李錫堤(2012)亦提出板岩區存在大規模山坡潛移為臺灣常見而一般瞭解不多或容易被誤解的困難地質條件之一，故需特別注意其工程地質問題。承上述所言，再加上近年來隨著全球氣候變遷而導致高強度、長延時之降雨頻率漸增，邊坡崩塌之可能性也隨之增加，一再突顯板岩邊坡崩塌發生潛勢之事前判斷與評估的重要性。但目前對於板岩邊坡之變形機制、可能運動過程等相關模擬技術未臻成熟，實有深入研究及發展之必要。為模擬板岩邊坡之變形行為，考量到利用分離元素法(discrete element method)進行模擬時，具可將滑動塊體從變形、破裂到分離的運動過程完整呈現之特性，故選用以分離元素法為理論基礎之模擬軟體 (PFC3D)進行模擬。 1.

(16) 有鑑於上述所言，本研究為符合 PDC3D 之球元素，並考量到板岩葉理發達之特性，故以研磨石膠結為片狀試體之方式來模擬板岩之變形行為。並設計出不同葉理角度與坡度之物理試驗，藉此觀察僅受重力影響下之岩體潛變行為，最後以 PFC3D 進行模擬，與物理試驗成果做驗證，確立此分析方法之可行性。. 1.2 研究目的本研究之主要目的在於： 1. 簡化板岩邊坡變形機制，應用室內物理模型試驗，求得試體僅於重力影響下之試體變形行為。 2. 利用以分離元素法為理論基礎之模擬軟體(PFC3D)，針對不同寬度與弱面間距之試體、坡度及葉理角度等條件下之潛變行為加以驗證，藉此確立此分析方法之可行性。 3. 不同力學行為下，僅在重力影響下之試體變形行為。. 1.3 研究方法及範疇本研究為求得試體之變形行為，簡化板岩邊坡變形機制後，設計一組室內物理模型試驗箱，針對不同寬度、弱面間距、坡度、葉理角度及乾燥或泡水弱化狀態下之試體進行試驗，了解僅在重力影響下之試體變形行為。. 2.

(17) 接下來採用以分離元素法為理論基礎之數值軟體 PFC3D (Particle Flow Code in 3 Dimensions)進行模擬。俟驗證物理試驗結果後，即確立模擬方法，便可沿用模擬物理試驗時之各參數或模型，探討參數之影響及不同坡度與葉理角度情況下之影響。除此之外，另以 Burger’s model 進行模擬，探討不同參數之影響。最後，選定一實際邊坡，進行現地模擬，進一步驗證模擬方法之可行性。. 3.

(18) 第二章文獻回顧本章分別針對邊坡變形型態及數值分析等兩個部分做基本介紹。. 2.1 臺灣板岩板岩，為一變質度低之變質岩。具有極密且良好的葉理狀岩理，且極易順其葉理方向裂開，稱為板劈理，為異向性岩石。位於臺灣本島的板岩層 (slate formation)約佔三分之一或山地之一半面積，大部分板岩分布區域如圖 1 所示(李, 2012)。根據台灣地質概論，以下針對主要由板岩所構成之地層(即標示於表 1 中粗體字之地層)進行描述(何, 1997)。. 圖 1 臺灣大部分板岩分布區域. 4.

(19) 表 1 臺灣雪山山脈帶北部及中央脊樑山脈帶地層分類及對比表地質時代. 雪山山脈帶北部. 中央脊樑山脈帶廬山層. 中新世. 蘇樂層禮觀層. 中新-漸新世. 澳底層大桶山層 (粗窟砂岩) 地層間斷. 漸新世. 乾溝層四稜砂岩. 漸新-始新世. 西村層. 畢祿山層(新高層). 始新世西村層此地層主要由葉理發達之深灰色板岩和千枚岩質板岩所構成，夾有暗灰色、中粗粒、堅硬的石英砂岩互層，這種互層常見於此層的下部。畢祿山層此地層主要由板岩和千枚岩所構成，但是在變質的泥質岩層中夾有較厚的變質砂岩層，有的為石灰質砂岩，有的為長石質砂岩，砂岩粒度由細粒至粗粒，層厚也有薄有厚。. 5.

(20) 廬山層此層大部分由黑色到深灰色的硬頁岩、板岩及千枚岩和深灰色的硬砂岩互層組成，含有零星散布的泥灰岩團塊，其全部厚度估計在數千公尺上下，主要分布在脊樑山脈帶中。. 2.2 邊坡變形型態為判斷板岩邊坡之露頭是否因重力而產生變形，本研究之方法為先找出板岩層面，再找出葉理角度不連續處，且位於邊坡折彎處之岩層易產生破碎，並可藉由植生判斷，上述所提及之破碎處有益於植生生長。藉由上述之判斷方法，為觀察板岩邊坡可能的變形型態，本研究至南部的板岩地區進行現地調查，其中發現雖同為順向坡邊坡，但有不同的邊坡型態，調查點位如圖 2 所示，圖 3 為高傾角葉理角度之邊坡，其邊坡之型態有上方邊坡因重力往下擠壓，因此產生 S 型的折彎變形(圖 3(a))或波浪狀的折彎褶皺(圖 3(b))，但於圖 3(c)中顯示，雖然同樣為高葉理角度之邊坡，但也有可能會沒產生變形；圖 4 為低傾角葉理角度之邊坡，其邊坡之型態有隨著岩層弱化，邊坡產生崩壞並滑落(圖 4(a))、或是上方岩層沿層面滑落，但下方岩層因底部的束限，故產生折彎變形(圖 4(b))。. 6.

(21) 圖 2 現地調查點位. (a) 1 號點位. (b) 2 號點位圖 3 高傾角葉理角度邊坡 7.

(22) (c) 圖 3 高傾角葉理角度邊坡(3 號點位). (a) 4 號點位. (b) 5 號點位圖 4 低傾角葉理角度邊坡. 8.

(23) 根據 Chigira(1992)對長時期受重力影響之岩體潛移作出的討論，總結出岩體受重力變形的潛變作用大致分為四大類型，如圖 5 所示，並將其分述如下： Ⅰ. Ⅱ. Ⅲ. Ⅳ. 圖 5 岩石潛變四大模態(Chigira, 1992) 1. 順向坡的折彎褶皺型(Ⅰ)：發生在坡度比葉理傾角平緩且為順向坡時，岩層受到底層基岩的束限無法活動，順葉理下滑時下邊坡將使坡趾處產生折彎摺皺的潛移現象。 2. 逆向坡的拖曳褶皺型(Ⅱ)：發生在坡向為逆向坡時，葉理位態呈現高角度，岩層中葉理產生上下層如同被拖曳但尚未切斷「S」型褶皺的潛移現象。. 9.

(24) 3. 高傾角葉理的折彎褶皺型(Ⅲ)：發生在葉理位態為高角度時，邊坡上層隨重力長期作用下，葉理當中有一折點，折點上葉理平行向下或向外產生折彎，剖面的葉理型態呈現「ㄑ」型褶皺的潛移現象。 4. 順向坡的斷層型(Ⅳ)：發生在葉理位態緩和且與坡面趨於平行，且邊坡上產生類似斷層之潛移面。邱坤豪(2000)調查中橫公路德基至梨山一帶的板岩層劈理位態，發現主要有三群劈理位態，其中兩種位態出露地區有特定坡向及坡度，與坡向及坡度的關連性密切，因此推論其劈理為偏離過的位態，是受到高山地區山坡大規模潛移作用的影響，並評估其符合上述 Chigira 所分類的第三型高傾角葉理的折彎褶皺型（Bending folds）潛移模式，如圖 6 所示。. 圖 6 佳陽地區劈理位態分布圖(邱, 2000). 10.

(25) 陳天健等人(2012)以臺灣南部高屏溪流域為研究區域，針對 42 處現地調查之大型崩塌地做現地調查，結果發現高屏溪流域之岩體滑動模式有四種主要破壞型態，分別為潛移、傾覆、大面積崩塌、構造，其中潛移佔 17 筆調查筆數，傾覆佔 6 筆調查筆數。圖 7 為此研究現地調查照片，其左圖符合 Chigira 所分類的第一型順向坡的折彎褶皺型；右圖符合第三型高傾角葉理的折彎褶皺型。. 圖 7 高屏溪現地調查成果圖(陳等, 2012) (左圖為露頭受重力作用產生潛變現象；右圖為折彎褶皺型傾覆) 李錫提(2012)提出臺灣常見而一般瞭解不多或容易被誤解的困難地質條件有三種，分別為：順向坡內部普遍會含有泥層的問題、斷層泥極為敏感的工程特性及其影響、板岩層區普遍已發生大規模山坡潛移產生之工程地質問題。板岩層的關鍵問題就是板岩受重力作用會產生撓曲的現象。圖 8 是大雪山林道邊坡所見之板岩重力變形情形。. 11.

(26) 圖 8 大雪山林道邊坡所見之板岩重力變形(李, 2012) 綜合以上現地邊坡變形行為成果，發現臺灣山地區域許多板岩層之變形型態符合 Chigira 所分類的第一型順向坡的折彎褶皺型及第三型高傾角葉理的折彎褶皺型。. 12.

(27) 2.3 PFC3D 理論介紹為了解塊體運動之力學行為，本研究採用以分離元素法(Discrete Element method, DEM)作為理論基礎之顆粒流分析軟體，名為 Particle Flow Code in 3 Dimensions (簡稱 PFC3D)，是由 Itasca Consulting Group, Inc. 所發展之軟體，其就空間維度有分為 PFC2D 以及 PFC3D 兩種，本研究使用的軟體為 PFC3D，底下根據 PFC3D 之手冊做一理論介紹。. 2.3.1 PFC3D 程式簡介 PFC3D 主要模擬元素有兩種，分別為球元素與牆元素，牆元素主要用於建立邊界，球元素主要用於透過圓形顆粒間之交互作用計算、分析，可模擬靜態或動態問題。除了主要功能之外，其它功能如下述： 1.. PFC3D 允許顆粒間鍵結以組合成任意形狀之物體，或是用指令將. 顆粒組合成塊體。 2.. 由鍵結所產生之塊體，當所承受之應力超過其鍵結強度時，允許. 鍵結產生斷裂，使得塊體破裂。 3.. 除了內建指令外，程式提供 FISH 語言供使用者撰寫，使用者可. 加入自行定義之模型以因應不同的需要。. 2.3.2 PFC3D 基本假設 PFC3D 內有以下六大基本假設： 1. 顆粒視為剛體。. 13.

(28) 2. 顆粒間的接觸面積趨近於零(點接觸)。 3. 顆粒間的接觸行為屬於軟接觸(soft contact)，且允許顆粒在接觸點上產生重疊。 4. 重疊量的大小與由力－位移定律所決定之接觸力有關，且重疊量小於顆粒粒徑。 5. 顆粒間的接觸可存在鍵結。 6. 所有顆粒均為圓球狀，且可利用指令將許多顆粒組合成任意形狀之物體。. 2.3.3 運算原理 PFC3D 在進行運算時需要：每個顆粒皆利用運動方程式重複計算、每個接觸力採用力-位移定律來得出、以及不斷更新牆的最新位置，其運算流程圖如圖 9 所示。關於力-位移定律與運動方程式之敘述茲將分別說明於後續小節中。. 圖 9 PFC3D 運算流程圖(PFC3D manual) 14.

(29) 2.3.3.1 力-位移定律在此將說明接觸之判斷，以及力量的計算。首先，接觸行為分別有顆粒與顆粒接觸(ball-ball)、顆粒與牆接觸(ball-wall)兩種，其示意圖如圖 10 及圖 11 所示。ni 表示單位法線向量，且兩者接觸行為於法線向量的定義上有所不同，顆粒與顆粒接觸之法線向量為兩球心之連線方向；球與牆接觸之法線向量為球心與牆最短距離方向之連線。. 圖 10 顆粒與顆粒接觸之關係圖(PFC3D manual). 15.

(30) 圖 11 顆粒與牆接觸之關係圖(PFC3D manual). 顆粒與顆粒接觸時(圖 10)，A 球與 B 球球心的位置向量分別為 xi A 與 xi  ，接觸平面上之單位法線向量 ni 為： B. xi   xi  ni  (ball-ball) d B. A. (2. 1). 兩球心之距離 d 為：. d  xi   xi B. A. .  x   x    x   x   (ball-ball) B. i. A. B. i. i. A. i. (2. 2). 顆粒與牆接觸時，其法線向量為球心與牆最短距離方向之連線，比如說目前有 AB 及 BC 兩面牆，如圖 12 所示，若球與牆接觸時球心位置位於第 2 或第 4 區則單位法線向量 ni 將垂直牆面，若球心位置位於第 1、 3、5 區則 ni 方向為球心與牆端點之連線方向。. 16.

(31) 圖 12 顆粒與牆接觸之法線向量舉例說明圖(PFC3D manual). 定義重疊量 Un 為法線方向之相對接觸位移為：  A  B  R  R  d U   b   R d n. (ball-ball) (ball-wall). (2. 3). 其中，R[Φ]為Φ球之半徑。接觸點的位置為：   A  A 1 n   xi   R  2 U  ni   C   xi    x b    R  b   1 U n  n  i  i  2 . (ball-ball) (2. 4) (ball-wall). 代表顆粒與顆粒接觸行為及顆粒與牆接觸行為之接觸力向量 Fi，可以分為正向力向量 Fi n 與剪力向量 Fi s ：. Fi n  K nU n ni. 17. (2. 5).

(32) 2.3.3.2 運動方程式剛體顆粒之運動方式由合力與合力矩決定，在顆粒中分別以質點的平移運動與顆粒的轉動表示。質點平移運動中，以 xi 、 xi 、 xi 分別表示位置、速度、加速度，而顆粒的轉動中，分別以 i 與 i 表示角速度與角加速度。運動方程式由兩種方程式所構成，分別為合力與平移運動之關係，以及合力矩與轉動之關係。平移運動方程式為：. Fi  m( xi  gi ) (translational motion). (2. 6). 其中， Fi 為合力；m 為顆粒總質量； g i 為重力加速度。轉動方程式為：. M i  Hi. (2. 7). 其中， M i 為合力矩， H i 為角動量。上述關係式皆為局部座標系統，座標系統連接到顆粒質量中心，若局部座標糸統指向顆粒之主軸時，則 2.7 式可簡化成 Euler 運動方程式：. M1 =I11 +  I3  I2 32 M 2 =I22 +  I1  I3 13. (2. 8). M 3 =I33 +  I2  I1 21. 其中，I1、I 2 及 I3 為顆粒主慣性矩；1、2 及 3 為主軸角加速度；M1 、 M 2 及 M 3 為合力矩於主軸方向分量。. 18.

(33) 假設質量均勻分佈於顆粒內部，對於半徑 R 之圓球顆粒，其質心視同球心。就三維顆粒而言，通過質量中心之任何局部座標系統皆為主軸系統，三方向之主慣性矩皆相同，故 2.7 式可簡化為：. 2  M i =Ii   mR 2  i (rotation motion) 5 . (2. 9). 藉由平移運動方程式(式 2.6)與轉動方程式(式 2.9)，並利用中間有限差分法進行數值計算，當 t  n t / 2 時，得出 xi 、 i ；當 t  n t / 2 時，得出 xi 、 xi 、 i 、 Fi 。在 t 時的平移與轉動加速度為：. 1 (t  t /2) xi  xi(t  t /2)   t 1 i(t ) = i(t  t /2)  i(t  t /2)  t xi(t ) =. (2. 10). 將 2.10 式分別代入 2.6 式與 2.9 式可得 t  t / 2 時的速度： ( t  t /2) i. x. i(t . ( t  t /2) i. =x. t /2). =i(t .  Fi (t )    gi  t  m   M i(t )  t /2)   t I  . (2. 11). 最後，使用 2.11 式來更新顆粒中心位置：. xi(t  t ) =xi(t )  xi(t . 19. t /2). t. (2. 12).

(34) 2.3.4 接觸組成模式(contact constitutive models) PFC3D 中提供三種接觸組成模式，分別為接觸勁度模式 (contact-stiffness model)、滑動模式(slip model)及鍵結模式(bonding model)，茲將各種分述於下。. 2.3.4.1 接觸勁度模式(contact-stiffness model) 接觸勁度提供接觸力與相對位移之間的彈性關係。正向勁度為一割線勁度(secant stiffness)，為一總正向力與總位移之關係。. Fin  K nU n ni. (2. 13). 剪力勁度為一個切線勁度(tangent stiffness)，表示著剪力方向力量增量跟其剪位移之關係。. Fi s  k s U s. (2. 14). PFC3D 有兩種不同的接觸模式，分別為線性接觸模式(Linear contact model)以及非線性接觸模式(Hertz-Mindlincontactmodel)，而不同的接觸模式有不同的勁度計算方式，本研究採用 Linear contact model，其計算方式如下所示。線性接觸模式(Linear contact model)的正向勁度與剪力勁度分別為：  n kn[ A]kn[ B ] K   kn[ A]  kn[ B ] Linear contact model  [ A] [ B ]  K s  ks ks  ks[ A]  ks[ B ] 20.

(35) 2.3.4.2 滑動模式(slip model) 滑動模式為透過摩擦係數 μ 跟正向力 Fin 之乘積，計算出接觸剪力最 s s 大值 Fmax ，以判斷顆粒是否產生滑動， Fmax 表示如下：. s Fmax = Fi n. (2. 15). s s 若 Fin >Fmax 時，即產生顆粒滑動現象，此時 Fin =Fmax ，如下所示：. . s Fis ←Fis Fmax / Fis. . (2. 16). 另外，若接觸之連結為接觸鍵結模式(contact-bond-model)時，則滑動模式不作用。. 2.3.4.3 鍵結模式(bonding model) PFC3D 中允許顆粒之間相互鍵結，共有兩種鍵結模式，分別為接觸鍵結模式(contact-bond model)以及平行鍵結模式(parallel-bond model)。接觸鍵結模式，為假設接觸點上有一組正向及剪力彈簧作用其中，此組彈簧分別提供正向及剪力強度，而此模式之剪力強度用以抵抗滑動行為，故滑動模式將不在此模式下作用。此模式之正向與剪向相對於位移關係圖如圖 13 所示。平行鍵結模式，為假設相鄰接觸之顆粒間，有一圓柱狀之連接塊體 (圖 14)，此塊體可以承受力量與彎矩，以模擬膠結材料之行為(Potyondy and Cundall, 2004)，若所承受之強度達鍵結強度時，則鍵結判定失效。. 21.

(36) (a). (b). 圖 13 接觸鍵結模式行為，(a)正向力-位移關係圖，(b)剪力-位移關係圖 (PFC3D manual). 圖 14 平行鍵結模式示意圖(PFC3D manual). 22.

(37) 平行鍵結由五個參數來定義，分別為正向力方向與剪力方向之勁度與強度，與兩球重疊處之半徑，以上五種參數於 PFC3D 中之指令分別為， pb_kn、pb_ks、pb_nstrength、pb_sstrength、pb_radius。平行鍵結之總力與彎矩分別以 Fi 與 M i 表示，且皆於接觸面上可分為 n. s. n. s. 正向力方向與剪力方向(圖 14)，即 Fi 、 Fi 、 M i 、 M i 。 n. s. Fi  Fi  Fi n i. Mi  M  M n. (2. 17). s i. n. 其中， Fi 與 M i 可表示為： n.    M n n  M n n. Fi  F j n j ni  F ni M. n i. n. j. j. i. (2. 18). i. 當鍵結形成時， Fi 與 M i 之初始值為零。後續之相對位移與旋轉增量轉變為彈性力與彎矩的增量。在 t 時彈性力增量之計算為：. . n. . n. Fi  k AU n ni s. s. Fi  k AU is with U i  Vi t. (2. 19). 彈性彎矩增量之計算為：. . n. . s.  Mi  k J  n ni s. n. Mi  k I is. (2. 20). with i   i[ B ]   i[ A]  t 其中，A 為重疊區域；J 為重疊區域之極慣性矩；I 為重疊區域之一個通過接觸點的軸且在 is 方向上的轉動慣量，可由以下式子得到：. 23.

(38) 2. A R 4 1 J  R 2 4 1 I R 4. (2. 21). 與平行鍵結關聯之新的力和力矩向量由舊的彈性力和力矩增量向量總和計算而成，新的力向量計算為： n. n. n. Fi  F ni  Fi. . s. s. (2. 22). s. Fi  Fi.   Fi. rot .2. 新的力矩向量計算為： n. n. n. Mi  M ni  Mi s.   s. s. Mi  Mi.   Mi. (2. 23). rot .2. 作用於鍵結附近之最大拉伸應力和剪應力的計算為： s. Mi. n.  max   max . F  R A I F. s i. A. . M. n. J. (2. 24). R. 若最大拉伸應力超過了正常的強度(  max   c )，或最大剪應力超過剪切強度( max   c )，則平行鍵結斷裂。 Fi[ A]  Fi[ A]  F i Fi[ B ]  Fi[ B ]  F i. M i[ A]  M i[ A]  eijk  x[jC ]  x[jA]  F k  M i. (2. 25). M i[ B ]  M i[ B ]  eijk  x[jC ]  x[jB ]  F k  M i Fi[  ] 與 M i[  ] 分別為球  的應力與彎矩總和，而 F i 與 M i 可由 2.17 求. 得。. 24.

(39) 2.3.5 Burger’s model 除 2.2.4 節中所提及之三種內建接觸組成模式外，PFC3D 另外提供由使用者選擇的模式，Burger’s model 便是其中之一。此種模式主要常被用於模擬潛變之行為，由 Kelvin model 與 Maxwell model 組成。在此對其模式之計算做簡單敘述。. 圖 15PFC3D 之 Burger’s model 示意圖(PFC3D manual) Burger’smodel 的總位移是由 Kelvin model 部分的位移( uk )與 Maxwell model 部分的位移( umK 、 umC )總和得出： u  uk  umK  umC. (2. 26). (2.26)式之一次微分與二次微分分別為： u  uk  umK  umC u  uk  umK  umC. (2. 27) (2. 28). Kelvin model 部分的接觸力 f 及其一次微分為： f   Kk uk  Ck uk f   Kk uk  Ck uk 25. (2. 29) (2. 30).

(40) 利用 Maxwell model 部分的勁度 K m 與黏度 Cm 為： f   K mumK. f   K mumK. f   K mumK f  CmumC f  CmumC. (2. 31) (2. 32) (2. 33) (2. 34) (2. 35). 利用(2.27)式至(2.35)式，整理其二次微分方程式為： C  1 1  Ck Cm CC f   k  Cm   f  Cmu  k m u  f  Kk Km Kk  Kk Km   Kk. (2. 36). 根據(2.29)式， uk .  K k uk  f Ck. (2. 37). 利用有限差分法的中心差分近似法對時間取導數，並採用 uk 與 f 的平均值： t 1 t ukt 1  ukt 1  K k  uk  uk  f t 1  f t      t Ck  2 2  . (2. 38). 因此， ukt 1 . 1 t t t 1 t  Bu  f  f    k A 2Ck . (2. 39). 其中， K k t 2Ck K t B 1 k 2Ck A 1. (2. 40) (2. 41). 就 Maxwell model 部分，位移與其一次微分為： um  umK  umC um  umK  umC 26. (2. 42) (2. 43).

(41) 將(2.32)式與(2.35)式代入(2.43)： um  . f f  K m Cm. (2. 44). umt 1  umt f t 1  f t f t 1  f t   t K m t 2Cm. (2. 45). 因此， u. t 1 m. t 1 t f t 1  f t t  f  f     umt Km 2Cm. (2. 46). 總位移與 Burger’s model 的一次微分為： u  uk  um u  uk  um. (2. 47) (2. 48). u t 1  u t  ukt 1  ukt  umt 1  umt. (2. 49). 接觸力 f t 1 之計算為： f t 1  . 1  t 1 B  u  u t  1   ukt  C A .  Df t  . (2. 50). 其中， t 1 t   2Ck A K m 2Cm t 1 t D   2Ck A K m 2Cm. C. 接觸力 f t 1 可由已知值 u t 1 、 u t 、 ukt 及 f t 計算得出。. 27. (2. 51).

(42) 2.4 數值模擬研究 Nichol 等人以 UDEC 程式為分析工具，模擬大規模邊坡脆性岩石及延性岩石之傾覆，探討其破壞模式並找出導致岩石變成脆性及延性的原因(Nichol et al., 2002)。 Preh 與 Poisel 以 UDEC 程式作為分析工具，模擬一種岩石邊坡之破壞型式，此破壞型式稱為“kink band slumping”，藉此探討岩石邊坡之破壞模式(Preh and Poisel, 2004)。李宏輝以 PFC2D 程式作為分析工具，對砂岩微觀組成與互層岩體力學行為進行模擬與分析，以確切了解砂岩與互層岩體在單壓試驗條件下之力學行為，探討其微觀機制與微觀性質(李, 2008)。唐昭榮等人及羅佳明等人利用不同時期數值高程模型與地形圖等資訊，搭配 PFC3D 程式針對 1999 年草嶺及 2009 年小林村之山崩事件，進行災前地形分析與災中山崩動態模擬，並藉由運動速度之監測，探討整個事件之運動過程與特性(唐等, 2009；羅等, 2011)。 Scholtes and Donze 利用分離元素法之理論基礎，並結合三維分析程式，探討破碎性岩體之漸進式破壞(Scholtes and Donze, 2012)。蔣志宏以分離元素法 PFC2D 探討為土石流體內顆粒之運動機制，包括顆粒間接觸力、土體位移、土體流動行為模式及粒間微觀參數影響，其分析成果中顯示，分析出來的土石流流動特性與實際土石流流動特性 28.

(43) 相當相近(蔣,2009)。有鑑於分離元素法(discrete element method)進行模擬時，具可將滑動塊體從變形、破裂到分離的運動過程完整呈現之特性，現今國內外已有許多利用此分析方法進行相關邊坡災害之數值模擬研究，故本研究選用以分離元素法為理論基礎之 PFC3D 模擬軟體進行模擬。. 29.

(44) 第三章研究方法為了解板岩邊坡變形行為，本研究擬利用以分離元素法為理論基礎之模擬軟體(PFC3D)模擬出其變形機制。且為符合 PDC3D 之球元素，並考量到板岩葉理發達之特性，設計出特有的室內物理模型試驗，並進行於不同寬度、弱面間距、坡度、葉理角度、乾燥或泡水弱化狀態下之室內物理模型試驗，再利用 PFC3D 進行模擬，進一步驗證其物理試驗結果。在此將針對本研究中物理模型實驗與 PFC3D 模擬分析分別做方法介紹。. 3.1 物理模型試驗本節針對室內物理模型試驗做介紹，其內容包含試體基本物理性質及力學行為、試驗儀器與試驗方法。. 3.1.1 物理試驗基本假設本研究之物理模型設計主要基本假設分述如下： 1. 假設邊坡運動期間，板岩層沿葉理分離，且本研究不外加其他外力。 2. 簡化邊坡坡度，採用 17 度、30 度、45 度等單一角度為主，即暫不考慮較複雜之坡形或坡面粗糙度。 3. 簡化葉理角度，採用順向 60 度、逆向 90 度與逆向 60 度，即暫不考慮其餘複雜葉理或岩層交錯問題。. 30.

(45) 3.1.2 試體基本物理性質及力學性質本研究選定 PFC3D 作為模擬工具，為與其球元素相吻合，採用研磨石(圖 16)，又因板岩具異向性之特性，將研磨石利用樹脂(由聚醋酸乙烯酯、水組成)膠結後，形成一片狀試體(圖 17)。研磨石之基本物理性質及力學性質如下：直徑 5 mm，密度 2.53 g/cm3，摩擦角 48.9˚(圖 19)，滾動摩擦角 6.97˚。片狀試體之摩擦角 48.59˚(圖 21)。樹脂之拉應力 7.9 MPa(圖 22)。. 圖 16 研磨石. 圖 17 片狀試體. 圖 18 研磨石之直接剪力實驗結果 31.

(46) 圖 19 研磨石之破壞包絡線. 圖 20 片狀試體之直接剪力實驗結果. 32.

(47) 圖 21 片狀試體之破壞包絡線. 圖 22 樹脂之拉應力. 33.

(48) 3.1.3 試驗儀器本研究所使用之主要實驗設備包含兩大類：(1)影像分析設備；(2)物理試驗設備。影像分析設備部分，採用 Canon EOS 400D 相機，其有效畫素為 1010 萬畫素；物理試驗設備部分，根據前人研究經驗並考量未來相關模型試驗之研究需求，自主研發可模擬各岩層間崩塌行為之試驗模型，並將其分為三大部分：源頭區、滑動區及堆積區(圖 23)。茲將物理試驗設備分別說明如下：. 圖 23 物理試驗模型 1. 源頭區：因為本實驗考量乾燥、材料泡水弱化兩種狀態，故製作兩種試驗箱，其比較圖如圖 24 所示。乾式試驗箱的長、寬、高分別為 200 mm、200 mm、500 mm，側邊為了方便觀測試體運動行為，故使用透明壓克力板；材料泡水弱化用試驗箱的長、寬、高分別為 200 mm、. 34.

(49) 220 mm、455 mm，與乾式試驗箱不同之處為，為了在試驗箱中加水，使試體產生弱化，前面會加放玻璃板，即為圖 24(b)方框處。 2. 滑動區、堆積區：考量到邊坡運動過程中，可能造成試體滑落並堆積之情形，因此分別設置滑動區與堆積區。. (a). (b). 圖 24 物理試驗箱，(a)乾燥狀態，(b)泡水弱化狀態. 3.1.4 試驗方法為瞭解邊坡變形行為，本試驗方法考量試體於乾燥與泡水弱化兩種情況，並分別選用四種參數(坡度、葉理角度、弱面間距及寬度)，其編號規則如下：D 代表乾燥狀態；W 代表泡水弱化狀態；T 代表弱面間距，1、 2 分別為兩弱面間試體層數；B 代表寬度，3、12 分別為顆粒數；S 代表坡度，17、30、45 分別為其角度；F 代表葉理角度，C60、O90、O60 分. 35.

(50) 別表示順向 60 度、逆向 90 度、逆向 60 度三種角度。舉例而言，於乾式試驗中，弱面間距 1 層、寬度 3 顆、坡度 17 度、葉理順向 60 度之試驗，其編號為 DT1B3-S17FC60。. 圖 25 各參數物理意義示意圖. 茲將試驗步驟分述如下： 1. 試體製備本試驗利用樹脂將研磨石彼此膠結，步驟分述如下。步驟一：將第一排研磨石排列於鐵尺與塑膠尺中間，排約 10 公分左右，接著在研磨石上擠上體積約為 220 mm3 之樹脂。步驟二：將塑膠尺以靠近鐵尺那端為基點，旋轉適當角度後，將第二排之研磨石沿著塑膠尺放下去，使其剛好卡在第一排兩顆研磨石間。步驟三：本試驗中試體的長度為 20 公分，故步驟一、二完成後，沿 36.

(51) 著鐵尺再重複步驟一與二，其完成圖如下圖所示。步驟四：將塑膠尺移開，並從第二排研磨石旁擠上樹脂，其體積約為 440 mm3，接著將第三排研磨石從旁邊往第二排兩顆研磨石中間置入，完成後即為長 20 cm、寬為 3 顆之片狀試體，如下圖所示。長 20 cm、寬為 12 顆之片狀試體的話，則是重複執行步驟四，直到所需寬度即可。步驟五：為清楚辨識各岩層間運動行為，本實驗於試體完成並風乾幾天後，進行染色的步驟。即於通風處利用紅色噴漆均勻噴在試體側邊，因僅噴在試體側邊，故已避免材料力學的改變，並放置約五分鐘使其完全陰乾，完成圖如下圖所示。. 步驟一. 步驟二. 步驟三 (長度 10 公分). 步驟三 (長度 20 公分). 圖 26 試體製備步驟. 37.

(52) 步驟四 (上膠). 步驟四 (排列). 步驟五 (染色前). 步驟五 (染色後). 圖 26 試體製備步驟 2. 選定變因選定執行乾燥或泡水弱化試驗後，選定上述試驗之四種變因，坡度、葉理角度與片狀試體之寬度與弱面間距後，調整試驗箱至適當位置，擺放相機至固定位置，且調整光圈、快門等條件至最佳狀態。 3. 排列試體將原片狀試體與染色後片狀試體依序排列於試驗箱中，以明顯區分各岩層，放置片狀試體直到試驗箱擺滿為止，如圖 27 所示。 4. 開始實驗將隔板放開，試體逐漸產生變化，並使用相機拍攝其運動行為，試驗規劃表如表 2 所示，其中 D 代表乾燥狀態；W 代表泡水弱化狀態；T 代表弱面間距；B 代表寬度；S 代表坡度；F 代表葉理角度。. 38.

(53) 圖 27 片狀試體依序排列情形. 表 2 試驗規劃表弱面間距、環境. 弱面間距、葉理角度. 環境. 寬度、坡度. 葉理角度寬度、坡度. D. T1B3S17. FC60、FO90、FO60. W. T1B3S17. FC60. D. T1B3S30. FC60、FO90、FO60. W. T1B3S30. FC60. D. T1B3S45. FC60、FO90、FO60. W. T1B3S45. FC60、FO90. D. T1B12S17 FC60、FO90、FO60. W. T1B12S17. FC60. D. T1B12S30 FC60、FO90、FO60. W. T1B12S30. FC60、FO90. D. T1B12S45 FC60、FO90、FO60. W. T1B12S45. FC60. D. T2B12S30. D. T2B12S45 FC60、FO90、FO60. FO90. 39.

(54) 3.2 數值模擬分析為探討板岩邊坡變形機制，本研究採用以分離元素法為理論基礎之三維分析軟體 PFC3D(Particle Flow Code in 3Dimensions)，相較於其他分離元素軟體，具有下列使用上的優點：(1)顆粒的位移量限制不受傳統位移諧和條件之約制，可模擬粒狀材料之流動與大變形問題；(2)由顆粒連結所模擬的塊體，當所受的應力超過其鍵結強度時，允許塊體可以發生材料破裂；(3) 可透過程式提供的指令將圓形顆粒組合成不規則之多邊形塊體；(4) 考量到邊坡運動過程中，易有消能、碰撞、破裂之互制行為。. 3.2.1 數值模擬物理模型試驗數值模擬物理模型試驗，係指物理試驗模型配置為基礎，作為數值模型及參數設定之簡化依據。數值模擬物理模型試驗主要目的在模擬物理模型試驗之成果，以檢核此種數值分析方法之可行性。由於 PFC3D 接受多種度量衡單位制度，因此在建置模型前，須事先從基本參數單位表中挑選一套單位，本研究所選擇之單位如下表所示。表 3 基本參數單位表 Length. m. Density. 103 kg/m3. Force. kN. Stress. kPa. Gravity. m/ sec2. Ball stiffness. kN /m. Parallel bondstiffness. kPa / m. 40.

(55) 本研究之物理實驗模型中運動區部分以牆元素進行模擬，試體部分以球元素加上鍵結來進行模擬，其模型建立方法與模擬項目分別說明如下，模型建立後之範例圖如圖 28 所示： 1. 牆元素：主要用來模擬運動區，並於模擬各項實驗時變化 17 度、30 度、45 度等三種地形坡度。 2. 球元素：將一個個直徑 0.5 cm 之球元素透過鍵結模式，連結成一片片長 20 cm、寬為 3 顆或 12 顆的片狀試體，再逐層擺放在試驗箱中，與實際物理模型實驗中板岩層之排列堆置方式相同。 3. 模擬項目：包括改變坡度(17 度、30 度、45 度)、葉理(順向 60 度、逆向 90 度、逆向 60 度)、寬度(3 顆、12 顆)、弱面間距(1 層、2 層)，並設定乾燥或泡水弱化狀態。. 圖 28 數值模擬模型範例圖. 41.

(56) 3.2.2 數值模擬物理模型試驗之參數設定根據前述物理試驗模型中研磨石之各基本物理性質及力學性質，設定成 PFC3D 中顆粒之各參數，其顆粒之參數設定呈現於下表之中。表 4 顆粒參數設定 Radius. 0.0025 m. Density. 2.53×103 kg/m3. Friction coefficient. 1.15. Gravity. 9.81 m/ sec2. Ball stiffness. kn=105kN /m；ks=4.167×104kN /m. Parallel bond stiffness. kn=105 kPa / m；ks=105 kPa / m pb_nstren=7.9 MPa pb_sstren=7.9 MPa. Parallel bond strength. 42.

(57) 第四章物理試驗暨數值模擬驗證成果本章將呈現物理模型試驗成果及數值模擬物理試驗成果，旨在利用數值模擬成果驗證物理試驗成果，藉此確立此分析方法之可行性。. 4.1 物理模型試驗為確立模擬方法之可行性，本研究設計出特有的物理模型試驗，於第三章時已說明此試驗方法與流程，在此將呈現其物理模型試驗成果。. 4.1.1 乾式試驗為了解試體僅受重力情況下之變形行為，於開始試驗後，本試驗便以定時拍攝來觀察物理試驗之變形行為，下圖為 DT1B3-S30FO90 之物理試驗過程，由圖 29 顯示，試體先產生開裂後，再往下折彎，形成ㄑ字型變形。五分鐘後，全部試體皆已變形，且無明顯變化，故視為該試驗結束。在此將以此組試驗做為基準，探討不同變因下之變形行為差異。. (a). (b) 圖 29 DT1B3-S30FO90 試驗過程，(a)試驗前，(b)經過 30 秒. 43.

(58) (c). (d) 圖 30 DT1B3-S30FO90 試驗過程，(c)1 分鐘後，(d)5 分鐘後. 1. 坡度坡度為 17 與 30 度之試體，其實驗約於 1 分鐘後即無明顯變化；而坡度為 45 度之試體，其實驗約於 30 秒後即無明顯變化。得出當坡度越陡，因受重力影響較遽，故其自開始變化至穩定現象之時間間隔較短。 (a). (b). (c). 圖 31 DT1B3-FO90，(a)坡度 17 度，(b)坡度 30 度，(c)坡度 45 度 2. 寬度進行物理試驗時，為防試體從側邊滑動，故於側邊利用保麗龍來提供較高摩擦力，寬度 3 顆情況下之試驗結果底部因受保麗龍的束制，故. 44.

(59) 底層無太大變形，但寬度 12 顆之試體因其自重較重，保麗龍的摩擦力不足以支撐，故底層之試體呈現往前傾倒之趨勢。. (a). (b). 圖 32 DT1-S30FO90 試驗成果，(a)寬度 3 顆，(b)寬度 12 顆. 3. 葉理角度由下圖顯示，當試體葉理角度為順向情況下，其整體並無明顯變形。但當試體葉理角度為逆向情況下，其葉理整體結構會先崩壞，導致片狀試體產生ㄑ字型折彎現象。 (a). (b). (c). 圖 33 DT1B3-S30 試驗成果，(a)葉理順向 60 度，(b)葉理逆向 90 度，(c) 葉理逆向 60 度. 45.

(60) 4. 弱面間距由源頭區上方往下看時，試體斷面圖如圖 34 所示，其中試體直徑為 b，則弱面間距 1 層之慣性矩為12 . 120 .  b4 64.  b4 64. ；弱面間距 2 層之慣性矩為. ，即弱面間距 2 層之慣性矩為弱面間距 1 層之 10 倍，因此弱面. 間距 2 層之試體難以產生變形，故弱面間距 2 層之試體雖處於高傾角葉理角度情況，但其試驗成果呈現無明顯變化的現象。. 圖 34 試體斷面圖. (a). (b). 圖 35 DB12-S30FO90 試驗成果，(a)弱面間距 1 層，(b)弱面間距 2 層. 46.

(61) 4.1.2 材料泡水弱化試驗為了解試體僅受重力情況下之變形行為，於開始試驗後，本試驗便以定時拍攝來觀察物理試驗之變形行為，由下圖顯示，於葉理順向 60 度情況下，試體因泡水讓膠軟化，會於底部先產生開裂(圖 36(b))，之後又因上方試體持續向下擠壓滑動，迫使下方泡在水中的軟化部分因受擠壓產生彎曲褶皺(圖 37(c))。經過 40 小時候，已有大部分試體產生彎曲褶皺，且產生變形之區域因為隨重力往下滑動，因此擠壓在隔板處，無法再有明顯變形，故判斷試體已呈現穩定，視為該試驗結束。在此將以此組試驗做為基準，探討不同變因下之變形行為差異。. (a). (b). 圖 36 WT1B3-S30FC60 試驗過程，(a)試驗前，(b)經過 1 小時. 47.

(62) (c). (d). 圖 37 WT1B3-S30FC60 試驗過程，(c)經過 2 小時，(d)經過 40 小時 1. 寬度於溼式試驗中，試體因泡水讓膠軟化，導致其強度降低，又受重力影響下，上方試體持續向下擠壓滑動，迫使下方泡在水中的軟化部分因受擠壓產生彎曲褶皺，寬度為 12 顆情況下之試體，因自重較重，故試體除於坡趾處產生凸肚變形，後方試體會隨著源頭區底板往下滑動。. (a). (b). 圖 38 WT1-S30FC60 試驗成果，(a)寬度 3 顆，(b)寬度 12 顆. 48.

(63) 2. 坡度隨著坡度 17 度與坡度 45 度情況下之試體，皆因寬度不足，後期試體會往側向翻倒，如圖 39、圖 40 所示。有鑑於此，圖 41 為寬度 12 顆各坡度情況下之試驗成果，由圖中顯示，坡度 30 與 45 度情況下之試驗，其變形程度較坡度 17 度情況明顯，但因試體會不斷往下滑動，後期將擠壓在隔板處，無法再有明顯變形，故無從比較坡度 30 與 45 度情況下試驗之變形程度。. (a). (b). 圖 39 WT1B3-S17FC60 試驗成果，(a)側翻前，(b)側翻後. (a). (b). 圖 40 WT1B3-S45FC60 試驗成果，(a)側翻前，(b)側翻後. 49.

(64) 圖 41 WT1B12-FC60 試驗成果(由左至右分別為坡度 17、30、45 度) 3. 葉理於葉理角度為逆向 90 度之泡水弱化試驗中，一開始試體呈現ㄑ字形的折彎現象，後又因泡水讓膠軟化使其強度降低，且受重力影響下，上方試體持續向下擠壓滑動，迫使下方泡在水中的軟化部分因受擠壓產生彎曲褶皺形，但多次試驗下發現，因前方有玻璃隔板阻擋，造成下方泡在水中之試體難以產生變形，而葉理逆向 60 度之實驗也有相同情形，故於溼式實驗中，僅執行葉理順向 60 度之實驗。. (a). (b). 圖 42 試驗成果，(a)WT1B3-S30FC60，(b)WT1B3-S30FO90. 50.

(65) 4.2 數值模擬物理試驗為了解利用 PFC3D 軟體模擬板岩變形之可行性，本研究根據前述之室內物理模型試驗參數設定數值模擬中各參數，對物理模型試驗成果做比對與驗證，在此將呈現數值模擬物理模型試驗之成果。. 4.2.1 乾式試驗由上一小節得出於乾式試驗中，葉理角度的不同，試體會有不同的變形行為，故在此將不同葉理角度分別作比對驗證。 1. 葉理角度順向 60 度葉理順向 60 度、寬度 3 及 12 顆、弱面間距 1 層與各坡度之數值模擬物理試驗成果如圖 43~圖 48 所示。由圖中顯示，數值模擬成果與物理試驗成果中的試體大部分皆呈現穩定狀態。其中於 DT1B3-S30FC60 中，試體從靠近背板處開始往下滑動，此行為與物理試驗有些微不同，評估此因為在進行物理試驗時，採用保麗龍來提高側壁摩擦力，因此在試驗期間，因保麗龍的摩擦而沒出現滑動現象。由此可判斷，於低葉理情況下，岩層趨於穩定狀態，但有可能產生滑動。. 51.

(66) 圖 43 DT1B3-S17FC60 模擬成果. 圖 44 DT1B3-S30FC60 模擬成果. 圖 45 DT1B3-S45FC60 模擬成果. 52.

(67) 圖 46 DT1B12-S17FC60 模擬成果. 圖 47 DT1B12-S30FC60 模擬成果. 圖 48 DT1B12-S45FC60 模擬成果. 53.

(68) 2. 葉理角度逆向 85 度葉理逆向 85 度、寬度 3 及 12 顆、弱面間距 1 層與各坡度之數值模擬物理試驗成果如圖 49~圖 54 所示。經過多次模擬後，本研究最終採用葉理逆向 85 度之數值模型來模擬葉理逆向 90 度情況下各坡度之物理試驗成果，並將最左邊一層試體與每層試體之最下面幾顆顆粒球使用 PFC3D 中的 fix 功能將其鎖住，以達到最佳模擬成果。圖中顯示，數值模擬坡度 17 度情況下之成果與物理試驗成果皆呈現倒ㄑ字型的變形，且試體折彎程度相當吻合。而數值模擬坡度 30 與 45 度情況下之成果與物理試驗成果雖同樣呈現倒ㄑ字型的變形，但試體的折彎程度不太相同，評估此因為在進行物理試驗時，為防試體從側邊產生滑動，故採用保麗龍來防止此情形發生，因此在進行數值模擬時，模型的設置上考量到此點，故將每層最底部幾顆顆粒球鎖住以更有效的模擬，但在高坡度情況下，保麗龍的摩擦力雖會造成試體變形時的阻力，但不至於阻止底部的研磨石產生移動，因此導致模擬成果與物理試驗成果有些許不同。. 54.

(69) 圖 49 DT1B3-S17FO85 模擬成果. 圖 50 DT1B3-S30FO85 模擬成果. 圖 51 DT1B3-S45FO85 模擬成果. 55.


(71) 3. 葉理角度逆向 60 度葉理逆向 60 度、寬度 3 顆、弱面間距 1 層與各坡度之數值模擬物理試驗成果如圖 55~圖 60 所示。由圖中顯示，數值模擬成果與物理試驗成果之變化皆呈現倒ㄑ字型變形，且試體的折彎程度相當吻合。. 圖 55 DT1B3-S17FO60 模擬成果. 圖 56 DT1B3-S30FO60 模擬成果. 57.


(73) 圖 60 DT1B12-S45FO60 模擬成果. 4.2.2 材料泡水弱化試驗 WT1B3-S30FC60 之數值模擬物理試驗成果如圖 61 所示。由圖中可以發現，模擬成果與物理試驗成果一樣會因泡水弱化後受上部試體往下擠壓，進而於底部產生凸肚行為的變形。但於物理試驗中，因需於源頭區及運動區間設置隔板，以阻擋水的流出，而進行數值模擬時，並未將其前方設置牆來阻擋其變形，故成果與物理試驗成果稍有不同。. 圖 61 WT1B3-S30FC60 之模擬成果. 59.

(74) 4.3 小結根據以上物理試驗成果與數值模擬成果比對成果後，證實採用以分離元素法為理論基礎的數值軟體 PFC3D (Particle Flow Code in 3 Dimensions)進行模擬，此種方法之可行性，在之後的章節將會沿用數值模擬物理試驗時所設定之各參數，進行不同坡度、不同葉理情況下之模擬。除此之外，另以數值模擬物理試驗時所設定之各參數為底，改變不同參數進行模擬，探討其導致出的影響。. 60.

(75) 第五章參數之影響有鑑於第四章數值模擬物理試驗之成果驗證成功，本研究將在此章中沿用其模擬物理試驗之模型與參數，並分別更改摩擦係數、球勁度、平行鍵結勁度與平行鍵結之鍵結強度，以探討不同參數情況下之變形行為影響。. 5.1 摩擦係數順向坡模型改變摩擦係數後之成果如表 5 所示，由成果圖顯示，當摩擦係數介於 0.4 至 0.86 時，在同一模擬步數情況下，隨著摩擦係數的降低，其滑動程度越大，當摩擦係數小於 0.4，其滑動程度已達穩定，無太大變化。而逆向坡模型改變摩擦係數後之成果如表 6 所示，由成果圖顯示逆向坡模型的變形行為與模擬物理試驗時一樣，產生試體往下折彎的現象。. 61.

(76) 表 5 順向坡模型改變摩擦係數 Friction coefficient=0.3. Friction coefficient=0.4. Friction coefficient=0.86. Friction coefficient=1.15. 62. Friction coefficient=0.58.

(77) 表 6 逆向坡模型改變摩擦係數 Friction coefficient=0.3. Friction coefficient=0.4. Friction coefficient=0.86. Friction coefficient=1.15. Friction coefficient=0.58. 為了探討不同參數下之影響，本研究先定義出整個試體的坡高 H 與試體長度 L，再將變形後之試體向下折彎的第一群變形曲率最大處連線，藉由此線之中點定義試體的變形高度 h、變形長度 l 以及變形角度 β，其示意圖如圖 62 所示，並利用坡高與變形高度之間的比值(h / H)、試體長度與變形長度之比值(l/L)與變形角度 β 來作為探討之依據。. 63.

(78) 圖 62 不同參數下影響行為之探討因子. 圖 63 逆向坡模型改變摩擦係數之 h/H、l/L 比值圖. 64.

(79) 圖 64 逆向坡模型改變摩擦係數之變形角度 β 比較圖其成果圖如圖 63 至圖 64 所示，由圖中顯示，試體長度與變形長度之比值(l/L)約位於 0.4 至 0.8 之間，由其趨勢線顯示，當摩擦係數約為 0.6 時，其比值最高。坡高與變形高度之比值(h / H)約位於 0 至 0.2 之間，由其趨勢線顯示，隨著摩擦係數的提高，其 h/H 比值愈高。變形角度 β 之變形角度約介於 30 度到 40 度左右。. 65.

(80) 5.2 球勁度順向坡模型改變球的勁度後之成果如表 7 所示，由成果圖顯示，在同一模擬步數情況下，當球的勁度為 5 104 kN/m 以下時，在模型的左下角部分會產生折彎，且折彎程度隨著勁度越小，折彎情況越厲害。而當勁度小於等於 5 103 kN/m，其折彎程度趨於相近。. 表 7 順向坡模型改變球的勁度 Kn=1103 kN/m. Kn=5 103 kN/m. Kn=5 104 kN/m. Kn=1105 kN/m. 66. Kn=2.5 104 kN/m.

(81) 逆向坡模型改變球的勁度後之成果如表 8 所示，由成果圖顯示在同一模擬步數情況下，逆向坡模型的變形行為與模擬物理試驗時一樣，產生試體往下折彎的現象。. 表 8 逆向坡模型改變球的勁度 Kn=1103 kN/m. Kn=5 103 kN/m. Kn=5 104 kN/m. Kn=1105 kN/m. 67. Kn=2.5 104 kN/m.

(82) 而同樣為了能探討不同參數下之影響，在此一樣根據 5.1 節所提及之探討依據做整理，並將成果呈現如圖 65 至圖 66 所示。試體長度與變形長度之比值(l/L)約位於 0.4 至 0.8 之間，由其趨勢線顯示，比值最大值為 0.78，且隨球勁度的提高，其 l/L 比值愈小。坡高與變形高度之比值(h / H) 約位於 0 至 0.5 之間，由其趨勢線顯示，當球勁度約為 7 104 kN/m 時，其比值最小。此外，根據圖 67 中變形角度之趨勢線，發現於球勁度約為 7 104 kN/m 時，其角度最大。. 圖 65 逆向坡模型改變球勁度之 h/H、l/L 比值圖. 68.

(83) 圖 66 逆向坡模型改變球勁度之變形角度 β 比較圖. 69.

(84) 5.3 平行鍵結勁度順向坡模型改變平行鍵結勁度後之成果如表 9 所示，由成果圖顯示，在同一模擬步數、不同平行鍵結勁度情況下之邊坡皆呈現穩定狀態，由此判斷平行鍵結勁度對順向坡模型之變形無太大影響力。. 表 9 順向坡模型改變平行鍵結勁度 pb_kn=1103 kPa/m. pb_kn=5 103 kPa/m. pb_kn=5 104 kPa/m. pb_kn=1105 kPa/m. 70. pb_kn=1104 kPa/m.

(85) 逆向坡模型改變平行鍵結勁度後之成果如表 10 所示，由成果圖顯示，在同一模擬步數情況下，當平行鍵結勁度等於或小於 1104 kPa/m 時，試體產生波浪狀的折彎現象，且隨著勁度的降低，其波浪狀的折彎越多。. 表 10 逆向坡模型改變平行鍵結勁度 pb_kn=1103 kPa/m. pb_kn=5 103 kPa/m. pb_kn=5 104 kPa/m. pb_kn=1105 kPa/m. 71. pb_kn=1104 kPa/m.

(86) 5.4 平行鍵結鍵結強度順向坡模型改變平行鍵結鍵結強度後之成果如表 11 所示，由成果圖顯示，在同一模擬步數情況下，當鍵結強度小於或等於 79 kPa 時，邊坡會從坡趾處產生碎裂，當鍵結強度大於 170 kPa 時，試體呈現穩定狀態。. 表 11 順向坡模型改變平行鍵結鍵結強度 pb_nstrengh=7.9 kPa. pb_nstrengh=79 kPa. pb_nstrengh=340 kPa. pb_nstrengh=790 kPa. 72. pb_nstrengh=170 kPa.

(87) 逆向坡模型改變平行鍵結鍵結強度後之成果如表 12 所示，由成果圖顯示，試體皆會由上往下產生折彎，但當鍵結強度小於或等於 170kPa 時，接近坡趾處之試體由折彎最顯處產生破碎。. 表 12 逆向坡模型改變平行鍵結鍵結強度 pb_nstrengh=7.9 kPa. pb_nstrengh=79 kPa. pb_nstrengh=340 kPa. pb_nstrengh=790 kPa. 73. pb_nstrengh=170 kPa.

(88) 第六章坡度及葉理角度之影響本章旨在了解試體於不同坡度、葉理角度情形下之影響為何，故在此將沿用模擬物理試驗時之參數，並統一模型尺寸與模擬步數等，但考量到模擬物理試驗成果時，針對葉理角度為順向 60~80 度、逆向 90 度與逆向 60~80 度這三種情況有些許的不同，其相異之處將於後續小節中分別敘述。. 6.1 模型建置為了解試體於不同坡度、葉理角度情形下之影響為何，本章將模型統一尺寸，設定 15 層試體，沿用模擬物理試驗時之參數，以及利用 PFC3D 中的 fix 指令將最右邊一層試體的底部五顆顆粒球鎖住，並統一模擬步數，以減少變因，得出於坡度與葉理角度之影響。但考量到模擬物理試驗成果時，針對葉理順向 60~80 度、葉理逆向 90 度與葉理逆向 60~80 度這三種情況有些許的不同，其相異之處分述如下。葉理順向 60~80 度情況下，會把最左邊一層試體利用 PFC3D 中的 fix 指令將其鎖住；逆向坡 90 度情況下，一樣將最左邊一層試體給予鎖住，且模擬時的葉理角度調為逆向 85 度；而葉理角度低於逆向 90 度時，將不把最左邊一層試體鎖住。整合上述幾點後所建置完成的模型，設定其僅受重力並跑相同步數後，整理之成果如表 13 所示，由表中顯示，乾式情況下且葉理角度位於 74.

(89) 順向 60~80 度時，試體呈現穩定狀態。根據模擬物理試驗時的經驗，當試體受到水的弱化後，會於坡趾處產生折彎變形，故本研究將葉理角度位於順向 60 度至 80 度之模型取試體坡高之 1/3 高度進行弱化，使其產生變形，並將其成果整理如表 14 所示。表 13 乾式情況模擬成果坡度葉理角度. 17. 30. 45. 順向 60 度. 順向 70 度. 順向 80 度. 逆向 85 度. 逆向 80 度. 75. 60.

(90) 坡度葉理角度. 17. 30. 45. 60. 逆向 70 度. 逆向 60 度. 表 14 坡度葉理角度. 17. 泡水弱化情況模擬成果 30. 45. 順向 60 度. 順向 70 度. 順向 80 度. 76. 60.

(91) 其中又為了能探討同坡度、不同葉理情形下之影響，以及同葉理、不同坡度情形下之趨勢為何，本研究定義整個試體的坡高 H、試體厚度 T、試體長度 L，再將變形後之試體向下折彎變形曲率最大處連線，藉由此線定義出變形厚度 t、變形高度 h、變形長度 l 以及變形角度 α 或 β，其示意圖如圖 67 所示，利用坡高與變形高度之間的比值(h / H)、試體深度與變形深度之比值(t/T)、試體長度與變形長度之比值(l/L)與變形角度(α 或 β)之趨勢線趨勢來作為探討的依據。. (a) 圖 67 探討因子示意圖，(a)泡水弱化情形，(b)乾式情形. 77.