Frederick Bloom Department of Mathematical Sciences,
Northern Illinois University,
DeKalb,IL60115Wave Propagation in Partially Saturated Soils1
1Introduction
In this article we will survey some of the methodology that has been developed in order to study the problem of wave propagation in soils or sediments;speci?cally,we are interested in the behav-ior of explosion induced shock waves propagating in partially saturated soils or sediments.Because waves in soils may be in-duced by many type of sources?e.g.,earthquakes,machinery,ex-plosions?the topic at hand has a long history and occupies a central position in the broad subject of soil dynamics;these dif-ferent sources of propagating waves in soils must often be treated in different ways.In particular for determining the plastic defor-mations of a soil,or partially saturated sediment,which are caused by an explosive detonation,it is not possible to use the rheological model employed in investigating the processes of soil compaction caused by impact or vibration;this is well understood in the lit-erature as being explained by the fact that the dynamic pressures and the velocities of deformation associated with an explosive detonation are many times greater than during an impact or vibra-tion.Thus,the rheological features of the sediment in the zone of plastic deformations are different with an explosion as opposed to a small dynamic disturbance.In addition,it was established,some time ago,that there is an appreciable difference,in terms of the response to loadings generated by explosive detonations,between soils saturated and not saturated with water;in the latter type of soil,or sediment,plastic deformations predominate under loading while the elastic deformations are small and are observed only when the load is removed.In saturated sediments,however,it is mainly elastic deformations which are originated even under very high pressures.These conclusions form the basis for the use of a “plastic gas”medium in modeling the behavior of a partially satu-rated sediment near the source of an explosive detonation in the sediment;such models,which were studied quite early in the Rus-sian literature?e.g.,Refs.?1–5??,have been brought into a satis-factory analytical state?Refs.?6–9??with respect to studying the dynamic compacting of a partially saturated soil or sediment in the neighborhood of a concentrated explosive source.The“plastic gas”models refer,in particular,to a medium which,upon being loaded,changes its density according to a de?nite law,but is such that when the load is removed?i.e.,when there is a marked re-duction in pressure?it retains the density obtained upon loading.If one has an interest primarily in the compaction of a partially satu-rated sediment in the immediate neighborhood of an explosive detonation,which is caused by the generation of a shock wave induced by that detonation,be it either spherical,cylindrical,or plane,then the appropriate approach would appear to be the one discussed in Sec.2.3,which is based on the use of a plastic gas continuum model.
An important point in dealing with the modeling of a partially saturated soil?or sediment?,and one that was explicitly recog-nized some time ago?e.g.,see Nelson,Baron,and Sandier?10?, Markle and Dass?11?,and Mandel?12??,is that many of the wave propagation problems for soils or sediments,which are associated with an explosive detonation and,therefore,concern themselves with waves generated by contained blasts,require a mathematical de?nition of the material behavior of the sediment over a wide range of pressures unless,of course,one is primarily interested in the compaction process near the source of the explosion.The observed material behavior in partially saturated soils or sedi-ments may range from the essentially hydrodynamic?uid behav-ior,at extremely high pressures?megabars?,in a neighborhood of the blast?the“plastic gas”zone?,to inelastic solid behavior in regions which are at intermediate pressures?of tens of kilobars or less?,to regions which behave,essentially,as an elastic medium at suf?ciently low pressures.
Our main focus in this article is on the shock wave induced compaction of a partially saturated soil,or sediment,which takes place in the immediate vicinity of an explosive detonation in such a sediment.Thus,in Sec.2we treat the problem of shock propa-gation in the plastic gas zone,beginning with the equations which govern the dynamic compacting of a sediment around a spherical source of explosion;an approximate solution method and numeri-cal formulation of this problem is also presented here along with some numerical results.Our work in Sec.2then continues by formulating the detonation problem in the plastic gas zone in space dimension n=1,2,3so as to be able to deal not only with spherical shock waves,but also with cylindrical shock waves and plane shock waves.Section2concludes by indicating how one may transition,mathematically,from the high pressure plastic gas zone to the moderate pressure elastic-plastic zone.
Once the induced wave?say,spherical shock-wave?front tran-sitions out of the plastic gas zone and into the elastic-plastic re-gime it is no longer a shock but is,rather,a propagating spherical wave;in Sec.3we review what the literature has to say about this aspect of the problem,focusing on models which correctly incor-porate dilatancy as a kinematical constraint.Section4takes the wave front into the low pressure elastic regime of the sediment and elucidates the well-known two-wave structure that exists in that zone.
2Shock Propagation in the Plastic Gas Zone
2.1Dynamic Compacting of Soil Around a Spherical Source of Explosion.In this section we analyze the spherical compacting of a partially saturated soil around the spherical source of an underground explosion.A slightly more general analysis,which applies to plane and cylindrical sources as well, will be presented in Sec.2.2.It is assumed in this section that at a high pressure of around107–109N/m2a soil undergoing loading behaves like an ideal multi-component?uid,while during unload-ing it is incompressible.The pressure-density graph which holds in this situation represents the behavior of what has been termed a plastic gas?1–5?.The basic equations which we present below are derived for a three-component soil and are valid for large defor-mations.While we limit ourselves in this section to an analysis of the solution of the large pressure shock wave problem,in Sec.2.3 we will consider the transition from the high pressure plastic gas zone to the moderate pressure elastic-plastic zone.
Large pressure caused by an explosion in a partially saturated soil is accompanied by plastic volumetric deformations and any suitable mathematical description of such an event must take into account the nonlinear physical properties of the soil as well as geometrical nonlinearities.The basic equations are derived below by employing Lyakhov’s constitutive equations for a three com-
1Based on a report prepared for the Mechanics of Materials Branch,CODE6386,
United States Naval Research Laboratory,Washington,DC20375,on Co-Op Job
Order63-5596-E4,Shock Wave Propagation in Gassy Sediments,1996.
Transmitted by Assoc.Editor S.Adali.
Applied Mechanics Reviews JULY2006,Vol.59/177
Copyright?2006by ASME
ponent soil consisting of sand grains?quartz?,water,and gas trapped in the pore space between the grains.As has been noted in Ref.?6?,in unsaturated soils under conditions of low and medium pressures?normal stresses?it is important to take into account shearing stresses so that both the deformation and stress?elds have a tensor character.But,in the range of very high pressures, shearing stresses have no signi?cant meaning in the process of the dynamic compacting of the soil.Thus,it may be assumed that under such high pressure conditions,during the loading process, the soil can?approximately?be treated as an ideal?uid;for a general discussion of the explosive compacting of soils the reader is referred to Refs.?6–8,33,34?.
2.1.1Basic Equations and Problem Formulation.We begin by considering a spherical cavity in a semi-in?nite half-space and suppose that a large pressure?eld is suddenly imposed at the boundary of the cavity?at time t=0?which then decreases mono-tonically.We introduce a spherical coordinate system X1=R,X2 =?,X3=?for the underformed state of the half-space in a neigh-borhood of the cavity while in the deformed state x1=r,x2=?, x3=?represents the position of a particle at time t?0.The mo-tion of a particle near the cavity will then be described by rela-tions of the form
r=r?R,t?,?=?,?=??2.1?
The deformation gradient associated with the motion in Eq.?2.1?may be decomposed as the composition of a rotation and a sym-metric deformation,i.e.,F=RU,and the corresponding Cauchy-Green tensor C=FF t?U2,or
C KL=g ij ?x i
?X k
?x j
?X L
?2.2?
where for the case of spherical symmetry the components of the metric tensor are given by
g11=1g22=r2g33=r2sin2??2.3?in the deformed con?guration,while in the undeformed con?gu-ration the metric tensor G ij has components
G11=1G22=R2G33=R2sin2??2.4?The nonzero components of the Cauchy-Green tensor C are then C11=??r?R
?2??r,R?2C22=r2C33=r2sin2??2.5?while the mixed components C L M?G KM C KL are
C11=?r,R?2C22=C33=r2
R2
?2.6?
The ratio of the initial density?0in the soil to the density?at the current time t is given by
?0?=?det C=?det?G KM C
KL
?=r
2
R2
r,R?r,R?0??2.7?
We assume that starting from time t=t*,the motion of the soil is incompressible so that,for t?t*,??R,t?=?u?R?.We now inte-grate Eq.?2.7?so as to obtain for some??t?:
r3?R,t?=3?R0R?0?u????2d?+?3?t?,t?t*?2.8?
and we note that r3?R0,t?=r03?t?,where R0is the initial radius of the spherical cavity.Thus,for the case of an incompressible mo-tion,it follows from Eq.?2.8?that
r3?R,t?=?3?t?+?3?R??2.9?We also compute,using Eq.?2.8?,that in the deformed con?gu-ration?i.e.,using the Euler description of the motion?the velocity ?eld has components
v r=
?r
?t=
1
3r2
d
dt
?3?t?
?2.10?
v?=v?=0
and,as a consequence,satis?es the incompressibility condition
??˙
?=
?v r
?r+2
v r
r
=0?2.11?Because the acceleration in an incompressible motion is given by
a r=
?v r
?t+v r
?v r
?r a?=a?=0?2.12?we have
a r=
1
3r2
d2
dt2
?3?t??2
9r5
d
dt
??3?t??2?2.13?
In this section we focus on the problem of propagation of a spheri-cal shock wave in the plastic gas zone around the source of an explosion;in subsequent sections we will consider,also,plane and cylindrical shock waves as well as the transition from the high pressure plastic gas zone to the moderate pressure elastic-plastic zone.We begin by elaborating the constitutive equation?i.e., equation of state?for the plastic gas.As has appeared very often ?especially in the Russian literature?in the high pressure range ?i.e.,107–109N/m2?there appears to be a consensus opinion that shearing stresses are almost insigni?cant in the process of dy-namic?global?compacting of a partially saturated soil,so that under these conditions,during the loading process,the soil tends to?approximately?behave like an ideal gas.For saturated soils, the equations of state?constitutive laws?appear to have been for-mulated?rst by Lyakhov?1?in1959;the basic ideas are described below.
We consider a three component?soilgrains,gas,water?volume element of soil of total mass M and we subscript masses and volumes of the individual components by k,where k=1corre-sponds to the soil matrix?sand grains?,k=2to the water in the pore space,and k=3to the gas in the pore space.We let the volume of the components at time t=0be denoted by V?k?0and the total volume of the element by V0;the corresponding volumes at any arbitrary time t?0are denoted by V?k??t?and V?t?,respec-tively.Thus
V?t?=?k=1
3
V?k??t?
?2.14?
V0=?k=1
3
V?k?0
Denoting by??k??t?and??t?the average density,respectively,of the individual components and the volume element at time t?0, and by??k?0and?0the corresponding values at t=0,we have
??k??t?=m k/V?k??t?,??t?=M/V?t?
?2.15?
??k?0=m k/V?k?0,?0=M/V0
where m k denotes the mass of the k th component in the volume element;it is assumed that both m k,k=1,2,3and M remain con-stant throughout the entire deformation process so that the mass
178/Vol.59,JULY2006Transactions of the ASME
fractions are also constant.If we now denote the volume fractions
at time t by v ?k ??t ?,and the initial volume fractions by v ?k ?0
?v ?k ??0?,k =1,2,3,then as the mass fractions of the elements are not changed by the deformation which occurs during an explosion we obtain
m k M =??k ??t ?V ?k ??t ???k ?V ?k ?
=??k ?0V ?k ?
?0V 0?2.16?
in which case
?0??t ?v ?k ??t ?=??k ?
??k ??t ?v ?k ?
0?2.17?
However,?k =13
v ?k ??t ?=1,so we have
?0
??t ?
=?k =1
3
??k ?
0?
?k ??t ?
v ?k ?
?2.18?
The relationship in Lyakhov ?1?which is adopted for constitu-tive behavior of the individual components ?see also Ref.?3??has
the form
p ?t ?
p 0=1+1?k
????k ?
??k ?
0??k
?1
?
?2.19?
1?k
=??k ?0c k 2?k p 0
where p ?t ?is the pressure at time t ?0,p 0=p ?0?,the ?k ,k =1,2,3,are the polytropic curve coef?cients of the individual components,and the c k ,k =1,2,3,are the speeds of sound propa-gation in the individual components ?which are computed for small disturbances at pressure p =p 0?.For the gas,which is trapped in the pore space of the sediment,we have
c 32=?3
p 0
??3?
T?3=1
?2.20?
Employing Eq.?2.18?in Eq.?2.19?we now obtain for the global constitutive relationship connecting p ?t ?and ??t ?
?0
??t ?
=?
k =1
3
??k
?
p ?t ?
p 0
?1?+1?
?1/?k
v ?k ?
0?2.21?
where the value of ?0may be computed from
?0=
?k =1
3
?
?k ?0v ?k ?
0?2.22?
The constitutive relation ?2.21?forms the basis for our discussion
of the compaction of a partially saturated soil ?sediment ?under the conditions of very high pressures which are present in the neigh-borhood of a spherical source of explosion in the soil.The param-eter values in Eq.?2.21?may be identi?ed by using the table of ?typical ?values for the components which is given below;this table has been reproduced from the similar table in Ref.?6?.In the model proposed by Lyakhov ?1?for the behavior of a partially saturated soil,it is assumed that Eq.?2.21?describes the loading process which occurs in an explosion;it is,however,also assumed that there exists a lower bound for the pressure,p =p L ,above which that constitutive relation may be used.If we go be-low p =p L then the compressibility of the sediment will depend,to a large extent,on the compressibility of the soil skeleton and the mechanism which governs the deformation will differ from that which is described by Eq.?2.21?.In the Lyakhov model ?1?,for p ?p L ,an experimentally derived function of the form p =p ???is used in place of Eq.?2.21?.We note,also,that an analysis of the present model of a plastic gas in a neighborhood of the source of
the explosion shows that the resulting relation p =p ???is very sensitive to changes in the parameter v 3:even a small gas content in the sediment yields considerable changes in the value of p for a given ?and,thus,has a major in?uence on the way in which a propagating shock front decays.In the present situation,dynamic loadings with pressures of the order of 5?109N/m 2can occur and Eq.?2.21?will be assumed to be applicable.
In Fig.1we depict the graph of the pressure p as a function of 1??0/?.As indicated,the loading process takes place along the curve OA which is given by ?a u subscript denotes the value of the indicated variable at the shock wave front ?
??p ?=1?
?0
?u
??ˉ???=p ?2.23?
As the pressure changes,abruptly,from p 0to p u ,so does the
density which goes from the value ?0to the value ?u =?u ?r ,t ?,which is the density of the soil at the shock wave front;however,the unloading process takes place along the line AB which corre-sponds to constant density:1??0/?u =const.The equations of state for the partially saturated soil ?sediment ?are,therefore,as follows:
?0
?2
p ˙=
???
?p p ˙if p ??ˉ???=0,p ˙?00if p ??ˉ
????0or p =?ˉ?p ?and p ˙?0?2.24?
Now,the equation of motion for the partially saturated soil,in Eulerian coordinates,is just
?a r =?
?p ?r
p =p ?r ,t ??2.25?
and Eq.?2.25?is subject to the initial conditions
p =p 0
?=?0v r =0
r =R at t =0?2.26?
as well as the boundary condition
p ?r ,t ??r =r 0?t ?=p 0?t ?
?2.27?
where r =r 0?t ?,r 0?0?=R 0,describes the motion of the spherical cavity;in fact r ?R ,t ??R =R 0=r 0?t ?.In Fig.2we depict the solution of the above initial-boundary value problem:a form for the pres-sure change p =p 0?t ?,at the boundary of the spherical void con-taining the explosive charge,is assumed in Eq.?2.27?,and the problem is formulated for a shock wave propagating into an un-disturbed region,D 0,ahead of the wave.In this undisturbed
re-Fig.1Pressure p variation as a function of 1??0/?…adopted
from Ref.?6?…
Applied Mechanics Reviews
JULY 2006,Vol.59/179
gion,D 0,we have ?=?0,p =p 0,and v r =0.Because radial dis-placements at the shock wave front are equal to zero,at the wave front we have r =R .If we let w ?t ?be the function which describes the motion of the shock wave front,i.e.,
r =R =w ?t ??t =w ?1?R ?
?2.28?
then at time t =0we have
w ?t =0=R 0?w ?1?R 0?=0
?2.29?The pressure,density,and velocity of material particles at the shock front are denoted by
?u =??r =w ?t ??,p u ?t ?=p ?r =w ?t ??
?2.30?
v u ?t ?=v r ?r =w ?t ???r 0
2w 2r 0
while
F u ?t ?=F ?r =w ?t ??
?2.31?
is the value of the deformation gradient F ?R ,t ?=?r ?R ,t ?/?R at the
shock wave front.Because r =R at a shock front,use of Eq.?2.7?yields
?0/?u ?t ?=F u ?t ??2.32?
This last result can also be derived from the principle of mass
conservation,which yields the following jump conditions at the wave front:
?F ?=??0/???F ?=F u ?1?2.33?
??0/??=?0/?u ?1
with,e.g.,?F ?=F ?r =w ?t ?+?F ?r =w ?t ??
We now seek to derive the kinematic continuity condition at the shock wave front;to do this we denote the shock wave propaga-tion velocity by ??t ?=dw /dt and we note that this is the propa-gation velocity in the Lagrangian and Euler descriptions because the medium into which the shock is propagating is undisturbed.At the shock wave front we have
r ?w ?t ?,t ?=R
?2.34?Using the de?nition of F ?R ,t ?,and the jump conditions in Eq.
?2.33?,we obtain
??t ?=v u ?t ?+F u ?t ???t ?
?2.35?from which it follows,by virtue of Eq.?2.32?,that
?u ?t ????t ??v u ?t ??=?0??t ??2.36?For the dynamic continuity condition ?e.g.,Mandel ?12??we have p u ?t ??p 0??u ?t ????t ??v u ?t ??v u ?t ?
?2.37?
and,therefore,Eq.?2.36?,with Eq.?2.37?,and the constitutive law ?2.24?,form a system of three equations for the four unknown quantities at the shock wave front:p u ,?u ,v u ,and ?.Using Eq.?2.36?,the relation ?2.37?may be simpli?ed so as to read
p u ?t ??p 0=?0?2?t ??1??0?u ?t ??
?2.38?
2.1.2The Equations in the Incompressible Zone D .We return to Eq.?2.8?and de?ne the function ??t ?by taking into account the fact that for r =r 0?t ?we have R =R 0;thus ??t ?=r 0?t ?and
r 3
?R ,t ?=3
?
R 0
R
?0?u ???
?2
d ?+r 03?t ??2.39?
Equation ?2.39?describes the general form of the motion of the
sediment in the incompressible zone after the explosion at time t =0;at this junction the density and the function r 0,which de-scribes the motion of the spherical cavity,are not explicitly known.Since Eq.?2.39?is valid at any point of the incompress-ible zone D it is valid at the shock wave front where v =R =w ?t ?;therefore,Eq.?2.39?implies that
w 3
?t ??3
?
R 0
w ?t ?
?0?u ???
?2
d ?=r 03?t ??2.40?
which,in turn,provides an additional equation for r 0?t ?.Differen-tiating Eq.?2.40?through with respect to t we ?nd that
r 02r
˙0=w 2??
1??0?u ?t ?
?
?2.41?
or,equivalently,
d dt ?r 03?t ??=?1??0?u ?t ??
d dt
?w 3?t ???2.42?
Equation ?2.42?has resulted from the assumption that the region behind the shock wave front is incompressible;it may also be obtained from the incompressibility condition,as follows:at time t ?0the volume of the incompressible zone D is ?see Fig.3?
V ?t ?=4
3
??w 3?t ??r 03?t ???2.43?
However,as a consequence of their increase in mass of zone D ?due to compaction ?
dV dm =
1?u ?t ?
?with dm =?0·4?w 2dw ?
?2.44?
so
dV dt =4?w 2
·?0?u ?t ?dw dt ?4?3·?0?u ?t ?d ?w 3?dt
?2.45?
However,in view of Eq.?2.43?
dV dt =4
3??d dt ?w 3??d dt ?r 0
3
??
?2.46?
and Eq.?2.42?now follows by combining Eqs.?2.45?and ?2.46?.
Equation ?2.36?and ?2.37?,the constitutive law p u =?ˉ
??u ?,and Eq.?2.41?form a system of four equations for the ?ve unknown functions p u ,?u ,v u ,w ?t ?,and r 0?t ?.If we employ Eq.?2.25?
,
Fig.2…r ,t …-phase plane …adopted from Ref.?6?…
180/Vol.59,JULY 2006
Transactions of the ASME
recall that the acceleration a r is given by Eq.?2.12?,coupled with Eq.?2.10?,and use,in Eq.?2.10?,the fact that??t?=r0?t?,then the equation of motion?2.25?assumes the form
??r,t??13r2d2dt2?r03?t???29r5?d dt?r03?t???2?=??p?r,t??t?2.47?
From?2.39?we may conclude that the density function?which appears in Eq.?2.47?is a function of the variable z=r3?r03?t?,in the zone D,where the incompressible motion occurs,i.e.,
??r,t?=?u?z?in D?2.48?We are now in a position to formulate the complete initial-boundary value problem which governs motion in the incompress-ible zone D between the expanding spherical cavity and the ex-panding shock wave front;our task is to?nd the functions?=??r,t?,p=p?r,t?in the domain
D=?r?r0?t??r?w?t??,?2.49?
as well as the functions r=r0?t?,r=w?t?which de?ne D.In D the functions?and p must satisfy the partial differential equations ??r,t??13r2d2dt2?r03?t???29r5?d dt?r03?t???2?=??p?r,t??r?2.50a?
?˙=0????r,t?
?t+
???r,t?
?r·
1
3r2
d
dt
?r03?t??=0?2.50b?
which follow from the equation of motion in D,and the incom-pressibility condition,while the functions r0and w satisfy an or-dinary differential equation,which relates the velocity of the boundary of the spherical cavity and the velocity of the shock wave,namely
r02?t?r˙0?t?=w2?t?w˙?t?F?w˙?t??
F?w˙?t??=1??0/?u?t??2.51a?where?see Eq.?2.38??
w˙2?t?=??ˉ??u?t???p0?/?0?1??0/?u?t???2.51b?
Along the curves r=r0?t?and r=w?t?,the functions??r,t?and p?r,t?satisfy the boundary conditions
p?r,t??r=r
?t?=p0?t?
p u?t??p0=?0?w˙?t??2?1??0/?u?t???2.52?
p u?t?=?ˉ??u?t??
where?u?t?=??r=w?t?,p u?t?=p?r=w?t?.Also,the functions r0?t?, w?t?must satisfy the initial conditions
r0?0?=w?0?=R0
r˙0?0?=vˉ?2.53?
If p0?0?is known then the values of w˙,r˙0,and?may be deter-mined at t=0,r=R0,i.e.,if we set
p0?0+?=pˉu w˙?0+?=?ˉ??R0,0+?=?ˉ,r˙0?0+?=vˉ
?2.54?it follows from Eq.?2.36?,p u=?ˉ??u?,and Eq.?2.38?that
?ˉu=?ˉ?1?pˉu?,?ˉ2=p
ˉ?p0
?0?1??0/?ˉu?
vˉ=?ˉ?
?0?ˉ
?ˉu
??ˉ?1??0/?ˉu?
?2.55?
We note that the extra boundary conditions in?2.52?serve to compensate for the lack of another ordinary differential equation for r0?t?and/or w?t?;if such an equation existed then only two of the three boundary conditions in?2.52?would suf?ce to determine ??r,t?and p?r,t?explicitly.
2.1.3Nondimensional Formulation of the Problem.For the problem at hand the following dimensionless variables are now introduced
x=
r
R0
?=a t
R0
??pˉu?ˉu·t R0?2.56a?y???=
r0???
R0
???x,??=?/?ˉu p??x,??=p/pˉu?2.56b?
p?0=p0/pˉu??0=?0/?ˉu?2.56c?
w????=w/R0??=w?????=??ˉu pˉu??2.56d?
where?=d/d?.In terms of the nondimensional variables speci-?ed in Eqs.?2.56a?–?2.56d?our basic equations and conditions now read as follows:
?i?The Dynamic Continuity Condition:
p?u?p?0=w??????2??0?1???0/??u??2.57?where
??u???=???x,???x=w????p?u???=p??x,???x=w?????2.58?
and
??u?0+?=1p?u?0+?=1y?0+?=1,w??0+?=1?2.59??ii?The State Equation for the
Sediment:
Fig.3The domain D of compaction…adopted from Ref.?6?…
Applied Mechanics Reviews JULY2006,Vol.59/181
??u =?ˉ?1?p
?u ?,p ?u =?ˉ??
?u ??2.60?
with ?ˉ
?1?1?=?ˉ?1?=1.If we denote the initial values of the dimen-sionless shock wave velocity and spherical cavity velocity by ?
?and v
?,respectively,i.e.,?
?=?ˉ
a
v
?=v
ˉa
?2.61?
so that
w
????=0+=??,y ??t =0+=v ?then
?
?2=1?p ?0?
?0?1???0?,
v
?=???1???0??2.62?
as a direct consequence of Eqs.?2.55?and ?2.57?
?iii ?The Incompressibility Condition:This takes the form
?
?=f ?x 3?y 3?????2.63?
?iv ?The Dynamic Equilibrium Condition:
In terms of the dimensionless variables introduced,above,it is easily shown that Eq.?2.50a ?assumes the form
?
?p ??x =???y
x
?2?
y ?+2
?y ??2y ?1??y
x
?3
??
?2.64?
?v ?The Relation between the Dimensionless Velocity of the
Spherical Cavity Boundary and the Shock Wave Velocity
In terms of the dimensionless variables,Eq.?2.41?becomes
y 2
y ?=w
?2
w ???1???0/??u ??2.65?
In the x ,?phase plane the solution has the form shown in Fig.
4;at point A in this ?gure we have
?
?u ?A ?=1p
?u ?A ?=1y ??A ?=v
?w
??A ?=???2.66?
We now refer to Fig.5;at the spherical cavity r =r 0?t ??i.e.,x
=y ????we assume that the explosive decompression in the spheri-cal cavity is polytropic with generalized polytropy exponent ?so that
p 0?t ??43?r 0
3?t ?
??
=p ˉu
?4
3
?R 0?
?
=const.?2.67?
In terms of the dimensionless variables we have introduced we may express Eq.?2.67?as
p
?0???y 3????=1,p
?0?t ?=p 0?t ?
p ˉu
?2.68?
and we note that
p
??x ,???x =y ????p ?0???=y ?3????,for ???*?2.69?
p
??x ,???x =y ????p ?0???=p ?0=const,for ???*
where ?*denotes the ?dimensionless ?time at which the explosive
gas pressure reaches the initial pressure p
?0.2.1.4An Approximate Solution Method .For the initial-boundary value problem delineated above we now outline an ap-proximate solution technique.Note,?rst of all,that the dimen-sionless density ?
?which appears in the dynamic equilibrium equation ?2.64?is,according to Eq.?2.63?,an unknown function of the argument x 3?y 3???.To avoid having to solve the nonlinear system of partial and ordinary differential equations that was for-mulated above,a simpli?ed system of equations will be introduced.
We begin by integrating Eq.?2.64?over the interval from x
=y ???to x =w
????so as to obtain:y ?=?2?
1?R 2R 1y 3?
?y ??2y +1R 1y
?3??2?y ?2?
?u ???R 1
?2.70?
where
R 1?w
????,??u ???,y ??????y ???
w
?????
?x 2
dx =?y ???
w
?????x 3?y 3???+1?n
x 2
dx
?2.71a ?
R 2?w
????,??u ???,y ??????
y ???
w
?????
?x 5
dx =?
y ???
w
?????x 3?y 3???+1?n
x 5
dx
?2.71b ?
and the exponent n ???is given by
n ???=
ln ?
?u ???ln ?w
?3????y 3???+1??2.72?
Although n ???is unde?ned at ?=0?which corresponds to point A in Fig.4?it is possible to compute the limit of n ???as ?→0,
i.e.,
Fig.4…x ,?…-phase plane in nondimensional variables
…adopted from Ref.?6
?…
Fig.5Pressure change as a function of time inside the spheri-cal cavity …adopted from Ref.?6?…
182/Vol.59,JULY 2006
Transactions of the ASME
lim ?→0n????n A=
1
3
??A?
??0??
?2.73?
The next step is to solve Eqs.?2.57?and?2.58?for y?and w??;the result of this exercise is the following:
w??=???u?p?0
??0?1???0/??u????u=?ˉ?1?p?u??2.74a?and
y?=w?2
y2?1??0?p?u?p?0?
?1???0??u????u=?ˉ?1?p?u??2.74b?
or
w??=
y2y?
w?2?1???0/??u?
?2.75?
where
1
??0
?p?u?p0??1???0??u????u=?ˉ?1?p?u?=?y?y2w??2?2.76?Equation?2.76?can be solved for p?u in the form
p?u=f??y?2y4w?2??2.77?
Once the function f?has been determined,the relationship between ??u and y?2y4/w?2may be obtained from the constitutive relation ??u=?ˉ?1?p?u?in the form
1???0
??u
=g??y?2y4w?2??2.78?
in which case Eq.?2.75?becomes
w??=Y2y?/w?2g??y?2y4/w?2??2.79?The relations in Eqs.?2.70?and?2.79?constitute the basic sys-tem of equations for the functions y and w??the non-dimensional forms for the positions,at time t?0,of the spherical cavity and the shock wave front?.In Eqs.?2.70?and?2.79?the coef?cients R1,R2and the function g?are complicated nonlinear functions of y,y?,and w?.
2.1.5Numerical Formulation of the Problem.To solve,nu-merically,the problem of an expanding spherical shock wave,in a partially saturated soil,based on the formulation presented here, we introduce the vector Y k???,k=1,...,5,with components Y1=y????Y2=y???Y3=w????Y4=p?u,Y5=??u
?2.80?so that the state of the sediment at time??0is completely deter-mined if we know all components of the vector Y k.Based on the de?nitions in Eq.?2.80?,the system of differential equations which have to be solved have the form
dY m
d?
=F m?U1,...,Y5?m=1,2,3?2.81?where
F1?Y k?=?2?1?R2R1Y23?Y12Y2+Y2?3??Y4R1Y22?2.82a?
F2?Y k?=Y1?2.82b?
F3?Y k?=
Y1Y22
Y32?1???0/Y4?
?2.82c?with
R1=?Y2Y3?x3?Y23+1?n x2dx
R2=?Y2Y3?x3?Y23+1?n x5dx?2.83?
n=ln Y4/ln?Y32?Y23+1??0
We assume that the value of the vector Y k??k?,k=1,...,5,is known at time?=?k so that the solution of the system?2.81?serves to determine the components of Y k,k=1,2,3at time?=?k+??k??k+1;then,knowing the components Y k??k+1?,k =1,2,3,the components Y k??k+1?for k=4,5may be found from Eqs.?2.75?and?2.76?,i.e.,
1
??0
?Y4???0??1???0Y5??Y5=?ˉ?1?Y4?=Y12Y24Y34?2.84?
while from the constitutive relation,Y5=?ˉ?1?Y4?.
2.1.6Numerical Results.For the numerical scheme outlined above,results have been generated in Ref.?6?which provide some insight into the following:?i?the maximum level of compacting of the partially saturated soil?sediment?in a neighborhood of the spherical cavity,?ii?the initial shock wave speed,?iii?the depen-dence of the density at the shock wave speed front on the position of the shockwave,?iv?the dependence of the radius of the spheri-cal cavity on the position of the shock wave,?v?and the radius of the compacted zone in the sediment?where it is assumed that the compacted zone is the region in the sediment where the soil den-sity,once the dynamic processes which accompany the explosion are over,is1%greater than the initial density?;we also?nd,in Ref.?6?,results related to?vi?the average level of compacting in the compacted zone where??1.01?0=?*.The results presented in Ref.?6?are for partially saturated soils in which the contents of water and gas?volume fraction?change in the range of0.1–0.2; we will now review some of the results that were generated by the numerical computations in Ref.?6?.
The?rst result assumes that the average pressure pˉu?Fig.5?generated by the explosive charges is known as well as the initial density of the soil?sediment?and the initial volume fractions of the soil components?water and gas?.A global density correspond-ing to pˉu may determined from Eq.?2.21?.The maximum density of the soil at the surface of the spherical cavity is reached at the initial time t=0,when r=R0,and is re?ected in the dimensionless parameter??0where
??0=?k=13
??k?1p?0?1?+1??1/?k v k0?2.85?
In Fig.6we display?using the data in Table1?a graph depicting the variation of??0with??0;from this graph,?ˉu can be determined provided that the soil properties and the pressure generated by the explosion are known.
An analysis of the basic equations,expressed in dimensionless form,shows that the?dimensionless?initial shock wave velocity ??=
??ˉ2/?
2,where?
2=pˉ
u
/?ˉu,the?isothermal?sound velocity for an ideal gas at pressure p=pˉu,is a function of p?0,v?k?0,?k,and ?k,i.e.,
Applied Mechanics Reviews JULY2006,Vol.59/183
??=???p ?0,v ?k ?0,?k ,?k ?
k =1,2,3?2.86?
Also,?
?is determined by the formula ??2=
1?p ?0??0?
1??p ?0?1/?3v ?3?0??
p ?0?2+p
?0?1??2??
1/?2
v ?2?
?
?
p
?0
?1+p
?0?1??1??1/?1
v ?1?
?
?1
?2.87?
whose graph,based on the data in Table 1,is depicted in Fig.7.
From Fig.7it is possible to determine the initial shock wave speed provided one knows the soil properties as well as the pres-
sure generated by the explosion.
The relationship between the dimensionless density p
?u and the dimensionless pressure p
?u is given by the following dimensionless form of Eq.?2.21?:
?
?=??0??
k =1
3v ?k ?
??k
??
?p
?0?1?+1?
?1/?k
?
?1
=?ˉ?1?p
???2.88?If we take into account Eq.?2.85?,it is easily seen that the relation
in Eq.?2.88?depends only on the dimensionless parameter p
?0and the constants v ?k ?0
,?k ,and ?k ,k =1,2,3;this relation is graphed in Fig.8,which has been obtained by using three values of the
parameter p
?0:p ?0=10?4,10?5,and 10?6for various water and gas contents in the sediment.From Fig.8,one may determine the pressures at the shock wave front.
As has already been indicated,density changes at the shock wave front may be calculated up to a maximum time t =t M at which the density ?u reaches a speci?ed value ?*;the points M ?in Figs.2and 4represent the position of the shock wave front at this
time when ?u =?*?i.e.,??u =?*/?ˉu ?,while the distance of the points
M ?from the source of the explosion determines the radius w M in
Fig.9and the dimensionless radius w
?M =w M /R 0.The limiting value of the radius of the spherical cavity at time t M is denoted by r 0
M in Fig.9while,in Fig.4,y M ?r 0M /R 0.Within the context of
the Fig.6Variation of ??0with ??0…adopted from Ref.?6?…1?v …2…
=0.099,v …3…0=0.001,2?v …2…0=0.095,v …3…v =0.005,3?v …2…0=0.09,v …3…
=0.01,4?v …2…0=0.199,v …3…0=0.001,5?v …2…0=0,v …3…0=0.005,6?v …2…
=0.19,v ...3 0
=0.01
Table 1Values of the components in a partially saturated soil …
adopted from Ref.?6?…
1-Quartz
2-Water 3-gas ?k
?0? 2.65·10?4kG s 2/cm
41·10?4kG s 2/cm 40.00125·10?4kG s/cm 4c ?k ?4500m/s
1500
m/s
330m/s ?k
5
7
1.4
Fig.7Graph of ?
?as a function of p ?0…adopted from Ref.?6?…Fig.8
p
?u as a function of ??u …adopted from Ref.?6?…Fig.9The radius w M …adopted from Ref.?6?…
184/Vol.59,JULY 2006Transactions of the ASME
theoretical model presented above the zone of compaction is de-?ned to be the spherical shell of inner radius r 0M
and outer radius
w M ?r 0M ?in dimensionless variables y M and w ?M ?y M ?;both y M and w ?M are functions of the parameter p ?0,the constants ?k ,v ?k ?0,?k ,and the polytropic exponent ?of the detonation prod-ucts.The relations between w ?M ,p ?0,and v ?k ?0
are shown in Figs.10
and 11;the calculations leading to these graphs are based on the constants in Table 1and a value of ?=3.One major result is the following:with the same average explosion pressure,the zone of the compacted soil ?or sediment ?for a dry soil is greater than that for a partially saturated soil.The results displayed in Figs.4and 5allow one to do a theoretical evaluation of the dimension of the zone of compaction of a soil ?sediment ?which has been subjected to the explosion of an embedded spherical charge;the analysis also leads to an indirect evaluation of the average level of com-paction.The average density ?m of the region y M ?y ?w
?M is given by
?m =?0??w M ?3?R 03?
?w M ?3??r 0M ?3=?0?w
?M ?3?1w
?M
??y M ?3?2.89?
An average compaction level z =?m /?0may then be computed
from the formula
z =?w ?M ?3?1?w
?M ?3??y M ?3?
1
1??y M /w ?M ?3?2.90?
In Fig.12we display the dependence of the average compaction level upon the dimensionless average explosion pressure p 0and the water and gas contents in the soil for the constants given in Table 1and a value of ?=3.Figure 13depicts the density de-crease at the shock wave front as a function of the position of the wave;the same ?gure also shows the increase in the radius of the
spherical cavity y as a function of the increasing distance w
?of the wave front from the cavity.At those points where ?u =?*?1.01?0,the level of water and gas contents appears to exert no in?uence while the change in the radius of the spherical cavity
appears to depend on the value of v ?2?0+v ?0?0
.In Fig.14we show the radius of the spherical cavity and the position of the shock wave as functions of the initial volume fraction of gas vapor in a pore.
The shape of the curves ?
?u ?w ??clearly shows that ?locally ?the level of soil compaction around the cavity may be very large and
that there is,close to w
?=1,a large gradient associated with the drop in density.Of course,various kinds of physical phenomena which take place in a neighborhood of the cavity are complex and are not entirely accounted for in the model described above,e.g.,soil degradation and water ?ltration which results from the nega-tive pressure in the ?nal phases of the decompression of the ex-plosive charges.However,the analysis presented here does allow for the determination of the dimension of the zone of weak com-paction as well as the level of compaction which results from
a
Fig.10
w M as a function of p
?…0……adopted from Ref.?6
?…Fig.11
y M as a function of p
?…0……adopted from Ref.?6
?…Fig.12The average compaction level as a function of the ex-plosion pressure …adopted from Ref.?6
?…
Fig.13Distribution of the density inside of the compacted
zone …adopted from Ref.?6?…
Applied Mechanics Reviews JULY 2006,Vol.59/185
single explosion of a spherical charge.
To summarize the results of the analysis of the model presented above,a single explosion of a spherical change in a partially satu-rated soil generates mechanical compaction of the soil in a spheri-cal ring of dimension
r 0M ?r ?w M
?2.91?
where the radius w M corresponds to the position of the shock
wave front at which ?=1.01?0and the parameter z =?m /?0,which is given by Eq.?2.90?,is a suitable measure of the degree of compaction resulting from the explosion in a neighborhood of the explosive charge.Within the spherical ring given by Eq.?2.91?,z determines an average ratio of soil density after the explosion to soil density prior to the explosion.The ?nal density in the spheri-cal cavity r ?r 0M
,after the explosion,cannot be determined by the method described above as the phenomena which occur in this zone during and after the explosion are extremely complex and dif?cult to formulate.The ?nal density is in?uenced by the pro-cess of soil degradation in a neighborhood of the cavity as well as by the process whereby the cavity is ?lled due to gravitational forces after decompression of the explosive charge is completed.Also there may be a signi?cant in?uence on the ?nal density in
the region r ?r 0
M
due to dynamic water ?ltration when a negative pressure exists in the cavity.
In Tables 2?a ?and 2?b ?,results are given for two different types of explosive charges and two different values R 0of the radius of the spherical cavity prior to the explosive detonation;as the actual value of p ˉu is unknown,these tables present results which corre-
sponds to assuming that the initial pressure is equal to the maxi-mum pressure at the front of the detonation wave.Tables 2?a ?and 2?b ?show an average ratio for the soil density after the explosion to the soil density in the compaction zone prior to the explosive detonation in the range of 1.019–1.089;in this range,the in?uence of the water content in the soil,and the computed value of p ˉu ,upon the ratio z appears to be insigni?cant but both factors do have a signi?cant in?uence on the radius of the zone of compac-tion,i.e.,for the case of an explosive charge,using the minimum value ?pressure ?generated by the dynamite,it is found that the radius of the zone of compaction varies in the range 0.72m ?w M ?5.49m while for the maximum value ?pressure ?generated by the dynamite charge it varies in the range 0.9m ?w M ?7.01m.
Remarks:In Ref.?8?the global constitutive relation ?2.21?be-tween p ?t ?and ??t ?was written explicitly in the form
?0??t ?
=v 10?9.32·10?5?p ?p 0?+1??1/5v 20
?3.11·10?3?p ?p 0?+1??1/7+v 30?10.285?p ?p 0?+1?
?1/1.4?2.92?
by using the values for the ?i ,c i ,and ?i ,i =1,2,3given in Table
1;as the p =p ???relation is very sensitive to a change in the parameter v 30,even a small change in the content of gas in the soil results in signi?cant changes in the pressure density relation.Sim-pli?ed equations of state were also introduced in Ref.?8?for those cases where the pressure values are contained in the interval 10MPa ?p ?103MPa:
?0??t ?
?v 10+v 20=1?v 30
?2.93?
and the authors indicate that for the range of pressures indicated above the error involved in using Eq.?2.93?does not exceed a few percent.In Fig.13,which is reproduced from Ref.?8?,we have depicted the distribution of the density within the compacted zone as predicted by the model represented by Eq.?2.92?,shown with the solid curves,and,also,as predicted by the simpli?ed model ?2.93?,shown with the dashed curves.It may be noted that ?i ?the
smaller the value of v 30the greater the propagation distance w
?max of the shock wave front and ?ii ?the larger the value of v 30the bigger the size the spherical cavity,i.e.,r
?0max
,is.Figure 14,which we have already referenced,is also taken from Nowacki and Gue-lin ?8?,and shows the radius of the spherical cavity and the
posi-
Fig.14Radius of the spherical cavity and the position of shock wave versus the initial volume fraction of the gas …adopted from Ref.?6?…
Table 2Compaction levels versus explosive pressure …adopted from Ref.?6?…
Explosive Average explosion pressure
Initial spherical cavity radius Water content in soil
Air
content in soil
Final
radius of cavity
Radius of compacted zone V olume of compacted zone Average compaction level
p ?*??MPa ?R 0?m ?v ?2?
0v ?3?
0r 0M ?m ?
w M ?m ?V p ?m 3?z ?a ?
Dynamite 505A 1100
0.150.0990.0010.320.72 1.50 1.0860.190.010.55 1.7120.12 1.0340.33
0.0990.0010.71 1.6115.90 1.0890.190.01 1.21 3.75214.2 1.034Granulated TNT 7111
0.150.0990.0010.380.98 3.66 1.0580.190.010.67 2.4963.71 1.0200.33
0.0990.0010.83 2.1538.99 1.0570.190.01 1.47 5.49678.4 1.019?b ?
Dynamite 505A 2200
0.150.0990.0010.360.90 2.86 1.0630.190.010.63 2.2144.3 1.0230.33
0.0990.0010.79 1.9830.4 1.0630.190.01 1.39 4.87471.9 1.023Granulated TNT
14222
0.150.0990.0010.45 1.227.20 1.0510.190.010.81 3.19133.43 1.0170.33
0.0990.0010.99 2.6876.68 1.0510.190.01 1.78
7.01
1420.7
1.021
186/Vol.59,JULY 2006Transactions of the ASME
tion of the shock wave front as functions of v30for both the full and simpli?ed equations of state;both positions are computed at the time t m when p0?t m??p lim=10MPa for the case of an initial pressure of100MPa.While there are signi?cant differences be-
tween the values of the cavity radius for different values of v30?at the same initial pressure p0?there are considerable differences in the dimensions of the compacted zones for different values of v30 and the same initial pressure p0.
2.2General Formulation of the Problem in Space Dimen-sion n.In the following we will generalize the analysis presented in Sec.2.1by employing a formulation due to Wlodarczyk?7?;in this work the theory of concentrated explosion is based on the assumption that all the explosive energy is produced in a region of ?essentially?zero volume which can be a point,a line,or a plane for the case of spherical,cylindrical,and plane symmetry,respec-tively,where in the two latter cases the energy is considered as a quantity per unit length of area.By idealizing the initial condition to be one of in?nite energy density at the center of the explosion, and neglecting the value of the pressure in front of the shock wave which results,closed-form solutions may be obtained which de-termine the characteristics of the nonstationary motion of a perfect gas due to the explosion of a concentrated charge.Although the explosive energy produced at time t=0is high?but?nite?the idealized model presented in Ref.?7?appears to model quite well the one-dimensional motion of saturated soil due to a sudden pro-duction of energy in a spherical or cylindrical region,or in a plane layer of the soil;the basic constitutive assumption employed in a neighborhood of the region in which the explosive charge is deto-nated is,once again,that of a plastic gas.
2.2.1Basic Equations and Assumptions.We are modeling the nonstationary motion of a mechanical system consisting of an inclusion characterized by spherical,cylindrical,or plane symme-try?a ball,an in?nite cylinder,or an in?nite layer,respectively?which is immersed in an in?nite space?lled by undisturbed,ho-mogeneous,water and gas saturated soil.It is assumed that a?nite amount E0of energy is produced over the entire region of the inclusion by the sudden detonation of an explosive charge.By r0 we denote the radius of a cylindrical?or spherical?charge or half the thickness of a plane charge.
The soil is modeled by an equation of state of the form
?0
??t?=f?v3
0??2.94?
where?0,??t?,and v30have the same interpretation as in Sec.2.1. The function f?v30?is,in turn,de?ned from a model in which ?0/??t?=??p?t??,i.e.,either a model of the type?2.19?,introduced by Lyakhov?1?,for which
?L?p?=v10??1?1c12?p?p0?+1??1/?1+v20??2?2c22?p?p0?+1??1/?2
+v30??3?3c32?p?p0?+1??1/?3?2.95?or that considered by Rakhmatulin?3?in which
?R?p?=v30?p0+p
p0+?p
+v20??2?2c22?p?p0?+1??1/?2
+v10??1?1c12?p?p0?+1??1/?3?2.96?
where?=??3+1?/??3?1?.The v i0,?i,?i,and c i all have the same meaning as in Sec.2.1;for typical values of these parameters we may again employ the data given in Table1.From the Rakhmatu-lin?3?model?2.96?,therefore,and the data in Table1,we have
?R?p?=v30?6p0+p p0+6p?v20?3.11·10?3?p?p0?+1??1/7
+v10?9.32·10?5?p?p0?+1??1/5?2.97?If we again consider the range of pressures,10MPa?p ?103MPa,then the full Lyakhov model may be approximated by the relation?2.93?,while
?R?p??1v30+v20+v10=1?5v30?f?v30??2.98?Using the laws of conservation of mass,momentum,and energy at the shock wave front generated by the detonation,and assuming that the shock wave is propagating into an undisturbed region of the soil ahead of the shock wave front,we have
v u=?1??0?u
??
p u?p0=?0v u??2.99?
e u?e0=
p u+p0
2?1?0?1?u
?
where e is the speci?c energy of the medium?sediment?,?,p,?,v have the same meanings as in Sec.2.1,and the subscripts“u”and “0”again indicate values at the wave front,and in front of the wavefront,respectively.The equation of state?2.94?is also satis-?ed at the shock wave front;the last relation in Eq.?2.99?may be employed to determine the energy losses at the shock wave front. For the continuous motion behind the shock wave front we have the following equations,which are the expressions of the local laws of conservation of mass and momentum
?
?r
?r+u?r,t??=?r r+u?r,t????1?0??2.100?and
?p
?r=??0
?r
r+u?r,t????1·?2u?t2?2.101?where?is a“symmetry”factor which has the values?=1for plane symmetry,?=2for cylindrical symmetry,and?=3for spherical symmetry.In order to establish a model for the motion of the gaseous products of a detonation,we begin by assuming that those products are not mixed,during the motion of the water saturated soil,with the plastic gas which is used to model the soil in a vicinity of the detonation.In Ref.?7?the initial values of the parameters which determine the state of the gaseous products of the explosion are?xed by using the theory of“immediate detona-tion”of the explosive charge;the use of such a theory is often justi?ed,in engineering practice,by noting that the contribution of the kinetic energy to the energy balance for the gaseous?explo-sion?products in a detonation wave is quite small.As has been indicated in Ref.?7?,if one considers the fraction of the potential and kinetic energies for the gaseous products of explosion in a detonation wave for the explosive trinitrotoluene,then the highest fraction of kinetic energy?8%?occurs in the case of spherical symmetry so that more than90%of the energy of the gaseous products of an explosion is potential energy?either elastic or ther-mal?.As a consequence,it may be assumed that particles of the gaseous products do not move so that the initial density of the gaseous products will be the same as the density?0e of the explo-sive.The initial pressure p0e is then given by
p0e=
1
2
p M=
1
2
·
?0e?2
k H+1
?2.102?
where p M is the pressure of the explosive gases at the Jouguet point,?is the velocity of the detonation?shock?wave,and k H is
the isentropic exponent of the explosive gases at the Jouguet point.If we then neglect wave phenomena in the explosive gases, and assume that they expand in an adiabatic manner,the mean density?e?t?and mean pressure p e?t?may be determined from
?e?t?=?0e?R0r0?t????2.103?
p e?t?=p0e?R0r0?t???k0
where R0is the initial“radius”of the cavity enclosing the explo-sive,r0?t?is the“radius”of the cavity at time t,r0?t?=R0 +u?R0,t?,and k0is the?correlated?polytropic exponent of the explosive gases.Examples of k0and k H for various explosives are given in Table3,which is reproduced from Ref.?7?.
From the conditions of continuity of the pressure and velocity at the?contact?boundary of the explosive gases we have
p?R0,t?=p e?t?v?R0,t?=v e?t??2.104?
where v=?u/?t.We also assume that the region ahead of the shock wave front is undisturbed so that
u?r,0?=v?r,0?=0?2.105?
Also,if we note that the simpli?ed models of both Lyakhov and Rakhmatulin,which we have delineated above,differ only in the value of the multiplier of v30,then we may employ a general equa-tion of state of the form
?0
??t?=1???const?2.106?
where?=v30for the Lyakhov model and?=?5/6?v30for the Ra-khmatulin model.
We now want to obtain the general solution of the problem posed above;that solution will then be specialized to the situa-tions in which?=1,2,or3.We?rst write the continuity equation ?2.100?in the form
1??
?r
?r+u??=
?0
?r
??1?2.107?
Integrating Eq.?2.107?we obtain
?r+u??=??R0r
??0??r??1dr+f1?t??2.108?
As R0+u?R0,t?=r0?t?we must have f1?t?=r0??t?so that Eq.?2.108?becomes
?r+u??=?1????r??R0??+r0??t??2.109?
We again let r=w?t?denote the position of the shock wave front at time t?0.Substituting r=w?t?into Eq.?2.109?,and using the fact that u?w?t?,t?=0,t?0,at the shock wave front,we have
w??t?=
??1
?R0
?+1
?r0
??t?=R
?+1
?
?r0??t??R0???2.110?
The velocity and acceleration of the soil particles are next ob-tained from Eq.?2.109?by repeated differentiation with respect to t,i.e.,
v?r,t??
?u?r,t?
?t=??1????r
??R
??+r
??t???1???/?r
??1?t?r˙
?t?
?2.111?and
a?r,t??
?2u?r,t?
?t2=??1????r
??R
??t????1???/?
?????1?r0??2?t?r˙02?t?+r0??1?t?r¨0?t??+?1???
???1????r??R0??+r0??t???1/2??/?r02???1??t?r˙02?t?
?2.112?If we now integrate Eq.?2.101?,and make use of the second condition in Eq.?2.103?,we obtain for the pressure the relation p?r,t?=??0?R0r
??ˉ?ˉ+u??ˉ,t????1u,tt??ˉ,t?d?ˉ+p e?t?
?2.113?However,from the relations?2.99?at the shock wave front,it follows that
v u?w?t?,t?=?w˙?t?
p u?w?t?,t?=p0+?0?w˙2?t??2.114?where w˙?t?=??t?.From the collection of relations delineated above it is easily seen that all of the quantities required to com-pletely characterize the solution of the problem at hand may be determined in terms of the function r0?t?and its derivatives r˙0?t?and r¨0?t?;the function r0?t?will be computed,below,for each type of symmetry?plane,cylindrical,and spherical?which may be of interest with respect to detonation of an explosive charge in a partially saturated soil.
2.2.2Dimension n=1:Plane Shock Waves.For plane shock waves?dim n=1?we employ?=1in the basic equations delin-eated above;this yields the following relations
u?r,t?=r0?t??R0???r?R0?
?u
?t
?r,t?=v?r,t?=r˙0?t??2.115?
?2u
?t2
?r,t?=a?r,t?=r¨0?t?
w?t?=R0+
1
?
?r0?t??R0??2.116?
w˙?t?=
1
?r
˙0?t??2.117?
v u?t?=v?w?t?,t?=r˙0?t??2.118?p u?t?=p0+?0?w˙2?t?=p0+
?0
?r
˙
2?t??2.119?and
Table3Key explosive parameters…adopted from Ref.?7?…
Explosive
?oe
kg/m3
p H
MPa
D
ms k H k0
Trinitrotoluene163021,0006930 2.73 3.00 TH36/64171729,5007980 2.71 3.00 Octogen189142,0009110 2.74 3.40 Pentryt177033,5008300 2.64 2.90 Nitromethane112812,5006280 2.56 2.73
p?r,t?=p e?t???0?r?R0?r¨0?t??2.120?In order to determine r0?t?,we substitute Eq.?2.120?into Eq.?2.119?using the fact that Eq.?2.120?is satis?ed at the shock wave front.After some algebraic manipulations,we arrive at the relation
?r0?t??R0?r¨0?t?+r˙02?t?=?
?0
?p e?t??p0??2.121?
or
d2
dt2??r0?t??R0?2
2?=??0?p e?t??p0??2.122?
Upon integrating Eq.?2.122?we obtain
r0?t?=R0+?2??0?0t?0??p e????p0?d?d??1/2?2.123?
which represents a closed form solution to the plane shock wave detonation problem if the pressure in the cavity is known as a function of time.A typical example in this case would be that of a pulse load which is generated by the explosion of a concentrated explosive charge;for such a case,p e?t?would be given in the form
p e?t??p0=?p m?1?t/???q0?t???
0t???
?2.124?
for some???0.
As in our earlier analysis of the concentrated spherical explo-sion in a gassy sediment,we now introduce the following dimen-sionless quantities for the problem associated with Eq.?2.124?
?=r
a0???0=R0
a0??
?a0=?p0/?0?,?=t?????,??=u?r???,t?
???
a0??
r0*???=r0?t??
a0??
w*???=
w?t????
a0??
?2.125?
r˙0*???=r˙0?t??
a0
w˙*???=
w˙?t????
a0
r¨0*???=
r¨0?t????
a0/??
P??,??=p?r???,t????
p0
P e???=
p e?t????
p0
P m=
p m
p0
In terms of the dimensionless variables in Eq.?2.125?,the rela-tions?2.115?–?2.120?take the form
U??,??=r0*?????0??????0?
U,???,???V??,??=r˙0*???=V u????2.126?
U,????,???A??,??=r¨0*???
w*???=?0+1
?
?r0*?????0?
P u???=1+??w˙*????2?2.127?
p??,??=P e????????0?r¨0*???
where
r0*???=?0+
?2?
q+1
·P m?1q+2??1???q+2?1?+?
?r˙0*???=
1
r0*?????0??q+1·P m?1??1???q+1?
??2.128?r¨0*???=
1
r0?????0
??P m?1???q??r˙0*????2?
for0???1,and
r0*???=?0+
??
q+2
P m?const
?2.129?
r˙0*???=r¨0*???=P e???=0
for??1.By extracting the limits as?→0,we determine that
r˙0*?0?=??P m,r¨0?0?=?13q??P m?2.130?
An analysis of the equations presented above is shown,in Ref.?7?,to lead to the following conclusions.
?i?Under the plane shock wave propagation conditions that are described above,the motion of the surface of the cavity and the shock wave front decay.
?ii?The maximum radius of the cavity is given by
?r0?max=R0+
2?
q+2
·
p m
p0
·a0???2.131?
?iii?The maximum penetration depth of the expanding shock wave front into the gassy sediment is given by
w max=R0+
2
q+2
·
p m
p0
a0???2.132?
From Eqs.?2.131?and?2.132?,it follows that the thickness of the compacted layer of sediment is given by
w max??r0?max=
2?1???
q+2
·
p m
p0
a0??=
2
q+2
·
?0
?·
p m
p0
·a0??
?2.133??iv?The value of the“overpressure”at the shock wave front, i.e.,?P u=P u?1,is given,by virtue of the?rst two equations in Eq.?2.127?,by
?P u=P m
2?q+1?
·
?1??1???q+1?2
1
q+2
??1???q+2?1?+?
?2.134?
so that?P u is strongly affected by the initial value of the pressure in the cavity,as well as by the manner in which the pressure in the cavity decays,i.e.,by the parameter q.In Fig.15we reproduce results of Wlodarczyk?7?,which depict the variation of?P u/P m as a function of?for various values of the parameter q.
?v?From Eqs.?2.127?and?2.128?it follows that the velocity of the medium increases in direct proportion to??=?1??0/?.The variation of V???,as given by Eq.?2.126?,as a function of?for various measures?of the air content of the soil and various values of the parameter q,is depicted in Fig.16.If q=0,so that the pressure in the cavity is constant,then the velocity also pre-serves the constant value
V??,???V u????r˙0*???=??P m?2.135?and
w *???=?0+
?
P m
?
?r 0*???=?0+??P m ?
?2.136?
In Fig.17we show the form of the wave trajectory for various
values of the coef?cient ?,and the exponent q ,in Eq.?2.124?;with increasing ?,i.e.,with increasing air content,the penetration depth of the shock wave front decreases very rapidly and the maximum penetration depth of the shock wave front into the sedi-ment is then determined by Eq.?2.132?.The variation of the func-tion U ???,as given by Eq.?2.126?,at the boundary of the cavity,i.e.,at ?=?0,is depicted in Fig.18for a few ?xed values of the parameters ?and q .
Now,suppose that only the initial detonation parameters of a particular explosive are known and that,in accord with Eq.?2.103?,for a plane shock wave with ?=1,
p e ?t ?=p m
?R 0r 0?t ?
?
k 0
p m =0.5p H ?2.137?
Substituting Eq.?2.137?into Eq.?2.123?we obtain for r 0?t ?,the radius of the cavity at time t ,a nonlinear V olterra integral equation of the second kind,i.e.,
r 0?t ?=R 0+
?
2??0
?0t
?0
?
?p m
?R 0r 0???
?k 0
?p 0?
d ?d ?
?
1/2
?2.138?
which can be solved either by an iteration technique,as in Wlo-darczyk ?7?,or numerically;the ?rst two approximations in the iterative scheme in Ref.?7?have the form
r 01?t ?=R 0+d 0t
d 0=
?
?·
p m ?p 0
?0
?2.139a ?
and r 02?t ?
=R 0+2?R 0
?
p m
?k 0?1?d 02
?0?1k 0?2
??
R 0R 0+d 0t
?k 0?2
?1
?
+d 0t R 0?
?p 02?0·t 2R 0
2?
1/2
?2.139b ?
It was also demonstrated,in Ref.?7?,that for p e ?t ?as given by Eq.?2.137?,Eq.?2.122?may,in fact,be solved in closed form as follows:we introduce the change of variables given by
x ?t ?=r 0?t ??R 0x ˙?t ?=r ˙0?t ??2.140?
y ?x ?=r ˙02?t ?=x
˙2dy
dx
=2x ¨Substituting Eqs.?2.137?and ?2.140?into Eq.?2.122?,we ob-tain
xy ??x ?+2y =2?
p m ?0?R 0
x +R 0
?
k 0
?2?
p 0?0
?2.141?
Introducing the new dependent variable z ?x ?=xy ?x ?,we may re-write Eq.?2.141?in the form
?xz ?x ???=2?
p m ?0·x ?R 0
x +R 0
?
k 0
?2?
p 0
?0
x ?2.142?
At the initial time ?which corresponds to x =0?r 0?0?=R 0,and the function y ?x ?,as given by Eq.?2.140?,assumes the value
y ?0?=r ˙02?0?=?
p m ?p 0
?0
?2.143?
at the surface of the cavity,i.e.,at r =r 0?t ?.The initial value of z ?x ?is then z ?0?=0.Thus,integrating Eq.?2.142?subject to z ?0?=0we ?nd that
xz ?x ?=x 2
y ?x ?=
?
2?R 0
2?k 0?2??k 0?1?·p m ?0?
?k 0?2?x R 0?R 0x +R 0
?
k 0?1
+
?R 0
x +R 0
?
k 0?2
?1
??
??
p 0?0
x 2?2.144?
for k 0=3,4,5,...,while for k 0=2we have
xz ?x ?=2?
p m ?0R 02?ln R 0+x R 0?x R 0+x
?
??p 0?0x 2?2.145?
and for k 0
=1,
Fig.15Variation of ?P u /P m with ?…adopted from Ref.?6
?…
Fig.16Variation of V with ?…adopted from Ref.?6?…
xz ?x ?=2?p m ?0R 0
2?x
R 0?ln R 0+x R 0
?
??p 0?0x 2?2.146?
It then follows that the velocity of the partially saturated soil
behind the shock front may be expressed,in view of Eqs.?2.140?,?2.144?–?2.146?,and ?2.115?,in the form
r
˙0?t ??v ?t ?=v u ?t ?=?
2??k 0?2??k 0?1?p m
?0?R 0r 0?t ??R 0
?2
??1
??k 0?2?
?
r 0?t ?
R 0
?1??R 0r 0?t ?
?k 0?1
?
?R 0r 0?t ?
?k 0?2
???
p 0
?0
?
1/2
?2.147?
for k 0=3,4,5,...,while for k 0=2,
r
˙0?t ?=?2?p m ?0?
R 0
r 0?t ??R 0
?2?
ln r 0?t ?R 0+R 0
r 0?t ?
?1??
p 0?0?
1/2
?2.148?
and,for k 0=1,
r
˙0?t ?=?2?p m ?0?
R 0
r 0?t ??R 0
?2
?
r 0?t ?
R 0?1?ln r 0?t ?R 0
???
p 0?0?
1/2
?2.149?
From Eq.?2.121?we next obtain the second derivative
r ¨0?t ?=1r 0?t ??R 0
?
?r ˙0
2?t ?+?p m ?0?R 0
p 0?t ??
k 0
??
p 0?0
?
?2.150?
At this point all the quantities that characterize the state and mo-tion of the sediment,including the motion of the shock wave front,have been determined as functions of the cavity thickness r 0?t ?.
From the relations ?2.140?it follows immediately that
r
˙0?t ?=?y ?r 0?t ??t ?0
?2.151?
so that
t =
?
R 0
r 0?t ?
d ?
?y ???
t ?0?2.152?
For k 0=3?which is quoted in Wlodarczyk ?7?as being a common
parameter value ?we obtain,as being a consequence of Eq.?2.147?,
r
˙0?t ?=???p m ?0?R 0r 0?t ?
?2
?p 0?0
?
1/2
t ?0?2.153?
or
Fig.17Form of the wave trajectory …adopted from Ref.?6
?…
Fig.18Variation of U with ?at the cavity boundary …adopted
from Ref.?6?…
??t =
?R 0
r 0?t ??p m ?0?R 0?
?2
?p 0?0
?
1/2
d ??2.154?Carrying out th
e integration indicated in Eq.?2.154?,we arrive at
the following result for k 0=3:for t ?0,
r 0?t ?=R 0
??p m
p 0????
p m
p 0
?1???
a 0t R 0?2?
1/2
?2.155?
where a 0=a ?R 0,0?.As a direct consequence of Eq.?2.155?,it
follows that
r 0?t ??R 0
?
p m p 0
t ?0?2.156?
From Eqs.?2.155?and ?2.156?we now easily ?nd that the maxi-mum “radius”of the cavity is
?r 0?max =
?
p m
?0
R 0?2.157?
and that this value is reached at time
t max =
R 0
a 0
?1??
p m
p 0
?1?
?2.158?
in which case the “maximum”distance traveled by the shock
wave front is
w max =?
1+
1
?
??p m /p o ?1??
R 0?2.159?
At the time t =t max ,the thickness of the compacted region in the partially saturated soil is given by
w max ??r 0?max =
?1?
?1???p m
p 0
?1?
R 0?2.160?
This completes the solution,in the plastic gas domain,for the
propagation of a plane shock wave in a partially saturated soil due to the detonation of an explosive charge in a cavity of “radius”R 0;we now consider the analogous problem for cylindrically symmet-ric shock waves which are generated by the detonation of a cylin-drical charge embedded in the sediment.
2.2.3Dimension n=2:Cylindrical Shock Waves .For a cylin-drical shock wave the relations ?2.109?,?2.111?,and ?2.112?,with ?=2,assume the form
u ?r ,t ?=??1????r 2
?
R 02?+r 02?t ??r ?2.161??u ?r ,t ??t =r 0?t ?r
˙0?t ???1????r 2?R 02?+r 0
2?t ??v ?r ,t ??2.162?
and
?2u ?r ,t ??t =?r 02?t ?r ˙02?t ????1?b ??r 2?R 02?+r 02?t ??
3+
r ˙02?t ?
+r 0?t ?r
¨0?t ???1?b ??r 2?R 02?+r 0
2?t ??a ?r ,t ??2.163?
while the position of the shock wave front at time t is given by Eq.
?2.110?with ?=2,i.e.,if
r =w ?t ?
t ?0
?2.164?
denotes the coordinate of the shock wave,then w ?t ?=
?
R 02+
1??r 0
2
?t ??R 02?t ?0?2.165?
The velocity of propagation associated with an expanding cylin-drical shock wave front is then
??t ??w
˙?t ?=r 0?t ?r
˙0?t ??
?
R 02+1?
?r 02?t ??R 02
??2.166?
which follows by direct differentiation of Eq.?2.165?.Using Eq.?2.166?the kinematic and dynamic compatability conditions ?2.99?,we ?nd for the velocity and pressure at the shock wave front the expressions
v u ?t ?=???t ??
r 0?t ?r ˙0?t ??
R 02+1?
?r 02
?t ??R 02??2.167?
and
p u ?t ??p 0+?0??2
?t ?=p 0+?0
r 02?t ?r ˙02?t ?
???1?R 02+r 02
?t ?
?2.168?
We now substitute Eqs.?2.161?and ?2.163?into expression ?2.113?for the pressure ?eld in the sediment so as to deduce that
p ?r ,t ?=p e ?t ??
?02?1???
??
?r ˙02?t ?+r 0?t ?r
¨0?t ???ln ?
?1????r 2?R 02?+r 02
?t ?
r 02?t ???
+?r 02?t ?r ˙02?t ??
?
?
1
?1????r 2
?R 02?+r 02?t ?
?
1
r 02?t ?
??
?2.169?
Once again,all those quantities which are necessary in order to characterize the state and motion of the sediment at any time t ?0,when the sediment in a vicinity of the detonation is in the plastic gas regime,have been expressed,explicitly,in terms of the function r 0?t ?,which governs the motion of the ?cylindrical ?cav-ity,and its derivatives;in this manner the problem has been re-duced to one of determining the function r 0?t ?.
To obtain an expression for r 0?t ?,we substitute Eq.?2.169?into Eq.?2.168?;this leads to a quasi-linear second order ordinary differential equation for r 0?t ?of the form
r ¨0?t ?+??r 0?t ??r ˙02?t ?=??r 0?t ??
?2.170?
We now set y =r ˙02?t ?,so that
dy dt =2r
¨0?t ?r ˙0?t ?Tdy
dr 0
=2r ¨0?t ?,?2.171?
in which case Eq.?2.170?may be rewritten as a linear ?rst order
equation for y =y ?r 0?:
dy
dr 0
+2??r 0?y =2??r 0??2.172?
The initial condition associated with Eq.?2.172?is
y 0?y ?r 0?=r ˙02?0?=v u 2?0?=??2
?0???
p m ?p 0
?0
?2.173?and it is not dif?cult to show that the functions ?and ?in Eq.?2.170?have the following forms:
??r 0?=1
r 0
+
??
ln
R 02+1?
?r 02?R 02
?r 0
2??1
?
?
r 0
R 0
2+1?
?r 02?R 02
?+
2?1???r 0???
1?R 0
2?r 0
2?
1r 0
??
?2.174?
and
??r 0?=
2?1???
?0
?
p m
?R 0r 0
?
2k 0
?p 0
r 0
???
ln
R 0
2+1?
?r 02?R 02
?r 0
2?
?1
?2.175?
If we now integrate Eq.?2.172?,and apply the initial condition
?2.173?,we obtain
r ˙0
2?y ?r 0?=exp ??2
?R 0
r 0
????d ???y 0
+2
?R 0r 0
????exp ?2
?R 0
?
???ˉ?d ?ˉ?d ?
?
?2.176?
For a prescribed value of the initial radius R 0of the cavity,Eq.
?2.176?yields an expression for the velocity r
˙0of the edge of the cavity and Eq.?2.170?then yields the value of r
¨0.Our last case is identical to the one considered in Sec.2.1,namely,an expanding spherical shock wave.
2.2.4Dimension n=3:Spherical Shock Waves .For a spherical shock wave,the relations ?2.109?,?2.111?,and ?2.112?,with ?=3,assume the form
u ?r ,t ?=??1????r 3?
R 03?
+r 03?t ?
?r
?2.177??u ?r ,t ??t =r 02?t ?r
˙0?t ??3
??1????r ?R 0?+r 0
?t ???v ?r ,t ??2.178?
?2u ?r ,t ??t 2=2r 0?t ??r ˙02?t ??+r 02
?t ?r
¨0?t ??3??1????r 3?R 03?+r 03?t ??
2?
2r 04?t ??r ˙02?t ??
?3
??1?
???r ?R 0?+r 0?t ??
?2.179?
while the position of the shock wave front at time t is given by Eq.
?2.110?with ?=3,i.e.,by
w ?t ?=
?
3
R 03+
1??r 0
3
?t ??R 03??2.180?
The velocity of propagation associated with an expanding spheri-cal shock wave front is then given by
??t ??w
˙?t ?=r 02?t ?r ˙0?t ??
?3
?
R 03+1?
?r 03?t ??R 03??
2
?2.181?
Using Eq.?2.181?,in the kinematic and dynamic compatability
conditions ?2.99?,we ?nd for the velocity and pressure at the shock wave front the following:v u ?t ?=???t ??
r 02?t ?r ˙0?t ??3
?
R 03+1?
?r 03?t ??R 03
??
2
?2.182?
and
p u ?t ??p 0+?0??2
?t ?=p 0+?0
r 04?t ?r ˙02?t ??
?3
?
R 0
3+1?
?r 03?t ??R 03
??
4
?2.183?
Substituting Eqs.?2.177?and ?2.179?into expression ?2.113?for
the pressure ?eld,in the partially saturated soil,we ?nd that p ?r ,t ?=p e ?t ?+
?01??
??1
?
3?1????r 3?R 03?+r 03
?t ?
?
1r 0?t ?
?
??2r 0?t ?r ˙02?t ?+r 02
?t ?r
¨0?t ???
?
?02?1???
?
1
?3
??1????r
3
?
R 03?+r 03?t ??
4
?
1r 04?t ?
?
r 04?t ?r ˙02?t ??2.184?
For the evolution of the radius of the spherical cavity r 0?t ?,we
again deduce an equation of the type ?2.170?,i.e.,
r ¨0?t ?+???r 0?t ??r ˙02?t ?=?
??r 0?t ???2.185?
which,as in the case of a cylindrical shock wave,may be reduced
to one of the form
dy
dr 0
+2???r 0?y =2?
??r 0??2.186?
by setting y =r ˙02?t ?.The initial condition for Eq.?2.186?assumes,
in this case,the form
y 0?y ?r 0?=?
p m ?p 0?0
p m =0.5p H
?2.187?
and a solution of the initial value problem is again obtained in the
form ?2.176?with ???·?,?
??·?replacing ??·?and ??·?.In lieu of ??r 0?and ??r 0?,i.e.,Eqs.?2.174?and ?2.175?,we now have,in Eq.?2.186?,
???r 0?=2r 0?r 0
22
·?
?
R 03+1??r 03
?R 03??
?4/3
?r 0?4?
R 03+1??r 03?R 03
??
?1/3
?r 0
?1?
???1???r 0
2??
R 03+
1??r 0
3
?R 03??
?1/3
?r 0
?1?
?2.188?
and
?
??r 0?=??1???
?0
p m
?R 0r 0
?
3k 0
?p 0
??
R 0
3+1??r 0
3?R 03??
?1/3
?r 0?1?
r 0
2?2.189?
This completes our discussion of the derivation of closed form solutions for the problem of the propagation,in a partially satu-rated soil,of ?pressure ?shock waves which are generated by the detonation of a charge possessing either plane,cylindrical,or
spherical symmetry.The analysis has turned upon the assumption that the mechanical properties of the sediment are approximated by a modi?ed three-component soil model of the type introduced by Lyakhov?1?and Rakhmatulin?3?in the zone where loading occurs,and by the plastic gas model of Rakhmatulin?3?in the zone where unloading occurs.For this problem closed analytical formulas have been obtained which determine both the state and motion of the?three-component?sediment consisting of gas?usu-ally,air?,water,and sand grains.It was noted,in Ref.?7?,that the equations which have been delineated here are suitable for engi-neering applications and that,in particular,the results can be used to estimate the degree and magnitude of the densi?cation zone of water saturated soil subjected to an explosive motion.Within the framework of the plastic gas model we have employed,the degree of densi?cation of the soil is determined by the air content,i.e., knowledge of the air content in water saturated soil enables one to estimate,without solving the full problem of a concentrated ex-plosion,the maximum degree of densi?cation of the sediment.
2.3Transitioning From the High Pressure Plastic Gas Zone to the Moderate Pressure Elastic-Plastic Zone
2.3.1Initial Conditions.At some time t=t M,when the pres-sure at the shock wave decreases to a value p lim,the propagation of the shock wave is assumed to be terminated and the wave is transformed into an ordinary spherical wave,cylindrical wave,or plane wave;as this wave propagates,the soil continues to be compacted but now undergoes,for t?t M,elastic-plastic deforma-tion in the region between r0?t?,the boundary of the cavity and w?t?,the boundary of the wave front.The following“initial”con-ditions then prevail at time t=t M:
u?r,t M?=0,v?r,t M?=v u?r?
?2.190???r,t M?=?M?r?,?rr?r,t M?=????r,t M?=p M?r?
where the functions v u?r?,?M?r?,and p M?r?are determined by the solution of the problem in the plastic gas zone at time t=t M.In Eq.?2.190?,?rr and???denote the non-zero radial and tangential components of the Cauchy stress tensor.In the zone where elastic-plastic deformation occurs,for t?t M,the stress components must be decomposed into their elastic and plastic parts and the results substituted into the equation of motion.
2.3.2Decomposition of the Stress Components and the Equa-tions of Motion.Recall that the coef?cient?indicates whether we are dealing with a plane??=1?wave,cylindrical??=2?wave,or spherical??=3?wave.For t?t M the Cauchy stress tensor?com-ponents?are decomposed as
?rr=?rr E??1EP???=???E??2EP?2.191?
where for?=1,2,or3the elastic parts?rr E and???E are given, respectively,by
?rr E=??+2???u
?r+???1??
u
r
?2.192a?
???E=??u
?r+???1???+??
u
r
?2.192b?
with?and?the Lamécoef?cients of the medium in the elastic range.In particular,for the case of spherical symmetry,
?rr E=??+2???u
?r+2?
u
r
?2.193a?
???E=?
?u
?r+2??+u?
u
r
?2.193b?
The equation of motion in the compacted zone,for t?t M,has the
form
?
?r?
?1+u r???1?rr?+?+1r??1+u r???1?rr
??1+?u?r
??1+u r???2????=?0?2u?t2?2.194?
so that,e.g.,for a spherical wave
?
?r?
?1+u r?2?rr?+2r??1+u r?2?rr??1+?u?r??1+u r?????
=?0
?2u
?t2
?2.195?
All the nonlinear terms which occur in the equation of motion
?2.194?appear,as a consequence of the decomposition?2.191?,
and the constitutive relation?2.192a?and?2.192b?,in the func-
tions?1EP and?2EP;using this observation,we may,upon substi-
tuting Eq.?2.191?in Eq.?2.194?,rewrite the equation of motion
for the general case in the form
?
?r
?rr E+
???1?
r
??rr E????E???0
?2u
?t2=
?
?r
?
1
EP
+
???1?
r??1EP??2EP
?
?2.196?
with the corresponding result for the speci?c case of a propagating
spherical wave obtained by setting?=3in Eq.?2.196?.The func-
tions?1EP and?2EP must be obtained by integrating,from time t M
up to time t,the rate constitutive relations for the radial and tan-
gential components of the elastic-plastic stress;this,in turn,en-
tails making a suitable choice for the plasticity theory which gov-
erns the behavior of partially saturated soils in the moderate
pressure zone.Any suitable elastic-plastic model for the deforma-
tion of the soil,once the pressure generated by an explosion has
decreased to the level p lim,must,of necessity,incorporate the
dilatant behavior of partially saturated soils in an appropriate
manner,i.e.,it must treat dilatation as a kinematical constraint on
allowable motions of the sediment.In Sec.3of this review we
present a brief summary of those results which are currently avail-
able for modeling the dilatancy of soils under explosive motions.
3Spherical Wave Propagation in the Elastic-Plastic
Zone:Models Incorporating Dilatancy
3.1Work in the Russian Literature
3.1.1Elastic-Plastic Dilatancy Models.The term“dilatancy,”
and the experimental discovery of that phenomena in relation to
the behavior of water-saturated sands,is credited to Reynolds?13?
in1885.The term dilatancy refers to the effect whereby shearing
?of a granular material?leads to volume changes.Dilatancy of
geomaterials is a characteristic feature of their mechanical behav-
ior which appears even during the initial loading of such a me-
dium.The purely kinematical nature of dilatancy is directly re-
lated to the fact that the distortion of a set of particles in contact
leads to essential changes in the void space between particles.As
noted in Ref.?14?,“dilatant proportionality of shear and volume
increments following a nonassociated?ow rule have shown ex-
tremely good agreement with granular material response.”In fact,
most experimentalists,who work in the general area of geomate-
rial response,point to the effect of dilatancy as playing the major
role in earthquake damage.
In Ref.?14?,which is representative of much of the work of the
Russian work in this area,a model of elastic-plastic soil response
has been constructed in accordance with guidelines provided by sets of triaxial experimental data.The basic model of dilantant plasticity for soils that is formulated in Ref.?14?follows earlier work of Nikolaevsky ?15?and Nikolaevsky,Syrnikov,and Shefter ?16?and is a simple generalization of the Prandtl-Reuss plasticity
theory in the sense that the kinematical condition e ˙ij p ?ij =0,ex-pressing plastic incompressibility,is replaced by the dilatancy
constraint
?e ?e ˙ij p ??e ˙ij p ?ij ????ij ,????
˙p ?=0?3.1?
with
?
˙p =1
?2
??e ˙11p ?e ˙22p ?+6e ˙12p e ˙12p
+...?1/2?3.2?
where ?is the dilatancy rate,?ij the Cauchy stress tensor,?the hardening parameter,e ˙ij p the plastic strain rate,and ?
˙p the plastic shear rate.If ?
˙p =0,then Eq.?3.1?reduces to the usual constraint of plastic incompressibility.For geomaterials,the Coulomb-Mohr limit condition is employed as the yield criteria,namely,
????ij ????????p ?Y =0
?3.3?
with
?=
1
2?2
???11??22?2
+6?12?12+...?1/2?3.4?
where ?is the shear stress,p the pressure,?the internal friction coef?cient,and Y the cohesion.For the linear isotropic case we have,in addition,the nonassociated ?ow rule
e
˙ij p =Q ijkl p ?kl =??ij +p ?ij ????H +p ??ij ?˙?3.5?
which involves the two functions ?and ?;in Eq.?3.5?,Y =?H .In
the elastic domain,the usual law applies,i.e.,
e ij
e =e ij ?e ij p =Q ijkl p ?kl ?3.6?
with Q ijkl e the elastic modulus tensor.For the hardening parameter ?,in Eq.?3.1?,one may,e.g.,use the ratio ?=e p /e *p
where e p
=e ij p ?ij and e *p is a critical value of e p for which ??e p ?=0.For river sands we display in Figs.19and 20,respectively,graphs of the internal friction coef?cient ?as a function of the dilatancy rate ?
and the hardening parameter ?=e p /e *
p
?which may also be shown to be related to a ratio of porosity to critical porosity ?.Some authors have also used ?=?p as the hardening parameter.3.1.2Dilatancy Under Explosive Motions .The complete dila-tant plastic model considered in Ref.?14?was developed ?Ni-kolaevsky ?15?,Rodionev et al.?17??for use in the study of un-derground explosions.Based on that model we show,in Fig.21,the explosive wave pressures at three times;the wave pro?les exhibit two local maxima,the ?rst of which is interpreted as the wave front while the second one is interpreted as indicating a maximum for the residual stresses around the spherical cavity containing the explosive charge.The main work done by the ex-plosive charge consists of compaction ?up to the critical level at which p =p lim ?after which dilatant loosening of the soil begins;the latter phenomena corresponds to a change in the sign of ?.In Fig.22,we depict the post-explosion density ?and the sound velocity in the radial and tangential directions for sands with an initial porosity of 25%.The external boundary for dilatant plastic ?ow turns out to be in good agreement with the external radius for the actual ?experimentally determined ?dilatant compaction.Ex-plosions in soft soils yielded sublimit plastic strains up to the nondimensional radius R =r ??c 2/W 1/3?=50,W being the explosive yield of the charge,with the conclusion that permeability changes in highly porous geomaterials are connected with small residual deformations of the pore space.On the other hand,the model predicts that explosions in saturated porous media yield nonmono-tonic changes in permeability.3.2
Work in the European Literature
3.2.1Dilatation as a Kinematical Constraint .The work of Houlsby ?18?appears to be typical of much of the European work in which the phenomena of dilatancy is introduced into thermo-mechanical plasticity models for the deformation of soils by treat-ing dilatancy as a kinematical constraint,i.e.,it is shown that the introduction of appropriate constraints into the plasticity model results in the phenomena of dilatancy.In Ref.?18?,relationships are also obtained between the angle of dilation associated with the soil and the angle of friction;furthermore,
by extending these
Fig.19Internal friction coef?cient ?as a function of dilatancy rate ?for river sands …adopted from
Ref.?14?…
Fig.20Internal friction coef?cient ?as a function of harden-ing parameter ?for compaction of loose river sands …adopted from
Ref.?14?…
Fig.21Explosion wave pressures for three moments of time …adopted from Ref.?14?…
concepts the author ?16?is able to introduce a simple measure of anisotropy into the model by using a further constraint equation and to subsequently link the dilation of the soil to the anisotropy.The ?nal model in Ref.?16?is capable of reproducing not only dilatancy and anisotropy but also embodies a critical state and exhibits the property of densi?cation under small strain cycles.3.2.2Thermomechanical Formulation .In Ref.?18?Houlsby has described a method whereby plasticity theories can be derived from a formulation,based on thermodynamic principles,which was ?rst developed by Ziegler ?19?;the basic approach focuses on the importance of the kinematic variables,i.e.,the strains e ij ,and a set of internal variables ?ij ,and involves a free energy function F in the form
F =F ?e ij ,?kl ?
?3.7?
The rate of dissipation per unit volume,D ,is assumed to be a
positive function both of the kinematic variables and of their rates
D =D ?e ij ,?kl ,e
˙mn ,?˙pq ??3.8?
while the stresses ?ij may be derived from the relation
?ij =?F
?e ij +?
?D ?e
˙ij ?3.9?
The internal variables may be eliminated by using
?F
??y
˙+??D ??˙ij =0
?3.10?
The quantity ?in Eqs.?3.9?and ?3.10?is a scalar multiplier;if D
is homogeneous,and of order n in the ?ij ,then ?=1/n .If n =1it may easily be shown that the stresses are rate independent;all the cases considered here will involve rate independent theories so that n =1,?=1.As an example,isotropic linear elasticity may be derived from the expressions
F =1
2Ke ii e jj +Ge y ?e ij
?D =0?3.11?
where e ij ?denotes the deviator of the strain tensor and K and G are the bulk and shear moduli;in this example,there are no internal variables and Eq.?3.10?is not required.
An elastic-plastic theory,with a von Mises yield surface,and an associated ?ow rule,can be derived by introducing the plastic
strains e ij p
as internal variables and de?ning
F =1
2Ke ii e jj +G ?e ij ??e ij p ??e ij ??e ij p
?
?3.12?D =
?
83
·c ·?e ij p e ij p ?3.13?
It is shown in Ref.?18?that Eqs.?3.12?and ?3.13?yield the following results:
?i ??ii =3Ke ii ?the volumetric response is elastic ?
?ii ?e ij ?=1/2G ?ij ?+e ij
p and,either ?iii ?e ˙ij p =0?elastic deviator behavior ?or
e
˙ij p
0?plastic strain occurs ?
?ij
??ij ?=83
c 2?the stresses lie on a von Mises yiel
d surfac
e ?and
?ij
?=?8/3ce ˙ij p ?e ˙kl p e
˙kl p ?plastic strain rates satisfy the von Mises
flow rule ?
3.2.3Some Simple Models Incorporating Dilatancy .The phe-nomena of dilatancy can be introduced into the formulation delin-eated above by making use of kinematic constraints as described by Ziegler ?19?.Either the model may be formulated in terms of all the kinematic variables,with the constraint equations then in-troduced,or the model may be expressed in terms of a reduced number of kinematic variables which are unconstrained;it is noted in Ref.?18?that the ?rst of these approaches is a more general and powerful technique than the second.In order to intro-duce constraints,relations ?3.8?and ?3.9?are supplemented by one or more constraint equations of the form
C ?i ?=C ?i ??e ij ,?kl ,e
˙mn ,?˙pq ?=0?3.14?
with each of the constraint functions C ?i ?a homogeneous function
of order one in the rates of strain and/or the rates associated with the internal variables.We then modify Eqs.?3.9?and ?3.10?so as to read
?ij =
?F
?e ij +??D ?e ˙ij +??i ?
?C ?i ??e
˙ij ?3.15?
and
?F
??ij +??D ??˙ij +??i ??C ?i ???
˙ij =0
?3.16?
where ??i ?is a multiplier which can be eliminated by use of the constraint equation.
We begin by writing down a simple frictional model using the following variables ?in terms of principal strains and stresses ?which are appropriate for a triaxial test:p =1
3??1+2?3?
q =?1??3?=?1+2?3
and ?=2
3??1??3?
?
3.17?
The frictional model follows from taking
Fig.22Density,sound velocities,and permeability near ex-plosion cavity in a dry high porous medium …adopted from Ref.
?14?…