Numerical Analysis of Blast Induced Fracturing of Hard Rocks - PowerPoint PPT Presentation

Hotel Caesar Park, Rio de Janeiro, Brazil August 2-6, 2004 Numerical Analysis of Blast Induced Fracturing of Hard Rocks Araken Lima,Celso Romanel and Deane Roehl Civil Engineering Department

Motivation • Rock blasting plan • Prediction of fracture extention • Formation of blocks

Rock blasting 1 p 1 Shock wave energy for driving the 1 detonation in the explosive 2 Kinetic and strain shock wave energy. Blasthole pressure p 2 3 Energy release during crack propagation 4 Noise, heat and other wasted energy 2 4 3 p 3 AB C D Volume of expanding gases or time of gas expansion Energy partitioning during rock blasting (Whittaker, Singh e Sun [1992]). Detonation energy = stress waves + gas pressure + others (temperature, flying rock fragments, air displacements) Explosives: TNT – high stress waves and low gas production ANFO – high gas production and low energy in stress waves shock wave propagation - µ s (dynamic) gas pressurization - ms (quasi-static)

Rock blasting Blast induced P-waves •crushing around the borehole •dense radial cracks (4-8 radii) •dominant cracks propagate Number of dominant cracks Ghosh e Daemen (1995) - 8 a 12 Song e Kim (1995) - 10 a 12

Numerical Model 2D Adaptive finite element model • Loading: blast induced stress waves • Rock material: homogeneous, isotropic linear elastic up to breakage • Crushed and dense fracture zones are neglected • Mixed mode fracture (P- and S-waves)

Fracture model •Mixed mode I-II fracture criteria (static loading) •Dynamic fracture criteria - mode I Grady e Lipkin (1980)

Fracture model •Crack propagation direction – maximum tensile tangential stress   2 2     K K 1 1   I I     θ = ± + 2 arctan 8     m   4 K 4 K     II II     Mode I K I = K IC , K II = 0 → θ m =0 º Mode II K I = 0 , K II = K IIC → θ m = -70,53º •Crack growth in a time interval ∆ a = c . ∆ t c crack growth velocity (smaller than half shear wave speed Grady & Kipp, 1979)

Finite element model Adaptive finite element mesh • Mesh generation : recursive spatial enumeration techniques quadtree inner domain boundary triangulation (Araújo 1999) • Maximum element size Wave reflection on element interface due to change in element size Efficient frequency transmission: Celep e Bazant [1983] and Mullen e Belytschko [1982]: 1/10 smallest wave length • State variable mapping { } [ ] { } n − ∆ t = ξ P η P t t u N ( , ) . q [ ] { { } } n − ∆ = ξ η t � P P t t � u N ( , ) . q [ ] { { } } n − ∆ t = ξ P η P t t � � � � u N ( , ) . q

Finite element model •Singular quarter point elements Henshell & Shaw (1975) Barsoum (1976) y Crack 3l/4 l/4 l/4 3l/4 tip •Finite element rosette 8 elements j x i •Stress intensity factors Displacement correlation technique (Shih, de Lorenzi and German, 1976) Modified crack closure method (Raju, 1987)

Finite element model Stress intensity factors Displacement correlation technique (Shih, de Lorenzi and German, 1976) y µ ⋅ π ( )   2 = ⋅ ⋅ ⋅ − − K   4 v v − I j 1 j 2 κ +  1  L j j-2 j+2 j+1 y x j-1 µ ⋅ π ( )   2 = ⋅ ⋅ ⋅ − − K   4 u u − II j 1 j 2 κ +  1  L x L L element rosette at crack tip

Finite element model Stress intensity factors Modified crack closure method (Raju, 1987) { } 1 [ ( ) ( ) ] [ ( ) ( ) ] = − ⋅ − + − + ⋅ − + − G F t v v t v v F t v v t v v ′ ′ ′ ′ I y 11 m m 12 l l y 21 m m 22 l l ⋅ δ 2 a i j { } 1 [ ( ) ( ) ] [ ( ) ( ) ] = − ⋅ − + − + ⋅ − + − G F t u u t u u F t u u t u u ′ ′ ′ ′ II x 11 m m 12 l l x 21 m m 22 l l ⋅ δ 2 a i j π 1 = π − t 6 20 = − 21 = t 6 3 t 12 11 2 2 22 = t 1 δ a Normal stress distribution at crack tip Nodal forces at singular elements

Finite element model Fracture closure control Contact - penalty formulation Master node Node to edge contact 1 g(u) 3 Linearized surface gap function Element boundary 2 ⋅ + ⋅ + A x B y C ( ) = g u 0 0 + 2 2 A B Fracture surface penetration Normal to contact forces [ ] + ⋅ = K K d R Equilibrium ( ) c ≤ 0 ⇒ = − g(u) p λ g u ( ) nc ∑ = λ ⋅ ∇ ⋅ ∇ T with K g g c j j = j 1

Dynamic model •Numerical integration of the equation of motion [ ] { } [ ] { } [ ] { } { ( ) } ⋅ + ⋅ + ⋅ = � � � M w C w K w P t Wilson θ implicit algorithm •Pressure pulse on hole wall- Duvall (1953) ( ) = ⋅ − α ⋅ − − β ⋅ t t p ( t ) p e e 0 rock constants α and β (Aimone, 1982).

Dynamic model Viscous dampers on the boundary- radiation condition σ = ⋅ ρ ⋅ ⋅ � a c w L σ = ⋅ ρ ⋅ ⋅ L � a c u explosion τ = ⋅ ρ ⋅ ⋅ � b c u T u τ = ⋅ ρ ⋅ ⋅ � b c w T Stress waves + + = � � � M u C u Ku R w τ = ⋅ ρ ⋅ ⋅ σ = ⋅ ρ ⋅ ⋅ T � L � b c u a c w velocity (m/s) time (s) time (s)

Dominant Cracks Detonation hole with 4, 8, 10, 12 and 16 dominant cracks Gosh/Daemem (1995) 8-12 cracks Song/Kim (1995) 10-12 cracks •hole size a 0 =5 cm •material properties: ρ = 28 Mg/m 3 , E = 41 GPa, ν = 0,25 K ID = 1,65 MPa ⋅ m 1/2 K IID = 1,03MPa ⋅ m 1/2 fracture propagation velocity: 1210 m/s •t=100 µs hole pressure 715, 8 MPa •finite element meshes 16 dominant fractures 4 dominant fractures

Dominant Cracks Principal stresses 4 cracks σ 1 σ 3 8 cracks σ 3 σ 1

Dominant Cracks Principal stresses 10 cracks σ 1 σ 3 12 cracks σ 3 σ 1

Dominant Cracks Principal stresses 16 cracks σ 1 σ 3

Examples • 1 Blasthole with 9 dominant cracks • 1 Blasthole with 8 dominant cracks • 2 Blastholes with 8 dominant cracks Granite: Borehole: E = 60 GPa radius = 2.54 cm ν = 0.25 ρ = 2.80 Mg/m 3 K ID = 1.65 MPa.m 1/2 penalty parameter 10 7 K IID =1.03 MPa.m 1/2 penetration tol 10 3 initial crack length 1/3 hole radius

Example 1: Blasthole with 9 dominant cracks Time increments ∆ t = 20 µ s to 50 µ s 40 cm free surface Blasthole 90 o 135 o 45 o 157.5 o 0 0 0 o 202.5 o 225 o 315 o 270 o

Blasthole with 9 dominant cracks t = 3 µ s 4811 nodes and 2292 elements σ 1 σ 1

Blasthole with 9 dominant cracks t = 17 µ s 4811 nodes and 2292 elements σ 1 σ 1

Blasthole with 9 dominant cracks t = 37 µ s 4783 nodes and 2260 elements σ 1 σ 1

Blasthole with 9 dominant cracks t = 182 µ s 6942 nodes and 3277 elements free surface σ 1

Blasthole with 9 dominant cracks t = 297 µ s 8266 nodes and 3891 elements free surface free surface σ 1

Blasthole with 9 dominant cracks t = 330 µ s 7127 nodes and 3325 elements

Blasthole with 9 dominant cracks t = 405 µ s 7777 nodes and 3637 elements

Blasthole with 9 dominant cracks t = 545 µ s 8995 nodes and 4199 elements

Blasthole with 9 dominant cracks t = 680 µ s 10047 nodes and 4687 elements

Blasthole with 9 dominant cracks t = 970 µ s 14154 nodes and 6594 elements Plate border cracks Experimental resuts on glass plate Porter(1970)

Example 2: Blasthole with 8 dominant cracks Time increments ∆ t = 50 µ s to 125 µ s 40 cm free surface Blasthole 90 o 135 o 45 o 0 0 180 o 0 o 225 o 315 o 270 o

Blasthole with 8 dominant cracks t = 15 µ s 3 4 2 5 1 6 8 7 σ 1 σ 1

Blasthole with 8 dominant cracks t = 523 µ s (without fracture interaction control) free surface free surface Stress σ 1 Horizontal displacement

Blasthole with 8 dominant cracks t = 523 µ s (with fracture interaction control – penalty method) free surface free surface Horizontal Stress σ 1 displacement

Blasthole with 8 dominant cracks t = 523 µ s No contact control Contact control

Example 3: 2 Blastholes with 8 dominant cracks Simultaneous blasting Blasthole 2 free surface 40 cm Blasthole 1 40 cm Crack propation speed: 1464 m/s

2 Blastholes with 8 dominant cracks t = 217,7 µ s t = 289,7 µ s 13 free surface 14 free surface 16 15 4 2 3 5 Crack pattern

2 Blastholes with 8 dominant cracks Adaptive mesh t = 538,2 µ s free surface

2 Blastholes with 8 dominant cracks Crack pattern t = 538,2 µ s free surface

Example 4: 2 Blastholes with 8 dominant cracks Sequential blasting Blasthole 2 free surface ∆ t=123.9 µ s 40 cm Blasthole 1 40 cm Crack propation speed: 1464 m/s

2 Blastholes with 8 dominant cracks Crack pattern t = 123,9 µ s Sequential blasting free surface 1 2 1 2

2 Blastholes with 8 dominant cracks 5 1 13 t=199,4 µs t=285,4 µs 2 12 1 13 14 t=311,9 µs t=439,9 µs

2 Blastholes with 8 dominant cracks t=285,4 µ s t=199,4 µ s t=439,9 µ s t=311,9 µ s