Hollow Sphere Subject to an Internal Blast
Note
To view this project in FLAC3D, use the menu command . Choose “Dynamic/SphericalPulse” and select “SphericalPulse.prj” to load. The main data files used are shown at the end of this example. The remaining data files can be found in the project.
This problem concerns the propagation of a spherical wave due to an impulsive pressure (explosion) in a sphere. In unbounded (i.e., infinite) media, two types of waves can exist: compression and shear waves. In this problem, the axisymmetric nature of the problem eliminates the shear wave. Therefore, only the solution for the compression wave need be sought. The problem provides a test of the dynamic capabilities of FLAC3D and is applicable to impact and explosion modeling.
The analytical solution for this problem, assuming that the material is elastic and isotropic, was derived by Blake (1952). The solution is based on the following governing equation:
where: 
\(C_p\) 
= compressional wave velocity; 
\(t\) 
= time; 

\(\phi\) 
= a potential function; and 

\(\bigtriangledown^2\) 
= Laplacian operator. 
Let \(p(t)\) be an impulse that jumps from zero to \(p_0\) at \(t\) = 0, and then decays exponentially with the time constant \(\alpha^{1}\). Thus, the pressure function can be defined by
A step function of the pressure (\(\alpha\) = 0) will be used for this problem. For such a pressure function, the potential function that satisfies the governing equation is
where: 
\(a\) 
= radius of the sphere; 
\(K\) 
= \({{1  \nu} \over {2 (1  2 \nu)}}\); 

\(\nu\) 
= Poisson’s ratio; 

\(r\) 
= radial coordinate; 

\(\alpha_0\) 
= \({{C_p} \over {2 a K}}\) = radiation damping constant; 

\(\tau\) 
= \(t  {{ r  a} \over {C_p}}\); and 

\(\omega_0\) 
= \({{c} \over {2 a K}}\ \sqrt{4 K  1}\) = natural frequency. 
The radial displacement can be found by differentiating the potential function with respect to radial distance.
A sphere embedded in an infinite, isotropic medium is simulated in two ways. Figures 1 and 2 show two different grids used for the simulation. The first grid is a onedimensional representation in which a narrow, tapering tube with roller boundaries on all sides approximates the condition of spherical symmetry. The second grid is fully threedimensional, but it is bounded by three symmetry planes (i.e., oneeighth of the problem is modeled). The radius of the inner boundary is assumed to be 10 m, and the outer boundary is located at a distance ten times that of the inner boundary for the 1D case.
For both grids, normal movement is prevented at the axes of symmetry, and a viscous condition is imposed on the outer boundary to absorb the wave. The material properties used for the problem are
shear modulus (\(G\)) 
1 × 10^{10} 
bulk modulus (\(K\)) 
1.665 × 10^{10} Pa 
density (\(\rho\)) 
1675 kg/m^{3} 
A pressure equal to 1000 Pa is applied at the inner boundary to simulate the blasting.
The radial displacement histories recorded up to 0.1 second at \(r\) = 2.05\(a\), 3.42\(a\), and 4.87\(a\) are given in Figures 3 and 4 for the onedimensional and threedimensional simulation, respectively. The delay of the response at locations far from the sphere can be noted in both cases. In both cases, FLAC3D is able to capture the response at peak and steady states. The fluctuation at late time is due to the fact that the radiated wave is not absorbed completely by the viscous boundary.
References
Blake, F. G. “Spherical Wave Propagation in Solid Media,” J. Acous. Soc. Am., 24(2), 211215 (1952).
Data Files
SphericalPulse1D.dat
model new
model largestrain off
model title "spherical cavity  pressure pulse  1D model"
model configure dynamic
; Generate zones, representing a 1D column on the sphere.
zone create brick point 0 (10,0,0) point 1 (100,0,0) ...
point 3 (9.8769,1.5643,0) point 6 (98.770,15.645,0) ...
point 2 (9.8769,0,1.5643) point 4 (98.770,0,15.645) ...
point 5 (9.7552,1.5643,1.5451) ...
point 7 (97.555,15.645,15.450) ...
size 30 1 1 ratio 1.01 1 1
; Assign constitutive model and properties
zone cmodel assign elastic
zone property shear=1e10 bulk=1.665e10 density=1675
; Name the boundary faces
zone face skin
; Apply boundary conditions
zone face apply velocitynormal 0 range group 'South' or 'North' or 'Top' or 'Bottom'
zone face apply stressnormal 1000 range group 'West'
zone face apply quiet range group 'East'
; Set damping
zone dynamic damping rayleigh 0.01 50 stiffness
;Take some histories
history interval 5
model history name='time' dynamic timetotal
zone history name='xd1' displacementx position (20.5,0,0)
zone history name='xd2' displacementx position (34.2,0,0)
zone history name='xd3' displacementx position (48.7,0,0)
; Run model for 0.1 seconds, then save
model solve timetotal=0.1
model save 'SphericalPulse1D'
SphericalPulse3D.dat
model new
model largestrain off
fish automaticcreate off
model title 'spherical cavity  pressure pulse  3D model'
model configure dynamic
; Create zones, in a 1/8 brick
zone create radialbrick edge 100 size 16,16,16,30 ...
rat 1 1 1 1.01 dimension 10,10,10
; A FISH function is used to convert from a cube to a sphere
fish define make_sphere
; Loop over all GPs and remap their coordinates:
loop foreach local p_gp gp.list
local p = gp.pos(p_gp)
local dist = math.mag(p)
local k = 10.0/dist
local a = p * k
local maxp = math.max(p>x,math.max(p>y,p>z))
k = 100.0/maxp
local b = p * k
local u = (maxp10.0)/90.0
gp.pos(p_gp) = a + (ba)*u
end_loop
end
[make_sphere]
; Assign constitutive model and properties
zone cmodel assign elastic
zone property shear=1e10 bulk=1.665e10 density=1675
; Name boundary faces
zone face skin
; Assign boundary conditions
zone face apply velocitynormal 0 range group 'West' or 'South' or 'Bottom2'
zone face apply quiet range group 'East' or 'North' or 'Top'
zone face apply stressnormal 1000 range group 'Bottom1
; Assign damping
zone dynamic damping rayleigh 0.01 50 stiffness
; Take some histories
history interval 5
model history name='time' dynamic timetotal
zone history name='xd1' displacementx position (20.5,0,0)
zone history name='xd2' displacementx position (34.2,0,0)
zone history name='xd3' displacementx position (48.7,0,0)
; Run model for 0.1 seconds, then save
model solve timetotal=0.1
model save 'SphericalPulse3D'
Was this helpful? ...  Itasca Software © 2022, Itasca  Updated: Mar 09, 2023 