Simulation of crack propagation in rock in plasma blasting technology

Digital Object Identiﬁer (DOI) 10.1007/s00193-002-0141-6

Shock Waves (2002) 12: 145–152

Simulation of crack propagation in rock

in plasma blasting technology

V. R. Ikkurthi, K. Tahiliani, S. Chaturvedi

Institute for Plasma Research, Bhat, Gandhinagar-382 428, Gujarat, India

Received 4 September 2001 / Accepted 12 March2002

Published online 11 June 2002 – c

Springer-Verlag 2002

Abstract. Plasma Blasting Technology (PBT) involves the production of a pulsed electrical discharge by

inserting a blasting probe in a water-ﬁlled cavity drilled in a rock, which produces shocks or pressure waves

in the water. These pulses then propagate into the rock, leading to fracture. In this paper, we present the

results of two-dimensional hydrodynamic simulations using the SHALE code to study crack propagation in

rock. Three separate issues have been examined. Firstly, assuming that a constant pressure P is maintained

in the cavity for a time τ, we have determined the P-τ curve that just cracks a given rock into at least

two large-sized parts. This study shows that there exists an optimal pressure level for cracking a given

rock-type and geometry. Secondly, we have varied the volume of water in which the initial energy E is

deposited, which corresponds to diﬀerent initial peak pressures P

peak

. We have determined the E-P

peak

curve that just breaks the rock into four large-sized parts. It is found that there must be an optimal P

peak

that lowers the energy consumption, but with acceptable probe damage. Thirdly, we have attempted to

identify the dominant mechanism of rock fracture. We also highlight some numerical errors that must be

kept in mind in suchsimulations.

Key words: Plasma blasting technology, Brittle fracture, Crack propagation

PACS: 46.50.+a, 47.11.+j, 47.40.−x, 52.50.Lp, 64.30.+t

1 Introduction

Plasma Blasting Technology (PBT) is emerging as an en-

vironment friendly alternative to the use of explosives in

rock-fracturing (Nantel and Kitzinger 1990). A schematic

of PBT is shown in Fig. 1. A bore hole is drilled in the

rock to nearly half its height and ﬁlled with water. A

blasting probe is dipped in the water and a high-energy

capacitor bank is discharged through it. The discharge

dumps energy in a small region which can be called the

“Hot Zone”, creating a steam/plasma bubble. This bub-

ble drives shocks or pressure pulses into the water. These

pulses propagate into the rock, leading to fracture.

The extent and nature of cracks would depend upon

the temporal variation of surface pressure on the cav-

ity walls. Earlier simulation work (Madhavan et al. 2000)

has shown that varying the spatial and temporal proﬁles

of electrical power deposition in the water can produce

wide variations in the surface pressures. These simulations

also show that rapid deposition of power produces shocks,

which yield generally higher surface pressure levels, albeit

Correspondence to: V.R. Ikkurthi

(e-mail: ramana@ipr.res.in)

An abridged version of this paper was presented at the 23rd

Int. Symposium on Shock Waves at Fort Worth, Texas, from

July 22 to 27, 2001.

000000

111111

000000

111111

0000

1111

0000

1111

000

111

BLASTING
PROBE

STEEL

TAMPER

WATER

ROCK

HOT−ZONE

Fig. 1. Schematic of plasma blasting arrangement

with rapid temporal variation. Slower deposition yields

lower, but longer-lasting pressure levels. Similarly, higher

pressure levels can be produced by dumping energy in a

smaller volume.

146

V.R. Ikkurthi et al.: Simulation of rock fracture

Faster energy deposition requires extra experimental

eﬀort, e.g., to reduce the inductance of the pulsed power

system. Also, dumping the energy in a smaller volume

requires blasting probes that can withstand higher pres-

sures. Hence, these can only be justiﬁed if they are more

eﬀective at rock fracturing. It is thus of interest to exam-

ine the eﬀectiveness of diﬀerent power deposition proﬁles

in rock fracturing.

We have carried out two-dimensional hydrodynamic

simulations using the SHALE code (Demuth et al. 1985) to

study rock fracture, using a crack growth and propagation

model.

2 Fracture model description

SHALE is an Arbitrary Lagrangian Eulerian Hydrody-

namic code developed at Los Alamos National Laboratory

(Demuth et al. 1985). The Bedded Crack Model (BCM)

(Demuth et al. 1985), already incorporated in SHALE,

has been used to model crack propagation. An outline of

BCM is given below for completeness.

BCM is a microphysical model and contains two lev-

els of description of the solids. Macroscopically, it sim-

ulates an equivalent elastic continuum, where the eﬀec-

tive moduli characterizing the equivalent continuum vary

with certain internal variables that are calculated from

the microstructure, characterized by a statistical ensem-

ble of penny shaped cracks. The internal variables upon

which the eﬀective moduli depend are calculated from the

state of the microstructure. The eﬀective moduli are used

to calculate the macroscopic stress ﬁeld from the strain,

which is determined external to the constitutive law. Fi-

nally, the macroscopic stress is used to calculate changes

in the microstructure.

The BCM implementation in SHALE is a limited one

and not the full model developed by Dienes (Dienes 1979),

in the sense that it does not handle the growth of cracks

with random orientations. SHALE-BCM was especially

developed for Oil Shale and does calculations only for

cracks parallel to the bedding planes and assumes no

cracks in other orientations. This is a good ﬁrst approx-

imation, since the bedding planes in Oil Shale are dom-

inant planes of weakness. Also, SHALE-BCM does not

deal with detailed physics of crack nucleation and coales-

cence. An enhanced version of BCM, BCM3, developed

at LANL (Adams et al. 1983), treats cracks in two or-

thogonal orientations besides the bedding planes. PBCM

(Margolin and Smith 1985), an enhanced version, includes

both compaction and brittle fracture. The brittle fracture

model, such as the one described in (Curran and Seaman

1987), is a very detailed one and deals with nucleation,

growth and coalescence of cracks with random orienta-

tions.

Crack growth, which takes place under certain condi-

tions, is an irreversible process and is governed by Grif-

ﬁth’s Growth Criterion (Griﬃth 1920). It states that a

crack will grow if the energy released exceeds the energy

required for the growth of the crack:

∂(W − S)

∂c

≥ 0

(1)

where ‘W’ is the Elastic Strain energy, ‘S’ is the increase

in the surface energy and ‘c’ is the crack radius.

All equations given in this section, and hence our sim-

ulations, are in a cylindrical coordinate system. When the

normal stress σ

is positive, i.e., the material is in ten-

sion, a crack whose normal points in the y-direction will

grow if

+ σ

2 (1 + ν)

≥

π T E

2 (1 − ν

) c

(2)

where ‘T’ is the energy required to create a new surface

and can be related to the critical strain energy release

rate, ν is the Poisson ratio of the initial matrix material,

‘E’ is Young’s modulus and σ

is the shear stress.

Equation (2) shows that in any stress ﬁeld, there is

a “critical crack size” found by assuming equality. Cracks

larger than the critical size are unstable and grow, whereas

smaller cracks are stable. The presence of shear stress

tends to decrease the critical crack size.

When σ

is negative, i.e., the material is in compres-

sion, the cracks tend to close. If this is taken into account,

the energy balance equation would have to include the

additional energy dissipated by friction between the crack

planes. The Griﬃth criterion would then take the form

∂(W − F − S)

∂c

≥ 0

(3)

where ‘F’ is the work done against friction (Margolin and

Smith 1985).

This means that the eﬀect of friction is to reduce the

eﬀective stress in Griﬃth’s criterion, thereby increasing

the critical crack size.

In the present form of SHALE, friction is not consid-

ered. Hence the frictional term in the above equation is

neglected. Also, SHALE-BCM only deals with the changes

in crack distribution due to the expansion of cracks and

the crack closing is not taken into account. Hence in com-

pression, the crack can grow if

≥

π T E

(1 − ν) c

(4)

The above theory applies to a single crack. The SHALE-

BCM model evolves a statistical ensemble of cracks, given

by the following initial distribution:

N = N

exp

−

c
¯c

(5)

where ‘N

’ is the total number density of cracks, ‘N’ the

number density of cracks having radius greater than ‘c’

and ‘¯c’ is a characteristic crack size.

At each timestep, the critical crack radius is computed

for each cell, using Eqs. (2) and (4). Cracks larger than

the critical size will grow during that timestep. The exact

V.R. Ikkurthi et al.: Simulation of rock fracture

147

000000

111111

000000

111111

5 cm

Steel

Water

= 1.0 gm/cc)

(

ROCK

(

= 2.0 gm/cc)

100 cm

50 cm

= 7.92 gm/cc)

Fig. 2. Schematic of computational region

distribution function of cracks is reconstructed at each

timestep.

The BCM evolves a dimensionless parameter γ, deﬁned

γ = N c

∞

(6)

γ plays the role of an internal variable. Physically, γ is

the cube of the ratio of crack size to crack spacing, and is

considered to be a macroscopic measure of the amount of

fracture. When γ is approximately equal to unity, cracks

are as big as they are far apart. This can be taken to

indicate fragmentation.

For the assumed crack orientation and distribution

the eﬀective components of the Compliance tensor C

yyyy

which is inverse of the elastic modulus tensor, is related

to the compliance of the matrix material C

yyyy

= C

yyyy

1 +

1 − ν

(7)

For more general crack distributions, ‘γ’ becomes a

two-index tensor. Now, the constitutive law takes the form

dε

ijkl

)

(8)

dσ

= C

−1

klij

dε

−

ijmn

(9)

3 Description of the method

Simulations have been done for Oil Shale, since necessary

data for the model is available in the literature (Demuth

et al. 1985). A schematic of the two-dimensional compu-

tational region is shown in Fig. 2. Since PBT is frequently

used for free-standing boulders weighing a few tonnes, we

have assumed a typical boulder size 1 m in diameter and 1

m long, having a density of 2 g/cc. This corresponds to a

mass of 1,570 kg. At the center of the rock, a 5 cm radius

bore hole is drilled to half the height, i.e., 0.5 m, and ﬁlled

with water. A steel tamper is placed on top of the hole,

which serves to conﬁne the high pressure inside the bore

during the discharge. We assume that energy is dumped

in the water over a time-scale that is short compared with

sound-transit times through the cavity. Numerically, this

is implemented by assigning a suitably high speciﬁc inter-

nal energy in a few computational cells located near the

bottom of the bore-hole.

We have used the BCM constitutive relation to model

the behaviour of oil shale (Demuth et al. 1985) and the

HOM Equation of State (EOS) for water and steel (Mader

1979) (see Appendix A).

We monitor the evolution of the system, particularly

the pressure-time waveforms at several points and the dis-

tribution of γ throughout the rock. We have performed

three studies, which are described in later sections.

In general, we assume that when γ exceeds unity in

a connected set of cells extending from the bore-hole to

a boundary, the rock would fracture along that surface.

From a practical point of view for PBT, the boulder should

be broken into 4–5 pieces of similar size. Keeping in mind

that our simulations are only two-dimensional, we have

deﬁned ‘rock fracturing’ as the condition when computa-

tional cells with γ >1 divide the rock into at least two

large-sized parts.

4 Results

Three diﬀerent issues have been examined in this study,

and are described below.

The four frames in Fig. 3 show, at diﬀerent times, the

distribution of cells in the rock that have γ >1. This is

done for a typical simulation with a uniform mesh having

40 and 80 cells in the radial and axial directions, respec-

tively. The cell width is 1.25 cm. At time t = 0, we start

with a speciﬁc internal energy of 0.337 × 10

12 erg

in one

cell located at the bottom of the water-ﬁlled cavity. All

the boundaries are treated as free boundaries.

We can see in Fig. 3 that γ exceeds unity in a con-

nected set of cells extending from the bore-hole to the

rock boundary in the radial direction. There is also con-

siderable cracking in the axial direction. Since the simu-

lation is two-dimensional (axisymmetric), this can be in-

terpreted to mean that the rock breaks into two or more

large-sized fragments of comparable size. Complete details

of the number of fragments and their sizes can only be

obtained from a 3-D simulation with a suitable fracture

model.

4.1 Pressure-time curve

Higher pressure levels in the water-ﬁlled cavity would

launch stronger pressure pulses, thereby speeding up the

process of crack growth. This means that a higher pres-

sure needs to be applied for a shorter time to produce

rock fracture. However, higher pressures come with their

own problems. For a given total electrical energy input,

148

V.R. Ikkurthi et al.: Simulation of rock fracture

Z (cm)

t = 0.1 msec)

0 25 50

R (cm)

100

Z (cm) Z (cm)

(t = 0.2 msec)

0 25 50

R (cm)

100

Z (cm) Z (cm) Z (cm)

(t = 0.5 msec)

0 25 50

R (cm)

100

Z (cm) Z (cm) Z (cm) Z (cm)

(t = 1.0 msec)

0 25 50

R (cm)

100

Fig. 3. Propagation of brittle crack in rock. Shading indicates

regions with γ >1

producing a higher pressure in the cavity requires depo-

sition of the energy into a smaller hot zone (Madhavan

et al. 2000). Now, energy deposition in water takes place

through the formation of an underwater arc between the

electrodes in the probe tip. Reducing the size of the hot

zone requires fabrication of a smaller probe tip, at the

same time retaining its ability to absorb high energies.

Furthermore, higher pressures in the hot zone also cause

greater damage to the probe tip. Of course, if the time

for which the higher pressure is applied is small enough,

the damage may be limited. It is, therefore, of interest to

determine the threshold P(t) curve for a given rock frac-

turing application.

With any practical arrangement for plasma blasting,

the pressure in the cavity ﬂuctuates in time, so there could

be an inﬁnite set of P(t) proﬁles that perform a given frac-

turing job. Our purpose is to gain physical insight, rather

than an exact determination of the optimum P(t). Hence

we simplify the problem by assuming that a constant pres-

sure ‘P’ is maintained in the cavity for a time τ, following

which the pressure instantaneously falls to zero. The simu-

lation must, of course, be continued well beyond τ to allow

time for the pressure waves to complete a few to-and-fro

transits between the bore-hole and the outer boundary.

We have assumed that all outer boundaries are held free.

The results of this study are shown in Fig. 4. As ex-

pected, the τ required to produce rock fracture decreases

with increasing P. However, the rate of decrease of τ be-

comes smaller at higher P, which means that there is little

incremental beneﬁt in generating very high pressures. It

would, in fact, be counter productive, since damage to the

blasting probe would be greater at higher pressures. At the

other extreme, τ increases rapidly as pressure decreases,

which shows that there is also a minimum pressure below

which cracking is not signiﬁcant.

This study shows that, from the point of view of mini-

mizing damage to the probe, there exists an optimal pres-

sure level for cracking a given rock-type and geometry.

4.2 P

peak

- E curve

Another quantity to be optimized is the volume of the

region in which the electrical arc deposits energy, i.e., the

‘Hot Zone’. For a given total energy E deposited in this

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

TIME (

sec)

PRESSURE (GPa)

Fig. 4. P-τ curve

200

250

300

350

400

450

500

550

600

650

0.5

1.5

2.5

3.5

ENERGY (kJ)

PEAK PRESSURE (GPa)

Fig. 5. P

peak

- E curve

region, a smaller volume of the Hot Zone implies a higher

peak initial pressure P

peak

in the vicinity of the probe tip.

The objective of the second problem is to determine the

curve in P

peak

− E space that just cracks the rock. The

boundary conditions for this study are the same as for the

P-τ study.

The results are shown in Fig. 5. We ﬁnd from the plot

that the minimum energy E required to produce rock frac-

ture decreases with increasing P

peak

. This means that the

system size and cost, which scales with the capacitor bank

energy, can be signiﬁcantly decreased by dumping the en-

ergy in a smaller region in water. On the other hand, be-

yond some P

peak

, there is only marginal reduction in E.

Higher P

peak

would also reduce the lifetime of the blast-

ing probe. Hence, there must be an optimal P

peak

that

lowers the energy consumption but with acceptable probe

damage. The optimal level obtained in our simulations is

E = 200 kJ for a 1,570 kg mass, i.e., ∼ 0.13 kJ/kg. This

matches well with the numbers reported in (Nantel and

Kitzinger 1990) for hard rock.

V.R. Ikkurthi et al.: Simulation of rock fracture

149

-1.2e+09

-1e+09

-8e+08

-6e+08

-4e+08

-2e+08

2e+08

4e+08

6e+08

2e-05

4e-05

6e-05

8e-05

0.0001

0.5

1.5

2.5

3.5

(dyne/cm

)

TIME (sec.)

(a)

-2e+08

-1.5e+08

-1e+08

-5e+07

5e+07

1e+08

1.5e+08

2e+08

5e-05

0.0001 0.00015 0.0002 0.00025 0.0003

0.05

0.1

0.15

0.2

0.25

0.3

0.35

0.4

0.45

0.5

(dyne/cm

)

TIME (sec.)

(b)

Fig. 6. Variation of σ

, σ

and γ in a a cell near the cavity,

b a cell at the boundary

4.3 Mechanism of failure

The pressure pulses or shocks originating in the Hot Zone

travel into the rock, producing compression. When they

reach the free boundaries of the boulder, they give rise to

tensile waves which travel inward.

The simulation yields the temporal evolution of the

crack parameter γ throughout the rock. The evolution of

γ in a computational cell depends upon the applied stress

level and the duration of its application. We have studied

the mechanism by which the rock fails, i.e., whether the

cracks are growing in tension or in compression. If the

material fails in compression, it is due to the shear stress

(σ

), while if it is in tension, it is due to a combination

of normal stress (σ

) and shear stress (σ

It is observed that regions near the boundaries fail

mainly in tension (σ

> 0), while regions near the water-

ﬁlled cavity normally fail in compression (σ

< 0). This

can be seen from Fig. 6, which shows crack evolution in

cells near the boundary and the bore-hole.

00000

11111

00000

11111

Steel

Water

ROCK

(FREE or FIXED BOUNDARY)

ROCK BASE

AXIS

FREE BOUNDARY

Fig. 7. Schematic of computational mesh showing boundary

conditions

4.4 Asymmetry in pressure levels between axial

and radial directions

It has been observed experimentally that a boulder with

a ﬂat base, resting on the ground, is more diﬃcult to frac-

ture as compared to one where the ground supports only

a small part of the base. This can be physically explained

as follows. When the shock breaks out of the unsupported

base (‘free boundary’), a rarefaction wave moves in. The

material goes into tension, leading to the growth of cracks.

If the base is well-supported on the ground, it approxi-

mately behaves like a ﬁxed boundary, so a rarefaction is

not produced. In the latter case, the cracks would have

to grow only under compression, which appears to require

higher energy input.

We have attempted to reproduce this eﬀect in sim-

ulations. This has been done by using a ﬁxed-boundary

condition instead of a free-boundary for one ﬂat surface,

i.e., the ‘base’ of the boulder, as shown in Fig. 7. How-

ever, we ﬁnd only a small change in the energy required

for fracture. This surprising result arises from a numerical

artiﬁce. Since this artiﬁce could also aﬀect other results in

2-D axisymmetric simulations, it is suﬃciently important

to merit a detailed discussion, which apppears below.

If the hot zone were idealized as a point, and in the ab-

sence of fracture, we would expect it to drive a spherically

diverging shock/pressure wave. Such a wave should pro-

duce identical pressure levels at a given distance from the

hot zone. In our simulations, the hot zone is cylindrical in

shape due to the mesh structure in SHALE. This seems

to produce generally higher pressure levels in the axial di-

rection in comparison with the radial direction. Due to

this asymmetry, regions ‘below’ the water-ﬁlled cavity fail

under compression during the ﬁrst passage of the shock.

Cells at the same axial location as the hot zone, but dis-

placed radially, experience a comparatively low pressure.

Hence these cells fail only after 2–3 back-and-forth tran-

sits of the shock. Since cells below the hot zone generally

fail during the ﬁrst transit, the boundary condition at the

bottom does not play a signiﬁcant role.

150

V.R. Ikkurthi et al.: Simulation of rock fracture

50 cm

100 cm

OTZONE

( 1 )

( 3 )

Fig. 8. Schematic showing marker cells (not to scale)

This asymmetry in pressure distribution can be due

to three eﬀects. The ﬁrst is the cylindrical shape of the

hot zone. The second is the presence of the water-ﬁlled

cavity. The third is that in the BCM model implemented

in SHALE, crack planes are assumed to be oriented normal

to the axial direction. All these factors contribute towards

destroying spherical symmetry. Of these, the second eﬀect

has a physical basis, while the ﬁrst and third eﬀects are

of purely numerical origin. The ﬁrst eﬀect is curious and

is therefore examined below.

To study the ﬁrst eﬀect, we monitored the pressure

histories, P(t), in the computational cells denoted 1 and

3 in Fig. 8. The cells are equidistant from the hot zone.

In order to eliminate the second and third eﬀects, two

runs have been performed, with the entire cylindrical do-

main ﬁlled with water. All boundaries are taken to be

ﬁxed. The initial density in all cells is taken equal to 1.001

g/cc, which is slightly greater than the normal density

of water. This results in an initial pressure ∼ 2.2 × 10

dyne/cm

in the marker cells.

Figure 9 shows the pressure history in marker cells 1

and 3. Considerable asymmetry can be observed. Clearly,

this asymmetry is caused by the cylindrical shape of the

hot zone.

Now, we expect that a pressure pulse starting from a

cylindrical hot zone should become progressively spheri-

cal as it moves away from the source. To test this idea, we

performed another simulation with a non-uniform mesh,

with the mesh size being a minimum in the hot zone, and

increasing in all directions as we move further away. This

means that, for the same cell numbers of the marker cells,

the size of the hot zone becomes smaller in relation to

its distance from the marker cell. In this run, the ratio

is chosen to be 0.01. Again, the simulation was run with

water ﬁlling the entire domain. The results are shown in

Fig. 10. In the above simulations with uniform and vari-

able mesh in water, the same amount of total input energy

2e+07

4e+07

6e+07

8e+07

1e+08

1.2e+08

1.4e+08

1.6e+08

0.0001

0.0002

0.0003

0.0004

Pressure (dyne/cm

)

Time (sec.)

(1)

(3)

Fig. 9. Asymmetry in axial and radial pressure levels in water

withuniform mesh

2e+07

4e+07

6e+07

8e+07

1e+08

1.2e+08

1.4e+08

0.0001

0.0002

0.0003

0.0004

Pressure (dyne/cm

)

Time (sec.)

(1)

(3)

Fig. 10. Axial and radial pressures in water witha variable

mesh

equal to 2.04 × 10

ergs is dumped in the hot zone, but

the volume of the hot zone in the two cases is 6.136 cm

and 0.201 cm

respectively. The higher speciﬁc internal

energy of the hot zone in the variable mesh case leads to a

higher initial pressure of 1.015 × 10

dyne/cm

, as com-

pared to 3.33 × 10

dyne/cm

in the uniform mesh case.

This diﬀerence in initial pressure levels in the hot zone im-

plies some diﬀerences of detail between the P(t) curves in

Figs. 9 and 10. This is indeed seen in the ﬁgures. However,

since the energy input is the same, the peak pressure at

the observation points is almost the same in both cases,

since the observation points are far from the hot zone.

The plots show better agreement than in Fig. 9, sup-

porting our claim that the cylindrical hot zone is to blame

for asymmetry.

This conclusion implies that using a spherical hot zone

would further improve spherical symmetry. To test this,

we have set up a hot zone which is approximately spher-

ical, and with a variable mesh as before. The volume of

the spherical hot zone is taken as 4.7 cm

. The speciﬁc

internal energy, for the same total input energy as in the

previous two cases, is 4.34×10

erg/g, which corresponds

to an initial pressure of 4.34 × 10

dyne/cm

in the hot-

V.R. Ikkurthi et al.: Simulation of rock fracture

151

2e+07

4e+07

6e+07

8e+07

1e+08

1.2e+08

0.0001

0.0002

0.0003

0.0004

Pressure (dyne/cm

)

Time (sec.)

(1)

(3)

Fig. 11. Axial and radial pressures in water when hot-zone is

spherical and with variable mesh

zone. Figure 11 shows the pressure levels. The plot shows

better agreement at the ﬁrst peak, as compared to Figs. 9

and 10.

As stated earlier, apart from the shape of the hot zone,

there are two other eﬀects which can cause asymmetry. A

study of these eﬀects is beyond the scope of this work.

In an experimental situation, the shape of the hot zone

will be determined by the design of the blasting probe

and the physics of underwater arc formation. These will

vary from problem to problem. Hence asymmetries in frac-

turing may arise due to physical reasons. The purpose of

this section is to highlight numerical problems which could

produce spurious asymmetry.

5 Limitations of the model

1. The version of the Bedded Crack Model (BCM) in-

corporated in the published version of SHALE does

not deal with detailed physics of crack nucleation and

growth. For example, it does not handle the growth

of cracks with random orientations. A more detailed

brittle fracture model must be incorporated, such as

the one described in (Curran and Seaman 1987).

2. A high-accuracy EOS ﬁt is available for water (Saul

and Wagner 1989). A comparison of this ﬁt with val-

ues from the HOM EOS shows signiﬁcant deviations,

even at low pressures. This EOS has not been used in

the present work since it is computationally too ex-

pensive to incorporate into a 2-D code. However, the

mismatch means that the pressure proﬁles, and frac-

ture performance, may change to some extent.

3. Experiments sometimes use a suspension of bentonite

in water (Nantel and Kitzinger 1990) in order to in-

crease its viscosity. EOS data must be generated for

such suspensions.

4. With any fracture model, the reliability of the simu-

lation depends upon the model parameters input to

the code. In the present work, these have been taken

from (Demuth et al. 1985). In a more detailed study,

these parameters should be obtained after matching

with experimentally-observable properties for the rock

material of interest. This matching should be done for

tensile as well as compressive properties.

5. One detail of interest to the blasting community is to

determine the precise pattern of rock fracture, such

as the number of large-sized fragments. This would

require a three-dimensional simulation, along with a

suitable fracture model.

6 Conclusions and future work

We have examined the eﬀectiveness of diﬀerent power de-

position proﬁles in rock fracturing, using two-dimensional

hydrodynamic simulations. There are three major conclu-

sions.

Firstly, assuming that a constant pressure P is main-

tained in the cavity for a time τ, we ﬁnd that there exists

an optimal pressure level for cracking a given rock-type

and geometry, such that damage to the blasting probe is

minimized. Secondly, we have varied the volume of water

in which the initial energy E is deposited, which corre-

sponds to diﬀerent initial peak pressures P

peak

. It is found

that there must be an optimal P

peak

that reduces the en-

ergy consumption, but with acceptable probe damage. For

comparable damage, a high stress for very short duration

takes less energy than longer pulses of lower amplitude.

This result contrasts with the usual assumption that dam-

age or other eﬀects are proportional to the total energy.

Thirdly, we ﬁnd that portions of the rock near the water-

ﬁlled cavity undergo fracture in compression due to shear

stress, while portions near the boundary fracture in ten-

sion due to a combination of normal and shear stresses.

We have also highlighted certain numerical errors that

can arise in these simulations.

In future, we plan to incorporate a more detailed nu-

cleation and growth model for cracks into our code, along

with more accurate equation-of-state data for materials

found in these experiments.

Appendix

A Equation of State

The HOM EOS has been used for water and steel -details

are available from (Mader 1979). The major features are

given below for the sake of completeness.

A.1 Material in compression

Here we have V

> V

, where V

and V

are the speciﬁc

volumes of the condensed material in the initial and cur-

rent states.

The experimental principal hugoniot is given by

= C + S Up

(10)

152

V.R. Ikkurthi et al.: Simulation of rock fracture

where U

and U

are the shock and particle velocities re-

spectively, and C and S are best-ﬁt coeﬃcients, speciﬁc to

the material.

The Rankine-Hugoniot relations then yield the hugo-

niot pressure as

− V

)

− S (V

− V

)]

(11)

The temperature on the hugoniot is given by

= F

lnV

(lnV

)

(lnV

)

(lnV

)

(12)

where F

, G

, H

, I

and J

are best-ﬁt coeﬃcients,

and the speciﬁc internal energy on the hugoniot by

− V

)

Mbar-cc

(13)

The required pressure P and temperature T are obtained

by adding corrections to the corresponding hugoniot val-

ues at the current density:

P = P

− I

)

(Mbar)

(14)

T = T

− I

)

(15)

where γ

is the Gruneisen coeﬃcient, I

the current spe-

ciﬁc internal energy and C

the speciﬁc heat capacity.

A.2 Material in rarefaction

For V

< V

, we have

P = P

−

3α

− 1

(16)

T = T

(17)

where α is the linear expansion coeﬃcient of thermal ex-

pansion.

Acknowledgements. We are grateful to one of the reviewers,

whose comments have helped in improving the quality of the

paper.

References

Adams TH, DemuthRB, Margolin LG, Nichols BD (1983) Sim-

ulation of Rock blasting withthe SHALE code, Report

LA-UR-83-399, Los Alamos National Laboratory

Curran DR, Seaman L (1987) Dynamic Failure of Solids,

Physics Reports, 147(5 and 6):253–388

DemuthRB, Margolin LG, Nichols BD, Adams TF, SmithBW

(1985) SHALE: A Computer Program for Solid Dynamics,

Report LA-10236, Los Alamos National Laboratory

Dienes JK (1979) Explosively Produced Fracture of Oil Shale,

Report LA-8104-PR, Los Alamos Scientiﬁc Laboratory An-

nual Report

Griﬃth AA (1920) The Phenomena of Rupture and Flow in

Solids, Philosophical Transactions of the Royal Society A,

221:137–154

Mader C (1979) Numerical Modeling of Detonations, Univer-

sity of California Press, Berkeley

Madhavan S, Doiphode P, Chaturvedi S (2000) Modeling of

Shock wave generation in water by electrical discharges.

IEEE Transactions on Plasma Science: Special issue on

Pulse Power, 28(5):1552–1557

Margolin LG, SmithBW (1985) A numerical model for simu-

lating dynamic processes in rocks, Report LA-UR-84-2947,

Los Alamos National Laboratory

Nantel JH, Kitzinger F (1990) Plasma Blasting Techniques.

Proc. 3rd International Symposium on Rock Fragmenta-

tion by Blasting, pp. 79–82

Saul A, Wagner W (1989) A Fundamental Equation for Wa-

ter covering the range from the melting line to 1273 K

at pressures upto 25000 MPa, J. Phys. Chem. Ref. Data,

18(4):1537–1562

Wyszukiwarka

Podobne podstrony:
Farina, A Pyramid Tracing vs Ray Tracing for the simulation of sound propagation in large rooms
0 Simulation of fatigue failure in a full composite wind turbine blade Shokrieh Rafiee 2006
(wydrukowane)simulation of pollutants migr in porous media
The Effect of DNS Delays on Worm Propagation in an IPv6 Internet
Simulation Of Heavy Metals Migration In Peat Deposits
Being Warren Buffett [A Classroom Simulation of Risk And Wealth When Investing In The Stock Market]
A Survey of Irreducible Complexity in Computer Simulations
Biological Models of Security for Virus Propagation in Computer Networks
Kang, Meng Sound Propagation in Micro Scale Urban Areas Simulation and Animation
Dance, Shield Modelling of sound ®elds in enclosed spaces with absorbent room surfaces

więcej podobnych podstron