Original Article

Tunnel and Underground Space. August 2021. 270-288
https://doi.org/10.7474/TUS.2021.31.4.270

ABSTRACT


MAIN

  • 1. 서 론

  • 2. DECOVALEX-TASK G

  • 3. 입자기반 개별요소모델을 통한 암석균열 모델링

  • 4. 벤치마크 모델링 결과 및 토의

  •   4.1 역학 조건 벤치마크 모델링

  •   4.2 수리-역학 조건 벤치마크 모델링

  • 5. 요약 및 결론

1. 서 론

암반은 불균질 재료로서 열적, 수리적, 역학적 거동은 암석의 구조적, 광물학적 특성뿐만 아니라 불연속면 조건에 의해 크게 좌우된다. 여러 외부요인에 의한 미세균열의 개시와 성장, 전파 특성은 암석 자체의 비선형 거동 및 취성파괴와 밀접한 관계가 있으며, 절리나 단층 등 거시적 규모의 불연속면은 암반 거동의 불안정성과 불확실성을 초래하는 결정적인 인자로 알려져 있다. 국내외적으로 고온, 고압, 고응력 조건의 심부 암반을 대상으로 하는 에너지 저장 관련기술, CO2 지중저장, 방사성폐기물 지층처분과 관련된 기술 수요가 증가하면서 열-수리-역학적으로 연계된 불연속 암반체의 복합거동이 중요한 연구주제로 인식되고 있으며, 많은 연구들에서 수치해석적 기법을 통해 이를 이해하고 평가하려는 시도가 계속되고 있다.

암반공학 분야에 사용되는 여러 수치해석기법 중 불연속면의 거동을 직접적인 방식으로 재현하는 방식을 불연속체 접근법(discontinuum approach)이라 통칭하며, 개별요소법(distinct element method, DEM)이 그 대표적인 예이다. 개별요소법에서는 암반을 개별적으로 거동하는 블록의 집합체(blocky system)로 가정하고, 블록과 그 경계에 구성식과 입력물성을 부여하여 그 상호작용을 해석한다. Itasca 사의 UDEC(Itasca, 2014)과 3DEC(Itasca, 2017)이 개별요소법을 대표하는 상용소프트웨어이며, 주로 암반구조물이나 암반사면 해석에서 불연속면의 영향을 평가하는 도구로 사용되어 왔다(Cundall, 1988, Gu and Huang, 2016, Li et al., 2019, Hu et al., 2020). 최근에는 이를 통해 암석의 구조적 특성과 미세균열의 전파특성을 미시적 관점에서 모델링하려는 연구들도 다수 보고되고 있다(Lan et al., 2010, Ghazvinian et al., 2014, Park et al., 2017, Fu et al., 2020, Wang et al., 2021). 한편, 유한요소법(finite element method, FEM)과 개별요소법을 병용하는 finite-discrete element method(FDEM)를 통해 암석의 불균질성과 구조적 특성, 파괴 메커니즘 등을 모사하려는 연구도 다수 발표되고 있다(Mahabadi et al., 2012, Lisjak et al., 2014, Wu et al., 2018). 이 방법은 연속체 기반의 유한요소를 통해 초기모델을 구성하되, 파괴기준에 도달했을 때 요소간 파괴와 분리를 허용함으로써 균열 개시와 전파, 새로운 개별 블록의 발생 등을 고려하는 방법이다.

본 연구에서는 입자기반 개별요소모델(grain-based distinct element model, GBDEM)을 통해 결정질 암석 균열의 역학적, 수리-역학적 거동을 모사할 수 있는 해석기법을 제시하였다. 이는 국제공동연구인 DECOVALEX-2023 Task G의 일환으로 수행된 연구로서, 본 논문에서는 1단계 연구주제로 수행중인 역학적 벤치마크 모델링과 수리-역학적 벤치마크 모델링의 해석결과를 소개하였다. 벤치마크 해석에서는 단일 균열을 포함한 암석 유사 재료를 사면체 개별블록의 집합을 통해 모델링하고, 3DEC을 이용하여 그 거동을 분석하였다. 경계응력과 압력, 균열의 거칠기와 경사에 따른 수치해석을 실시하여 균열의 수직 및 전단방향 거동을 살펴보았으며, 이를 해석해(analytical solution)와 비교함으로써 해석기법을 검증하고자 하였다.

2. DECOVALEX-TASK G

DECOVALEX 프로젝트는 방사성폐기물의 지층처분과 관련하여 수치해석 코드를 개발하고 검증하기 위한 국제공동연구로 1992년에 스웨덴의 SKI(Swedish Nuclear Power Inspectorate)와 미국의 LBNL(Lawrence Berkeley National Laboratory)을 중심으로 시작되어, 현재까지 지층처분을 고려하는 대부분의 국가가 참여하고 있다(Birkholzer et al., 2019). DECOVALEX 프로젝트에서는 4년을 주기로 처분 관련 THMC 문제에 대한 새로운 연구 주제(Task)들을 선정하여 각 주제에 대한 공동연구가 진행되어 왔으며, 현재 2020년 4월에 착수된 DECOVALEX-2023(2020-2023)에서 7가지 연구 주제에 대해 공동연구를 수행하고 있다(Birkholzer and Bond, 2020, Park et al., 2020).

Task G는 독일 Freiberg University of Mining and Technology (TUBAF)/Helmholtz Centre for Environmental Research (UFZ), 독일 DyanFrax 사, 스코틀랜드 Edinburgh 대학교가 공동제안한 연구 주제로, 결정질 암반 내 균열의 생성과 성장 메커니즘 및 균열에서 발생하는 열-수리-역학적 복합거동을 해석하기 위한 수치해석기법을 개발하는 데에 최종목표가 있다. 공동연구에 참여하는 연구팀들은 모두 벤치마크 모델링(benchmark modeling)과 블라인드 모델링(blind prediction)을 통해 각 팀의 해석모델을 개발, 검증하게 된다. 이후, 해석모델을 실험실 실험과 현장조건 모델링에 적용하여 연구팀간의 결과 비교와 토의를 거쳐 해석모델을 단계적으로 개선하게 된다.

Park et al.(2020)은 Task G의 상세한 소개와 함께 역학 조건 벤치마크 모델링의 해석결과를 보고한 바 있다. 주요 연구내용은 2차원 암석 시험편 내부에 포함된 단일 균열의 수직 및 전단거동을 모델링하고, Pollard and Segall(1987)의 해석해와 비교함으로써 해석기법을 검증하는 것이다. 그들은 3DEC을 이용하여 강체 블록(rigid block)으로 구성된 해석모델의 역학적 거동을 모델링하였으며, 균열에 유도된 응력을 바탕으로 구한 전단 및 수직변위가 해석해와 거의 동일하게 나타남을 확인하였다. 그러나 균열에 작용하는 응력이 이론해에 비해 낮게 계산되었고, 도메인 크기에 따른 경계효과(boundary effect)가 있는 것으로 나타나 추가연구의 필요성을 시사하였다. 위 연구에서는 강체 블록을 가정하였으므로 해석도메인의 크기에 비해 블록의 크기가 충분히 작지 않으면 응력 전달을 정교하게 모사하기 어렵고 포아송비를 구현하기 어렵다는 단점이 있었다. 본 연구에서는 후속연구로서 강체 대신 변형가능한 블록(deformable block)을 사용하여 해석모델을 개선하고, 역학조건 벤치마크 모델링의 해석결과를 업데이트하였다. 추가적으로 균열 거칠기와 압력 조건에 따른 역학적, 수리-역학적 거동을 해석하고 결과를 논의하였다.

3. 입자기반 개별요소모델을 통한 암석균열 모델링

3DEC에서는 해석대상을 개별 블록(block)의 집합체로 모델링하게 되며, 블록과 블록 간 접촉(contact)에서의 거동을 수치적으로 계산하게 된다. 외재적 시간전진 알고리즘(explicit time marching algorithm)에 따라, 접촉점에서는 힘-변위 법칙과 변위로부터 접촉력(contact force)을 계산하고, 블록에서는 Newton의 법칙과 접촉력의 합력으로부터 개별 블록의 운동을 결정하는 과정이 순차적으로 반복된다. 강체 블록(rigid block)이나 변형가능 블록(deformable block)을 가정할 수 있으며, 변형가능 블록의 경우 작은 사면체 요소(zone)로 분할되어 유한차분법(finite different method, FDM)과 거의 동일한 방식에 의해 사면체 요소의 응력과 변형이 해석된다. 해석대상의 특성과 목적에 따라 블록과 접촉점의 구성모델(탄성, 소성)과 입력변수를 설정할 수 있어, 연속체 접근법들에 비해 불연속면에서 발생하는 전단파괴나 인장파괴를 보다 직접적으로 모사할 수 있다는 장점이 있다.

입자기반 개별요소모델은 상기한 접근법을 암석 내 존재하는 미세한 균열에 적용하는 모델이라고 할 수 있다. 즉, 암반이 아닌 암석 자체를 광물과 같이 작은 입자들의 결합체로 표현한 뒤, 점진적인 접촉력의 전달과 결합 손실(미세 균열)을 모사하는 기법이다. 이러한 접근법을 사용하는 대표적인 모델로 원형 또는 구형의 강체 입자를 사용하는 입자결합모델(bonded-particle model, Potyondy and Cundall, 2004)이 있으며, 이는 Itasca 사의 PFC(Itasca, 2021)로 상용화되어 여러 지반, 암반공학 문제에 널리 사용되고 있다(Potyondy and Hazzard, 2008, Park and Song, 2009, Asadi et al., 2012, Yoon et al., 2015, Castro-Filgueira et al., 2020, Ma et al., 2021, Li and Rao, 2021). 그러나 위 모델은 암석의 낮은 인장강도와 내부마찰각을 구현하는 데 한계가 있는 것으로 알려져 있다. 즉, 하나의 미시변수 조합을 통해 암석의 일반적인 인장강도/압축강도 비(보통 1/10 ~ 1/20)를 재현하기 어렵고, 삼축압축 조건에서 구속압을 증가시키더라도 압축강도가 크게 증가하지 않는다. 이는 접촉력이 수직결합력을 초과하여 국부적 파괴가 발생하면 입자가 회전모멘트에 저항하지 못하고 거시적 파괴를 촉발하게 되므로, 접촉의 전단결합력을 아무리 증가시키더라도 재료의 압축강도 증가를 재현할 수 없기 때문이다(Potyondy, 2015). 결국 이러한 한계점은 입자간의 물리적 결합을 원형 또는 구형 입자의 접촉을 통해 모사하려는 데에서 기인한다. 이를 극복하기 위하여 Clump logic(Cho et al., 2007), Smooth joint model(Potyondy, 2010), Flat joint model(Potyondy, 2012, Potyondy, 2013) 등 입자간 맞물림(interlocking)과 접촉의 부분파쇄 등을 고려할 수 있는 모델링 기법들이 제안되어 왔으나 이들은 대부분 상기한 단점을 보완하는 순기능보다는 모델작성과 미시변수결정 과정에 더 많은 변수와 복잡성을 추가한다는 한계가 있다(Park et al., 2017).

Fig. 1에서와 같이, 압축하중을 받고 있는 암석의 인장파괴 메커니즘을 고려할 때 다각형이나 다면체 입자를 사용하는 개별요소모델은 원형 입자를 사용하는 경우에 대한 좋은 대안이 될 수 있다. 최근 십년 사이 여러 연구에서 다각형, 다면체, 랜덤 Voronoi 모델을 활용한 개별요소모델을 암석의 역학적 거동을 재현하고자 하는 시도가 활발히 이루어져 왔다(Lan et al., 2010, Ghazvinian et al., 2014, Park et al., 2017, Fu et al., 2020, Wang et al., 2021). 위 연구들에서는 광물 스케일의 다각형 또는 다면체 입자를 사용하는 개별요소모델이 입자의 맞물림(interlocking), 접촉력과 모멘트 저항, 균열의 개시와 성장과정을 직접적으로 모방함으로써 실제 암석에서 관찰되는 파괴 메커니즘, 이방성, 불균질성, 취성거동(파괴후거동, 인장강도/압축강도 비) 등을 효과적으로 재현할 수 있음을 보고하였다.

/media/sites/ksrm/2021-031-04/N0120310404/images/ksrm_31_04_04_F1.jpg
Fig. 1

Grain-based distinct element model: (a) micro-mechanisms for compression-induced tensile stress at the contacts between the particles (polygonal grain vs. circular grain, after Ghazvinkin et al., 2014) and (b) cylindrical model in 3DEC

4. 벤치마크 모델링 결과 및 토의

4.1 역학 조건 벤치마크 모델링

4.1.1 해석개요

역학 조건 벤치마크 모델링의 주요 연구내용은 암석 경계에 작용하는 압축응력에 의해 내포된 균열에 발생하는 역학적 거동을 모델링하는 것이다. 주어진 해석영역의 크기는 0.5 m×0.5 m의 2차원 평면이며, 암석 내부에 0.17 m 길이의 단일 균열이 존재하는 것으로 가정한다. 균열의 경우, 거칠기에 따라 편평한 면(planar fracture)과 거친 면(rough fracture)이 모두 고려되며, 편평한 면은 수평면으로부터 0°, 30°, 45°, 60°, 90°의 5가지 경사각을, 거친 면의 경우 60° 경사각을 갖는다. 거친 면의 모델링을 위해서는 실제 암석 절리면에서 측정된 세 가지 거칠기 프로파일이 사용된다. 모델링 시 암석은 탄성재료로서 암석 파괴와 균열의 성장은 고려하지 않으며, 균열의 강도(전단강도, 인장강도)는 0인 것으로 간주한다. Table 1은 역학 조건 벤치마크 모델링의 입력물성 및 해석조건을 요약한 것이다.

Table 1.

Mechanics benchmark model condition

Modeling condition and parameter Value
Domain size (m2) 0.5 × 0.5
Fracture length (m) 0.17
Boundary x-stress, σx (MPa) 5 (compression)
Boundary y-stress, σy (MPa) 10 (compression)
Initial stress (MPa) 0
Rock Young’s modulus (GPa) 50
Poisson’s ratio 0.26
Fracture Cohesion (MPa) 0
Friction angle (°) 0
Tensile strength (MPa) 0
Planar fracture inclination angle, α (°) 0, 30, 45, 60, and 90
Rough fracture inclination angle, α (°) 60
Roughness profile 3 different profiles

본 연구에서는 3차원 모델인 3DEC을 사용하였으므로, 2차원 문제를 모사하기 위해 평면외(out-of-plane) 방향으로 평면변형률 조건을 설정하였다. 즉, 0.5 m×0.5 m × 0.005 m 크기의 얇은 육면체 형태의 해석도메인을 작성하되, 해석 시에서는 두께 방향을 변형을 구속하여 2차원 해석조건을 구현하였다. 3DEC에서는 제공하는 균열 생성 기능을 사용하면, 해석도메인을 가로지르는 무한한 균열(infinite fracture)만을 모델링할 수 있으나, 개별 블록의 집합체를 사용하면 인접한 블록 간 경계면(contact 요소)에 균열 물성을 부여함으로써 유한한 균열(finite fracture)을 모사하는 것이 가능하다. 그리드 작성을 위하여 오픈소스 소프트웨어인 Gmsh(Geuzaine and Remacle, 2009)를 사용하였으며, 균열 프로파일의 좌표에 상응하는 위치에서 입자 간 경계가 연속적으로 발생할 수 있도록 함으로써 단일 균열을 포함한 사면체 해석그리드를 생성하였다. 3DEC에서는 내장함수(FISH)를 이용해 Gmsh 그리드 정보를 개별 블록(edge, vertex)의 정보로 변환하여 개별 블록을 생성하였으며, 블록 생성이 완료된 후 프로파일 좌표에 해당하는 접촉점(contact)에 균열 물성을 부여하였다. Fig. 2는 위 과정을 통해 3DEC에서 작성한 해석모델을 보여주는 것으로 입자(블록)의 개수는 총 75,000개 내외이며, 평균 부피를 이용해 환산한 등가 직경(equivalent spherical diameter, de)은 약 3.2 mm이다.

/media/sites/ksrm/2021-031-04/N0120310404/images/ksrm_31_04_04_F2.jpg
Fig. 2

Model generation for mechanics benchmark simulation

벤치마크 모델링의 검증을 위한 Pollard and Segall(1987)의 해석해(analytical solution)는 탄성재료 내부에 존재하는 강도가 없는 균열을 대상으로 한다. 이에 따르면 2차원 평면변형률 조건에서 응력변화에 의해 유도되는 균열 수직변위와 전단변위는 다음과 같이 표현된다. 이때 변위는 균열을 이루는 두 표면(fracture surfaces)의 상대변위를 의미하며, 수직변위의 경우 인장응력이 작용하는 경우에만 유효하다.

(1)
un=σn4(1-nu2)Ea2-r2
(2)
us=τ4(1-ν2)Ea2-r2

여기서 △un과 △us는 균열 수직변위와 전단변위, △σn, △τ는 균열 수직응력과 전단응력의 변화, a는 균열길이의 절반, r은 균열 중심으로부터의 거리를 의미한다.

4.1.2 등가연속체 접근법에 의한 미시물성 산정

일반적으로 입자기반 개별요소모델의 입력물성을 선택하기 위해서는 연속체 모델과는 달리 미시물성 조정 과정(micro-parameter calibration)이 필요하다. 위 모델은 전체 재료의 물성이 아닌 미시물성(입자나 접촉의 물성)을 대상으로 하며, 일반적인 암석역학, 암반공학 문제에서 이들은 미지값이거나 측정이 매우 어렵기 때문이다. 보통 해석모델이 기지의 실험실 시험(압축, 인장, 전단시험 등) 결과를 재현할 때까지 미시물성을 반복적으로 개선, 적용하여 최적의 조합을 결정한 뒤, 이를 입력변수를 하여 본 모델링(main modeling)을 수행한다.

본 연구에서는 해석대상의 탄성거동을 구현하기 위하여 위 시행착오법(trial and error)을 거치지 않고, 등가연속체 접근법(equivalent continuum approach)을 적용해 있는 미시물성을 결정할 수 있는 방법을 제시하였다. 등가연속체 접근법이란 불연속 암반에서 절리 밀도가 높고, 관심영역의 규모에 비해 절리간격이 상당히 작을 때, 암반을 등가의 연속체로 가정하여 변형 특성을 기술하는 방식이다. 이에 따르면 암반의 등가탄성계수(Em)와 등가전단탄성계수(Gm)는 무결암과 절리 조건을 반영해 다음과 같이 수식화할 수 있다.

(3)
1Em=1Er+1Kns
(4)
1Gm=1Gr+1Kss

여기서 ErGr은 각각 무결암의 탄성계수와 전단탄성계수, KnKs는 절리의 수직강성과 전단강성, s는 절리의 간격이다.

입자기반모델에 위 개념을 적용하면 식 (3)식 (4)의 좌항은 모델 전체의 변형, 우항은 블록 및 접촉의 변형에 상응한다. 입자기반모델의 탄성계수와 포아송비를 각각 Em, νm, 개별블록의 탄성계수와 포아송비를 각각 Eb, νb라 하고, Eb=nEm, νb=νm의 관계를 가정하면, 접촉의 수직강성과 전단강성을 식 (3)식 (4)를 변환해 계산가능하며, 식 (5)식 (6)으로 수식화할 수 있다.

(5)
Kn=ns(n-1)E
(6)
Ks=n2s(n-1)(1+ν)E

여기서 n은 모사하고자 하는 전체모델의 탄성계수(Em)에 대한 개별블록의 탄성계수(Eb)의 비로 항상 1보다 크며, s는 절리간격으로서 입자의 등가직경 de와 동일한 것으로 가정하였다.

본 연구에서는 전체시료의 변형에 대한 블록의 변형과 접촉의 변형이 기여도가 각각 99:1, 50:50, 1:99인 세 가지 경우(Case 1, Case 2, Case 3)를 가정하여 상기 접근법에 의해 탄성거동과 관련된 미시물성을 산정하였다. 한편, 각 해석 케이스의 모델에 대해 일축압축시험(최대 압축응력 50 MPa)을 수치적으로 수행하여 각 해석케이스가 목표한 탄성정수(Em=50 GPa, νm=0.26)를 재현할 수 있는지 확인하였다. Table 2는 입력물성 및 해석결과를 정리한 표로 Eν는 일축압축시험에서 얻어진 응력-변형률 곡선으로부터 구한 입자기반모델의 탄성정수이다. Fig. 3은 Case 1의 일축압축시험 해석결과를 도시한 것으로 응력-변형률 관계 및 압축시험 종료 시 발생한 변위벡터를 보여준다.

Table 2의 해석결과를 살펴보면, 블록과 접촉에 적용한 입력물성의 상당한 차이에도 불구하고 세 해석케이스에서 유사한 탄성거동(탄성정수 E, ν)을 보이는 것을 확인할 수 있다. 이는 등가연속체 접근법을 적용한 미시물성 산정방법이 합리적인 결과를 도출함을 의미한다. 세 해석 케이스 중 블록의 탄성계수가 상대적으로 낮고, 접촉의 강성이 크다고 가정한 Case 1의 결과가 목표값(Em, νm)과 가장 작은 오차를 보였고, 접촉의 기여도를 크게 가정할수록 오차가 커지는 결과가 나타났다. 이러한 경향성은 식 (5)식 (6)에서 s를 등가직경과 동일하게 가정한 데에서 기인하는 것으로 추정되며, 향후 개별블록의 형태에 따른 s 산정 방식에 대한 추가연구가 필요할 것으로 판단된다.

Table 2.

Determination of micro-parameters of GBDEM by equivalent continuum approach

Parameter Case 1 Case 2 Case 3
Input parameter 1Eb:1Kns 99 : 1 50 : 50 1 : 99
n=EbEm 10099 2 100
Block Young’s modulus, Eb (GPa) 50.5 100 5000
Block Poisson’s ratio, νb (assumed to be νm) 0.26 0.26 0.26
Contact normal stiffness, Kn (GPa/m) 1.56×106 3.13×104 1.58×104
Contact shear stiffness, Ks (GPa/m) 6.20×105 1.24×104 6.26×103
Joint spacing, s(m) (assumed to be de) 0.0032 0.0032 0.0032
Simulation result GBDEM model Young’s modulus, E (GPa) 49.6 48.7 48.4
GBDEM model Poisson’s ratio, ν 0.26 0.25 0.25

/media/sites/ksrm/2021-031-04/N0120310404/images/ksrm_31_04_04_F3.jpg
Fig. 3

Results of uniaxial compression test on GBDEM (Case 1 in Table 2): stress-strain curve (unit of stress: Pa) and displacement vector (unit of displacement: m)

벤치마크 모델링의 본 해석(main modeling) 시 Case 1의 블록 탄성계수, 블록 포아송비, 접촉 수직강성, 접촉 전단강성을 적용해 암석 시험편을 생성하였으며, 접촉에서의 인장강도와 전단강도는 매우 크게 산정하여 탄성재료로서의 재료 특성을 모사하였다. 균열(편평한 면, 거친 면) 프로파일 좌표에 해당하는 접촉에는 인장강도와 전단강도를 모두 0으로 설정하였으며, 전단강성과 수직강성은 별도의 조정 없이 암석 시험편과 동일한 것으로 가정하였다. 해석모델의 초기조건은 응력이 작용하지 않는 상태이며, XY 평면의 도메인 경계에는 일정 응력조건(σx = 5 MPa, σy = 10 MPa)을 설정하고, 평형상태에 이를 때까지 역학적 거동해석을 수행하였다.

Park et al.(2020)의 연구에서는 Pollard and Segall(1987)의 해석해를 통해 환산된 수직강성과 전단강성을 균열 접촉에 부여하고, 잔류인장강도와 잔류전단강도를 임의로 설정해 파괴 상태에 이른 이후에도 강성을 통해 변위가 제어되도록 하였다. 이러한 방법을 사용하면 균열의 응력에 비례하는 선형적 변위를 유도하여 해석해와 가까운 결과를 도출할 수 있지만, 해석해가 없는 임의문제에 적용하기에는 어렵다는 단점이 있다. 반면, 본 연구에서는 균열의 강도가 0이므로 재하와 동시에 파괴상태로 간주되며, 균열 강성이 아닌 암석의 탄성거동에 의해 균열 변위가 발생하게 된다.

4.1.3 해석결과

가) 편평한 균열

Fig. 4Fig. 5는 편평한 균열의 경사각(α)에 따른 해석결과로 전자는 암석의 변위크기 분포와 균열의 최대전단변위를 보여주며, 후자는 암석의 전단응력(σxy) 분포와 균열 끝단에서의 최대, 최소 전단응력를 보여준다.

전반적인 변위 분포를 살펴보면, 모든 해석케이스에서 도메인 중심에 대해 대칭적인 결과를 보였으며, 압축력이 작용하는 경계조건으로 인해 경계에서 해석 도메인의 중심으로 갈수록 작은 변위가 계산되었다. 경사각이 0°와 90°인 경우에는 균열이 주응력 방향으로 놓이게 되므로 전체 도메인의 변형에 균열의 존재가 미치는 영향은 크지 않았다. 균열전단변위는 0.1 μm 이하의 미미한 수준이며, 전단력에 의해 균열을 따라 발생한 유의미한 변위가 아닌 평형해석 과정에서 블록들의 움직임으로 인해 랜덤하게 발생한 변위로 판단된다. 경사각이 30°, 45°, 60°인 경우에는 응력 조건(σy>σx)으로 인해 균열 상부면이 하부면을 타고 미끄러지는 전단방향을 보였다. 최대전단변위는 균열의 중심에서 발생하였으며, 경사각이 45°인 경우 최대 16.6 μm로 경사각이 30°와 60°에 비해 높은 전단변위를 보였다.

/media/sites/ksrm/2021-031-04/N0120310404/images/ksrm_31_04_04_F4.jpg
Fig. 4

Rock displacement magnitude (unit: m) and maximum fracture shear displacement (at fracture center): (a) α = 0°, (b) α = 30°, (c) α = 45°, (d) α = 60°, (e) α = 90°

Fig. 5에서의 전단응력 σxy는 XY 평면에서의 전단응력을 의미하므로 균열의 전단응력과는 다른 의미이다. 하지만 이를 통해 균열의 미끄러짐으로 인해 유도되는 균열 주변부에서의 응력 집중과 방향성 등을 살펴볼 수 있다. 전단응력 분포를 살펴보면 균열의 양 끝단 주변에서 응력 집중이 발생하는 것을 확인할 수 있는데, 이는 균열을 따라 미끄러짐이 발생하지만 암석은 파괴되지 않고 끝단에서 변위가 구속되기 때문이다. 좌측 끝단의 상부면(또는 우측 끝단의 하부면)에서는 음의 전단응력이 좌측 끝단의 하부면(또는 우측 끝단의 상부면)에서는 양의 전단응력이 관찰되는데, 이는 균열 양 끝단 주위에서 전단방향을 따라 각각 인장력과 압축력이 작용하고 있음을 간접적으로 보여준다.

/media/sites/ksrm/2021-031-04/N0120310404/images/ksrm_31_04_04_F5.jpg
Fig. 5

Rock shear stress (σxy) distribution (unit: Pa): (a) α = 0°, (b) α = 30°, (c) α = 45°, (d) α = 60°, and (e) α = 90°

수치해석을 통해 얻어진 균열 전단변위를 해석해와 비교하기 위해서는 임의 방향의 균열에 작용하는 이론적인 응력을 계산할 필요가 있으며, 2차원 요소의 응력해석을 통해 기울어진 균열에 작용하는 수직응력과 전단응력을 수식화하면 다음과 같다.

(7)
σn=σx+σy2+σx-σy2cos2θ
(8)
τ=-σx-σy2sin2θ

여기서 σnτ는 균열에 수직응력과 전단응력, σxσy는 X, Y 방향 수직응력, θ는 기울어진 균열의 법선 벡터가 X축과 반시계방향으로 이루는 각도이다.

균열의 경사각(α=90-θ)에 따라 식 (8)에 의해 계산한 전단응력을 식 (2)에 대입하면, 경계응력으로 인해 균열에 유도되는 전단변위 해석해를 구할 수 있다. Table 3은 경사에 따른 이론적인 균열 전단응력, 해석해에 따른 균열의 전단변위, 수치해석을 통해 구해진 균열의 전단변위를 정리한 것이다. 표에서 전단변위는 모두 균열 중심에서 얻어진 값이다.

Table 3.

Theoretical shear stress and shear displacement of inclined rock fracture

Fracture inclination angle,
(α=90-θ)
Theoretical fracture shear stress
(MPa)
Theoretical fracture shear
displacement at center (μm)
Numerically obtained fracture shear
displacement at center (μm)
30° 2.17 13.73 14.04
45° 2.50 15.85 16.57
60° 2.17 13.73 14.28

Fig. 6은 해석해와 수치모델에서 구해진 균열 전단변위를 비교하여 도시한 것이다. 그림에서 수치해석 결과(scatter plot)는 균열면 위의 접촉(contact)에서 계산된 결과이며, 실선은 해석해이다. 각 그림에서 r은 균열 중심으로부터의 거리를 나타내며, 균열 중심을 기준으로 좌표가 1, 4분면에 있으면 양의 부호를, 2, 3분면에 있으면 음의 부호를 부여하였다. 비교를 통해 알 수 있듯이, 균열의 경사각이 45°일 때 수치해석모델에서 전단변위가 약간 크게 계산된 것을 제외하면 모든 해석케이스에서 거의 일치하는 결과를 보였다. 이는 선행연구(Park et al., 2020)에서 보고한 균열 유도응력 불일치, 도메인 경계효과 문제를 상당 부분 개선한 결과로서, 본 연구의 새로운 해석모델이 암석에 내포된 균열의 전단거동을 매우 합리적으로 모사할 수 있음을 보여준다.

/media/sites/ksrm/2021-031-04/N0120310404/images/ksrm_31_04_04_F6.jpg
Fig. 6

Comparison of analytically and numerically calculated fracture shear displacements: (a) α = 30°, (b) α = 45°, and (c) α = 60°

나) 거친 균열

거친 균열에 대한 벤치마크 모델링에서는 각 참여팀의 해석모델이 균열의 기하학적 특성으로 인한 거칠기 효과를 모사할 수 있는지를 검토한다. 거친 균열의 모델링을 위해 동일한 암석 절리면에서 채취한 세 개의 2차원 거칠기 프로파일을 사용된다. 본 연구에서는 각 프로파일에 대하여 2, 5, 10 mm의 측정간격(해상도)를 적용하여 총 9가지의 해석 모델을 작성하였다.

Fig. 7은 1번 프로파일에 대해 서로 다른 해상도를 적용하여 생성한 거친 균열 모델을 보여준다. 해석모델의 등가직경이 3.2 mm이므로, 2 mm의 해상도를 적용한 경우 균열 주위의 그리드 형태가 불규칙하고 블록 자체의 그리드 밀도도 다소 높게 생성된다. 그러나 시료 전체의 변형이 대부분 입자의 변형에 의해 제어되는 미시물성 조합(Table 2의 Case 1)을 사용하였으므로, 블록 그리드 밀도가 해석 결과에 미치는 영향은 크지 않았을 것으로 판단된다.

/media/sites/ksrm/2021-031-04/N0120310404/images/ksrm_31_04_04_F7.jpg
Fig. 7

Representation of rough fracture (profile 1) with different resolutions of (a) 2 mm (b) 5 mm and (c) 10 mm

Fig. 8Fig. 7에 나타낸 세 가지 모델(profile 1, 세 가지 해상도)의 해석 결과로서 균열 전단변위의 크기를 color scale로 표현한 것이다. 편평한 면의 해석결과(Fig. 4, Fig. 6)와 달리, 거친 균열에서는 모두 좌측 상단에서 최대전단변위를 보였으며 비대칭적인 결과를 보였다. 2 mm 해상도의 균열 모델에서는 최대전단변위가 7.1 μm로, 편평한 균열(α=60°) 모델의 1/2 수준에 지나지 않았으며, 프로파일의 해상도가 커질수록 전단변위는 증가하고 편평한 면과 유사한 해석 결과를 보였다.

/media/sites/ksrm/2021-031-04/N0120310404/images/ksrm_31_04_04_F8.jpg
Fig. 8

Shear displacement of rough fracture - profile 1 (unit: m): roughness resolutions of (a) 2 mm (b) 5 mm and (c) 10 mm

Fig. 9는 9가지 거친 균열 모델의 전단변위를 해상도에 따라 분류해 도시한 것으로 각 그림에는 편평한 균열 모델에 대한 해석해를 함께 도시하였다. 수치해석 결과(scatter plot)는 거친 균열을 따르는 접촉점(contact)에서 계산된 값이며, 실선은 60° 경사의 편평한 면에 대한 해석해이다. 모든 해석 케이스에서 거친 균열은 편평한 면에 비해 낮은 전단변위를 보였다. 세 프로파일에서 모두, 해상도가 증가함에 따라 전단변위가 증가하는 경향을 보였고 결과의 비대칭성은 감소하였다. 프로파일 간격이 2 mm일 때에는 세 균열 모델이 다소 상이한 결과를 나타내었으나, 간격이 증가함에 따라 거의 일치하는 경향을 보였고 10 mm 해상도를 적용한 경우 편평한 면과 거의 유사한 전단변위가 계산되었다.

/media/sites/ksrm/2021-031-04/N0120310404/images/ksrm_31_04_04_F9.jpg
Fig. 9

Comparison of numerically estimated rough fracture shear displacements (roughness resolutions of (a) 2 mm, (b) 5 mm, and (c) 10 mm) and analytical solution for planar fracture

Fig. 10은 프로파일 1을 이용한 해석모델의 균열 주변 전단응력(σxy) 분포를 보여주는 것으로, 비교를 위해 60° 경사각의 편평한 균열 모델의 결과를 함께 표시하였다. 프로파일의 해상도를 2 mm로 적용한 경우, asperity를 따라 응력집중이 균열면 전체에 걸쳐 발생하지만, 상대적으로 소규모의 응력집중영역이 불연속적으로 발달하였다. 해상도가 증가할수록 응력의 공간적 연속성이 커지며, 응력집중이 발생하는 영역의 크기도 증가하는 경향성을 보였다. 또한, 해상도가 10 mm인 경우 편평한 면과 거의 유사한 응력분포를 나타내었다.

/media/sites/ksrm/2021-031-04/N0120310404/images/ksrm_31_04_04_F10.jpg
Fig. 10

Shear stress (σxy) of rock around fracture (unit: Pa): rough fracture models (profile 1) with different resolutions of (a) 2 mm (b) 5 mm and (c) 10 mm and (d) planar fracture model

위에서 살펴본 해석결과는 모두 입자기반해석모델을 통해 구성된 균열 모델이 거칠기의 영향을 적절하게 모사하고 있음을 보여준다. 거칠기를 갖는 균열에서는 asperity의 맞물림으로 인해 편평한 균열에 비해 전단저항이 증가하고 전단변위가 감소하게 된다. 일반적인 암석역학에서 암석 절리면의 거칠기를 전단강도와 전단강성에 반영하는 것도 이와 같은 이유이다. 거칠기 프로파일의 측정간격이 감소할수록 전단변위가 감소하고 프로파일간의 결과 차이가 증가하는 경향성 역시 거칠기 효과가 적절히 반영됨을 시사한다. 자연적 연속성을 갖는 암석 표면의 특성을 일정 간격을 통해 측정하고 재현한다는 것은, 그 간격보다 작은 규모에서 asperity가 갖는 기하학적, 역학적 특성을 배제함을 전제한다. 즉, 작은 간격의 프로파일을 적용할수록 거칠기가 더 큰 균열모델을 생성하게 되며, Fig. 8, Fig. 9, Fig. 10의 결과에 이러한 거칠기 영향이 적절히 반영되어 있음을 확인할 수 있다.

4.2 수리-역학 조건 벤치마크 모델링

4.2.1 해석개요

수리-역학 조건의 벤치마크 모델링은 각 참여팀의 해석모델이 유체 주입 압력으로 인한 균열의 개방(fracture opening)을 적절히 모사할 수 있는지 검증하는 데에 목적이 있다. 해석모델의 형상과 물성은 역학 조건 벤치마크 모델과 동일하며, 균열의 경우 30° 경사각을 갖는 편평한 면으로 가정한다. 경계조건으로서 1) HM-1: σx= 0, σy = 0 조건과 2) HM-2: σx= 5 MPa, σx= 10 MPa 조건 두 가지가 고려되며, 균열에는 2 MPa의 일정한 유체 주입 압력이 작용하는 조건을 가정한다. 수리-역학 해석 시 암석은 불투수성으로 유동은 균열을 통해서만 발생한다.

균열의 개방은 균열에 작용하는 수직응력과 관련이 있다. 식 (7)에 따르면, 두 번째 응력조건에서는 균열 수직응력이 압축 방향으로 작용하게 되므로 균열의 개방을 기대하기 어렵다. 따라서 본 연구에서는 Task G에서 주어진 조건 이외에 추가적인 해석조건(HM-3)을 설정하여 두 번째 응력조건에서 균열의 미끄러짐과 개방을 동시에 재현하고자 하였다.

Table 4는 본 연구에서 고려한 세 가지 해석조건(HM-1, HM-2, HM-3)을 요약한 것이다. 추가 해석조건인 HM-3에서는 HM-2와 동일한 응력을 적용하되, 10.75 MPa의 압력(P)을 적용하여, 이론적으로 균열에 2 MPa의 인장응력이 작용할 수 있도록 설정하였다. 표에는 식 (7)식 (8)에 의해 계산한 균열 수직응력과 전단응력을 함께 정리하였다. 이론적으로 HM-3 모델에서 균열에 작용하는 유효수직응력(인장응력)과 전단응력은 각각 HM-1 모델의 수직응력, HM-2 모델의 전단응력과 동일하게 된다.

Table 4.

Hydro-mechanics benchmark model condition

Simulation case Boundary stress and pressure (MPa) Theoretical fracture stress (MPa)
σx σy P σn τ
HM-1 0 0 2 -2 0
HM-2 5 10 2 8.75 2.17
HM-3 5 10 10.75 -2 2.17

3DEC 내에서는 암석 블록 내 유동(matrix flow), 균열 내 유동(fracture flow), 블록과 균열 간 유동(leak off 등) 해석이 모두 가능하다. 그러나 수리-역학 벤치마크 모델링은 일정한 유체 압력으로 인한 균열의 변위를 해석하는 데에 초점이 맞추어져, 유동해석을 수행하지 않고, 압력을 고려한 유효응력 해석을 수행하였다. 수리역학 벤치마크 모델링에서 초기응력은 경계응력과 동일한 것으로 가정하였고, 유체의 압력에 따른 균열의 수직, 전단방향 거동을 살펴보았다.

4.2.2 해석결과

Fig. 11은 각 해석모델의 변위 분포를 보여주는 것으로 암석시료 내 변위 크기 분포도와 균열 주변 변위 벡터를 확대하여 도시한 것이다. 전반적인 결과를 살펴보면, 균열 중심에서의 변위가 가장 크며, 대칭적인 분포를 나타내었다. 역학 조건 벤치모델링에서는 초기응력이 0인 상태에서 경계응력의 재하로 인한 균열의 거동을 살펴보았으므로, 경계에서 모델의 중심으로 갈수록 작은 변위를 나타내었다. 그러나 수리-역학 벤치마크 모델에서는 경계응력 조건과 동일한 초기응력 상태에서 유체 주입 압력에 의한 균열 거동을 살펴보았으므로, 도메인 중심에서 가장 큰 변위를 보이게 된다. 한편, 균열의 변위의 크기는 상부면과 하부면의 상대변위를 의미하므로 Fig. 11에서 균열 주변에서 관찰되는 변위의 2배와 동일하게 된다.

HM-1 모델에서는 초기응력이 0이므로 유체 압력에 상응하는 2 MPa의 인장응력이 균열의 수직방향으로 작용하고, 전단응력이 작용하지 않으므로 균열 개방 현상만이 관찰된다. HM-2 모델에서는 균열의 수직방향으로 압축응력이 작용하므로 균열이 오히려 닫히게 되고(fracture closing), 전단응력에 의해 균열 미끄러짐(fracture slip)이 발생한다. 한편, HM-3 모델에서는 균열에 높은 압력을 적용하여 인장력을 유도함으로써 균열 개방과 미끄러짐이 동시에 유도되는 것을 확인할 수 있다.

균열 수직변위의 해석해(식 (1))를 유체 압력을 고려하여 보다 일반적인 형태로 표현하면 식 (9)와 같이 수식화할 수 있다(Parker, 1981). 유체 압력은 전단응력에 영향을 미치지 않으므로 전단변위에 대한 해석해는 식 (2)와 동일하다.

(9)
un=|σn-P|4(1-ν2)Ea2-r2

Fig. 12는 각 해석모델의 수치해석 결과(균열 수직변위, 전단변위)를 해석해(식 (9), 식 (2))와 비교한 그림으로, 수치해석 모델에서 수직변위가 약간 크게 계산된 것을 제외하면 모든 해석 케이스에서 수직변위과 전단변위가 해석해와 거의 일치하는 것으로 나타났다.

/media/sites/ksrm/2021-031-04/N0120310404/images/ksrm_31_04_04_F11.jpg
Fig. 11

Contour plot of displacement (at block center, unit: m) and maximum fracture shear displacement (at fracture center, unit: μm): (a) HM-1, (b) HM-2, and (c) HM-3

/media/sites/ksrm/2021-031-04/N0120310404/images/ksrm_31_04_04_F12.jpg
Fig. 12

Comparison of analytically and numerically calculated fracture shear displacements: (a) HM-1, (b) HM-2 and (c) HM-3

5. 요약 및 결론

본 논문에서는 국제공동연구인 DECOVALEX-2023 프로젝트 Task G의 일환으로 수행된 벤치마크 모델링 해석 결과를 소개하였다. 벤치마크 모델링의 주요내용은 암석 내부에 존재하는 단일 균열의 역학적, 수리-역학적 거동을 모사하기 위한 수치해석 기법을 개발하고, 이를 해석해와 비교함으로써 해석기법을 검증하는 것이다. 본 연구에서는 사면체 입자 모델과 3DEC를 이용한 입자기반 개별요소모델을 이용하여 암석 균열을 모델링할 수 있는 해석기법을 개발하였으며, 이 과정에서 등가연속체 접근법을 적용한 새로운 미시물성 산정기법을 제안하였다. 개발된 수치모델은 경계응력에 따른 균열의 미끄러짐(fracture slip)과 유체 주입 압력에 따른 균열의 개방(fracture opening), 균열 경사에 따른 응력 분포, 거칠기로 인한 전단저항의 증가와 전단변위의 구속 등을 합리적으로 재현할 수 있는 것으로 나타났으며, 모든 해석케이스에서 모델의 거동은 해석해와 거의 일치하는 거동을 보였다.

본 연구에서 개발 중인 입자기반 개별요소모델은 암석 내 존재하는 균열의 성장, 새로운 균열의 개시와 전파 특성 등을 모사하는 데에 특화된 강점이 있다. 그러나 벤치마크 모델링에서는 암석 균열을 탄성재료에 포함된 강도가 없는 불연속면으로 간주하였으므로, 실제 암석 균열의 거동 특성을 반영한 심도 깊은 논의가 이루어지지는 못하였다. 현재 국제공동연구가 초기 연구개발단계에 있고, 열-역학 조건 벤치마크 모델링, 실내시험 모델링, 현장조건 시뮬레이션 등 체계적인 연구 계획이 수립되어 있어, 본 연구에서 제시한 해석모델을 단계적으로 개선, 개발하고 다양한 문제에 적용할 수 있을 것으로 기대된다.

Acknowledgements

The authors appreciate and thank the DECOVALEX-2023 Funding Organizations, Andra, BASE, BGE, BGR, CAS, CNSC, COVRA, US DOE, ENRESA, ENSI, JAEA, KAERI, NWMO, RWM, SÚRAO, SSM and Taipower for their financial and technical support of the work described in this paper. The statements made in the paper are, however, solely those of the authors and do not necessarily reflect those of the Funding Organizations.

This research was supported by the Basic Research Project of the Korea Institute of Geoscience and Mineral Resources (GP2020-010) funded by the Ministry of Science and ICT, Korea. This work was also supported by the Institute for Korea Spent Nuclear Fuel (iKSNF) and National Research Foundation of Korea (NRF) grant (2021M2E1A1085193) funded by the Ministry of Science and ICT, Korea.

References

1
Asadi, M.S., Rasouli, and V., Barla, G., 2012, A bonded particle model simulation of shear strength and asperity degradation for rough rock fractures. Rock Mech. Rock Eng. 45(5), 649-675. 10.1007/s00603-012-0231-4
2
Birkholzer, J.T. and Bond, A.E., 2020, DECOVALEX-2023 Introduction. Presented in DECOVALEX-2023 2nd Workshop, November 16-20 2020, Virtual conference.
3
Birkholzer, J.T., Tsang, C.F., Bond, A.E., Hudson, J.A., Jing, L., and Stephansson, O., 2019, 25 years of DECOVALEX-Scientific advances and lessons learned from an international research collaboration in coupled subsurface processes. International Journal of Rock Mechanics and Mining Sciences 122, 103995. 10.1016/j.ijrmms.2019.03.015
4
Castro-Filgueira, U., Alejano, L.R., and Ivars, D.M., 2020, Particle flow code simulation of intact and fissured granitic rock samples. Journal of Rock Mechanics and Geotechnical Engineering 12(5), 960-974. 10.1016/j.jrmge.2020.01.005
5
Cho, N., Martin, C.D., and Sego, D.S., 2007, A clumped particle model for rock. Int. J. Rock Mech. Min. Sci. 44, 997-1010. 10.1016/j.ijrmms.2007.02.002
6
Cundall, P.A., 1988, Formulation of a three-dimensional distinct element model-Part I. A scheme to detect and represent contacts in a system composed of many polyhedral blocks, In International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts 25(3), 107-116. 10.1016/0148-9062(88)92293-0
7
Fu, T.F., Xu, T., Heap, M.J., Meredith, P.G., and Mitchell, T.M., 2020, Mesoscopic time-dependent behavior of rocks based on three-dimensional discrete element grain-based model. Computers and Geotechnics 121, 103472. 10.1016/j.compgeo.2020.103472
8
Geuzaine, C. and Remacle, J.F., 2009, Gmsh: A 3‐D finite element mesh generator with built‐in pre‐and post‐processing facilities. International journal for numerical methods in engineering 79(11), 1309-1331. 10.1002/nme.2579
9
Ghazvinian, E., Diederichs, M.S., and Quey, R., 2014, 3D random Voronoi grain-based models for simulation of brittle rock damage and fabric-guided micro-fracturing. Journal of Rock Mechanics and Geotechnical Engineering 6(6), 506-521. 10.1016/j.jrmge.2014.09.001
10
Gu, D. and Huang, D., 2016, A complex rock topple-rock slide failure of an anaclinal rock slope in the Wu Gorge, Yangtze River, China. Engineering Geology 208. 165-180. 10.1016/j.enggeo.2016.04.037
11
Hu, Z., Xu, N., Li, B., Xu, Y., Xu, J., and Wang, K., 2020, Stability Analysis of the Arch Crown of a Large-Scale Underground Powerhouse During Excavation. Rock Mechanics and Rock Engineering 53(6), 1-9. 10.1007/s00603-020-02077-4
12
Itasca Consulting Group Inc., 2014, UDEC (Universal Distinct Element Code) version 6.0. Minneapolis: Itasca CGI.
13
Itasca Consulting Group Inc., 2017, 3DEC (3 Dimensional Distinct Element Code) version 5.2. Minneapolis: Itasca CGI.
14
Itasca Consulting Group Inc., 2021, PFC (Particl Flow Code) online manual. Available at https://www.itascacg.com/software/PFC.
15
Lan, H.X., Martin, C.D., and Hu, B., 2010, Effect of heterogeneity of brittle rock on micromechanical extensile behavior during compression loading. Journal of Geophysical Research: Solid Earth 115, B1. 10.1029/2009JB006496
16
Li, A., Dai, F., Xu, N., Gu, G., and Hu, Z., 2019, Analysis of a complex flexural toppling failure of large underground caverns in layered rock masses. Rock Mechanics and Rock Engineering 52(9), 3157-3181. 10.1007/s00603-019-01760-5
17
Li, Z. and Rao, Q.H., 2021, Quantitative determination of PFC3D microscopic parameters. Journal of Central South University, 28(3), 911-925. 10.1007/s11771-021-4653-6
18
Lisjak, A., Tatone, B.S., Grasselli, G., and Vietor, T., 2014, Numerical modelling of the anisotropic mechanical behaviour of opalinus clay at the laboratory-scale using FEM/DEM. Rock mechanics and rock engineering 47(1), 187-206. 10.1007/s00603-012-0354-7
19
Ma, Q., Tan, Y.L., Liu, X.S., Zhao, Z.H., and Fan, D.Y., 2021, Mechanical and energy characteristics of coal-rock composite sample with different height ratios: a numerical study based on particle flow code. Environmental Earth Sciences 80(8), 1-14. 10.1007/s12665-021-09453-5
20
Mahabadi, O.K., Randall, N.X., Zong, Z., and Grasselli, G., 2012, A novel approach for micro-scale characterization and modeling of geomaterials incorporating actual material heterogeneity. Geophysical Research Letters 39(1), 1-6. 10.1029/2011GL050411
21
Park, J.W., Park, C.H., Yoon, J., and Lee, C., 2020, Grain-based distinct element modelling of the mechanical behavior of a single fracture embedded in rock: DECOVALEX-2023 Task G (Benchmark simulation). Tunnel & Underground Space 30(6), 573-590.
22
Park, J.W., Park, C., Song, J.W., Park, E.S., and Song, J.J., 2017, Polygonal grain-based distinct element modeling for mechanical behavior of brittle rock. International Journal for Numerical and Analytical Methods in Geomechanics 41(6), 880-898. 10.1002/nag.2634
23
Park, J.W. and Song, J.J., 2009, Numerical simulation of a direct shear test on a rock joint using a bonded-particle model. Int. J. Rock Mech. Min. Sci. 46, 131-1328. 10.1016/j.ijrmms.2009.03.007
24
Parker, A.P., 1981. The mechanics of fracture and fatigue: An introduction. London and New York, E. & F. N. Spon, Ltd., pp. 174.
25
Pollard, D.D. and Segall, P., 1987, Theoretical displacements and stresses near fractures in rock: with applications to faults, joints, veins, dikes, and solution surfaces. In: Fracture mechanics of rock. pp. 277-347. 10.1016/B978-0-12-066266-1.50013-2
26
Potyondy, D.O., 2010, A grain-based model for rock: approaching the true microstructure, In: Proceedings of Bergmekanikk I Norden. Kongsberg, Norway.
27
Potyondy, D.O., 2012, A flat-jointed bonded-particle material for hard rock, In: Proceedings of 46th U.S. Rock mechanics/ geomechanics symposium, Chicago, USA, ARMA 12-501.
28
Potyondy, D.O., 2013, PFC3D flat joint contact model version1, Itasca Consulting Group, Minneapolis, Technical Memorandum ICG7234-L.
29
Potyondy, D.O., 2015, The bonded-particle model as a tool for rock mechanics research and application: current trends and future directions. Geosystem Eng. 18, pp. 1-28. 10.1080/12269328.2014.998346
30
Potyondy, D.O. and Cundall, P.A., 2004, A bonded-particle model for rock, Int. J. Rock Mech. Min. Sci. 41(8), 1329-1364. 10.1016/j.ijrmms.2004.09.011
31
Potyondy, D.O. and Hazzard, J.F., 2008, Effects of stress and induced cracking on the static and dynamic moduli of rock, In: Proceedings of First International FLAC/DEM Symposium, Minneapolis, USA, pp. 147-156.
32
Wang, Z., Wang, T., Wu, S. and Hao, Y., 2021, Investigation of microcracking behaviors in brittle rock using polygonal grain‐based distinct method. International Journal for Numerical and Analytical Methods in Geomechanics. 10.1002/nag.3246
33
Wu, Z., Ma, L., and Fan, L., 2018, Investigation of the characteristics of rock fracture process zone using coupled FEM/DEM method. Engineering Fracture Mechanics 200, 355-374. 10.1016/j.engfracmech.2018.08.015
34
Yoon, J., Zimmermann, G., Zang, A., and Stephansson, O., 2015, Discrete element modeling of fluid injection-induced seismicity and activation of nearby fault 1. Can. Geotech. J. 52(10), 1457-1465. 10.1139/cgj-2014-0435
페이지 상단으로 이동하기