## Abstract

This paper investigates the robustness against localized impacts of elastic spherical shells pre-loaded under uniform external pressure. We subjected a pre-loaded spherical shell that is clamped at its equator to axisymmetric blast-like impacts applied to its polar region. The resulting axisymmetric dynamic response is computed for increasing amplitudes of the blast. Both perfect shells and shells with axisymmetric geometric imperfections are analyzed. The impact energy threshold causing buckling is identified and compared with the energy barrier that exists between the buckled and unbuckled static equilibrium states of the energy landscape associated with the pre-loaded pressure. The extent to which the impact energy of the threshold blast exceeds the energy barrier depends on the details of its shape and width. Targeted blasts that approximately replicate the size and shape of the energy barrier buckling mode defined in the paper have an energy threshold that is only modestly larger than the energy barrier. An extensive study is carried out for more realistic Gaussian-shaped blasts revealing that the buckling threshold energy for these blasts is typically in the range of at least 10–40% above the energy barrier, depending on the pressure pre-load and the blast width. The energy discrepancy between the buckling threshold and energy barrier is due to elastic waves spreading outward from the impact and dissipation associated with the numerical integration scheme. Buckling is confined to the vicinity of the pole such that, if the shell is not shallow, the buckling thresholds are not strongly dependent on the location of the clamping boundary, as illustrated for a shell clamped halfway between the pole and the equator.

## 1 Introduction

### 1.1 Static Behavior and Energy Barriers.

Figure 1 presents the energy barrier against axisymmetric buckling for a perfect elastic spherical shell subject to uniform external pressure for three sets of boundary conditions: a full sphere buckling symmetrically about its equator, a hemisphere clamped at its equator, and a deep spherical cap clamped at $45deg$ from the equator. In Fig. 1 and throughout this paper, elastic spherical shells are considered with radius R, thickness h, Young’s modulus E, and Poisson’s ratio ν. The reference values
$pC=2Eh23(1−ν2)R2andΔVC=4π(1−ν)R2h3(1−ν2)$
(1)
are the pressure and the volume decrease of the full perfect spherical shell at bifurcation, the so-called classical buckling condition. In Fig. 1, the normalized buckling energy barrier, Eb/[(pCΔVC/2)Ch/R], with $C=3/[(1−ν)1−ν2]$, is plotted as a function of the normalized pressure applied to the shell, p/pC, on the vertical axis. Each of the curves in Fig. 1 has been computed assuming axisymmetric behavior about the axis through the North–South poles of the sphere. For thin shells (e.g., R/h ≥ 50), the post-buckling mode is confined to the vicinity of the pole, as will be seen in results to follow, and has almost no interaction with the buckle on the opposite pole in the case of the full sphere or with the clamped boundary in the other cases. For thin shells, the dimensionless energy barrier plots in Fig. 1 are essentially independent of R/h and ν as well as the clamping condition for p/pC ≤ 0.85. The factor (pCΔVC/2)Ch/R used to normalize the energy barrier is the product of the elastic energy in the full spherical shell at the classical bifurcation pressure, pCΔVC/2, and the small term Ch/R proportional to the thickness to radius ratio of the shell. Thus, it is evident from Fig. 1 that the energy barrier is a small fraction of the total elastic energy stored in the shell except possibly at very low applied pressures. It is also evident in Fig. 1 that there is negligible difference between the energy barrier (per dimple buckle) for the clamped shells and the full shell in the range p/pC ≤ 0.85.

Based on Huang’s [7] analysis of the elastic buckling of deep spherical caps, non-axisymmetric quasi-static buckled states of the perfect spherical shell clamped at the equator are expected to generate lower saddle energies than those for axisymmetric states for applied pressures above about p/pC ≅ 0.85. However, for pressures below about p/pC ≅ 0.85, there is strong evidence that the quasi-static saddle buckled states having the lowest energy are associated with axisymmetric dimple-like modes, as will be discussed in the sequel. Because the study in this paper is limited to axisymmetric behavior, we will primarily focus our attention for the perfect shell on the range of applied pressures 0.25 ≤ p/pC ≤ 0.85.

The axisymmetric quasi-static pre-buckling and post-buckling normal deflections of the perfect shell clamped at the equator are presented in Fig. 2 for the full range of pressures and paired with the plot of pole deflection. A plot of pressure versus change in volume is also included. In the two plots on the left in Fig. 2, the upper part of each curve (w/h > −0.5) is the stable response characterizing the unbuckled state (or “node” state in the terminology of nonlinear dynamics), while the lower part of the curve (w/h < −0.5) is the unstable post-buckled state (saddle point). The static buckling pressure of the shell clamped at the equator is slightly reduced to pfold = 0.956pC compared with the perfect full spherical shell. Static buckling (assuming a rotationally symmetric response) occurs at the point p = pfold, with wpole/h ≈ −0.5, connecting the stable and unstable responses (see Fig. 2 middle panel, maximum of p/pC). The left panel shows the radial displacement distribution, w(θ)/h, for all equilibria. For saddles, the deformations are strongly localized for pressures at some distance from pfold (below about p = 0.85pC). The deformation accommodating the clamping is concentrated near the equator and is separate from the buckling dimple for saddles at pressures p away from pfold, explaining why the energy barrier for buckling in Fig. 1 is insensitive to the clamping condition in this range.

The right panel of Fig. 2 shows the stable and unstable branches in the (ΔV, p)-plane. In this plane, the fold is very sharp (at the top right) such that stable and unstable branches are close together. Showing the static response in this plane illustrates that the energy barrier Eb between the stable and unstable equilibrium states is very small: the energy difference Eb for a particular p equals the area under the curve above p, reflecting the trend seen in Fig. 1. Figure 3 shows the spatial profile of a stable and unstable equilibrium at p = 0.6pC again.

### 1.2 Outline of the Paper.

The energy barrier plays a central role in providing an understanding of the robustness of the pre-loaded shell against disturbances such as local impacts, and it enables a quantitative means of rationalizing the energies of disturbances required to buckle the shell. To back up this statement, an extensive study has been conducted in this paper of perfect and imperfect elastic spherical shells that are clamped at the equator and subject to blast-like impacts of various amplitudes and shapes in the vicinity of the pole. The quasi-static saddle buckling mode at pressures below pfold has the form of a dimple with axial symmetry about its center. Except at applied pressures just below pfold, published simulations and experimental observations indicate that the buckling behavior is dominantly axially symmetry about the center of the emerging dimple buckle. Non-axisymmetric features only develop deep into the post-buckling response, well beyond the range of relevance to the considerations in this paper [2]. For these reasons, an axisymmetric dynamic analysis of the hemispherical shell captures the essential aspects of the dynamic buckling process. For the imperfect shells investigated here, the initial geometric imperfection is also assumed to have axial symmetry. The pre-loaded shell is initially at rest and then subject to a suddenly imposed localized initial velocity distribution that simulates the blast-like impact. Two basic axisymmetric disturbances are considered: one in which the shape targets the initial velocity distribution directly toward the buckled state associated with the saddle equilibrium state, and the other in which the normal velocity component has a Gaussian shape centered at the pole. The amplitude of the blast impact, as measured by the initial kinetic energy imparted to the shell, is systematically increased until it attains a threshold large enough to buckle the shell. Plots of the threshold blast energy for buckling normalized by the energy barrier as a function of the applied pressure provide a clear picture of the influence of width and shape of this type of disturbance.

In brief, the outline of the paper is as follows. The shell model is introduced in the remainder of this section. The protocol of the numerical experiments for studying the blast load is presented in Sec. 2 with further details discussed in context later. The method is used to produce results for blasts tailored in shape to replicate the quasi-static saddle mode in Sec. 3, thereby targeting the saddle in the energy landscape and giving rise to buckling from disturbance energies only slightly above the energy barrier. Trends for the more realistic Gaussian-shaped blasts are presented in Sec. 4 where it is also demonstrated that the blast energy threshold for buckling depends only weakly on the meridional clamping angle, assuming the shell is not a shallow cap. Selected trends for the blast threshold are presented in Sec. 5 for hemispherical shells clamped at the equator and having a Gaussian-shaped initial geometric imperfection at the pole. Additional implications of the study are discussed in Sec. 6.

### 1.3 The Shell Model.

As previously noted, the model assumes a response that is rotationally symmetric around the vertical North–South pole axis with θ as the meridional angle with θ = 0 at the equator and θ = π/2 at the pole. The model is based on the small strain-moderate rotation theory of Sanders [9] and Koiter [10]. If the shell is clamped at the equator, then the lower limit is θ0 = 0 while θ0 = π/4 for 45 deg clamping. The upper limit θ1 in the numerical integration process is slightly below π/2 because the coordinates are singular at the pole (π/2 − θ1 ≈ 1.68 × 10−3π in all simulations). The model corresponds to the situation where a small rigid disc is inserted and “welded” at the North pole. Throughout all following sections, we will typically display all quantities in a non-dimensionalized form. The applied pressure p is measured in multiples of the critical (buckling) pressure of the perfect fully symmetric shell, pC. Dimensionless time is τ = t/T0, where T0 is the period of volume oscillations of the sinusoidal spherically symmetric vibration mode of the perfect shell
$T0=2(1−ν)ρEπR$
(2)
with ρ as the density of the material.
The radial and tangential displacement components of the shell, w(θ, t) and u(θ, t), are functions of θ and t, and the initial stress-free displacement of middle surface when the shell is imperfect is wimp(θ). Following the non-dimensionalization used in Refs. [11,12], the dimensionless counterparts are (W, U, Wimp) = (w, u, wimp)/R. Derivatives of the dependent variables with respect to angle θ are denoted by a prime (e.g., W′) and with respect to τ by a dot (e.g., $W˙$). The non-zero stretching strains of the shell middle surface are $(εθ,εω)$. The non-zero curvature changes of the middle surface are $(Kθ,Kω)$ with dimensionless bending strains defined as $(κθ,κω)=R(Kθ,Kω)$. The dimensionless strain–displacement relations are
$εθ=W+U′+12φ2−φWimp′,εω=W−Utanθκθ=φ′,κω=−φtanθ$
(3)
with φ = −W′ + U. The non-zero resultant in-plane stresses, $(Nθ,Nω)$, and bending stresses, $(Mθ,Mω)$, have dimensionless components $(nθ,nω)=(Nθ,Nω)/Eh$ and $(mθ,mω)=(Mθ,Mω)R/D$, with D = Eh3[12(1 − ν2)] as the bending stiffness. The dimensionless stress–strain relations are
$(nθ,nω)=(εθ+νεω,εω+νεθ)/(1−ν2),(mθ,mω)=(κθ+νκω,κω+νκθ)$
(4)
In dimensionless form, expressions for the volume decrease, the potential energy of the time-independent applied uniform external pressure, the kinetic energy, and the strain energy of the shell are, respectively
$ΔV(τ)2πR3=−∫θ0π/2Wcosθdθ,PE(τ)2πR3=p¯∫θ0π/2WcosθdθKE(τ)2πD=12∫θ0π/2[W˙2+U˙2]cosθdθSE(τ)2πD=12∫θ0π/2[(κθ2+κω2+2νκθκω)+12(R/h)2(εθ2+εω2+2νεθεω)]cosθdθ$
where $p¯=pR3/D$ is the dimensionless pressure, which is positive pointing inward. In evaluating the integrals, the integrands are extrapolated from values of θ greater than θ1 to the pole. The inertia of any fluid medium inside and outside the shell is neglected. For the problems considered in this paper, apart from dissipation arising from the numerical method discussed in the sequel, the formulations are conservative such that SE + PE + KE does not change once the shell is set in motion. The energy barrier at a fixed pressure p is defined in terms of the static equilibrium solutions at the node and saddle to be Eb = (SE + PE)s − (SE + PE)n with the subscripts s denoting saddle and n denoting node.

## 2 Spatial Shape of Initial Conditions and Buckling Threshold Trajectories

The numerical experiment subjects the clamped shell which is initially at rest and subject to a uniform pressure p to a locally applied “blast,” represented by an initial axisymmetric velocity distribution. The initial deformation will always be taken to be the stable equilibrium state at pressure p (the node), as shown in Fig. 2. We consider two blast shapes: targeted shapes and Gaussian shapes. By shape, we are referring to the spatial distribution of the velocity imposed on the shell at t = 0. The targeted shapes are selected to buckle the shell with as small initial kinetic energy as possible—these involve initial velocity components that are unlikely to be produced by a real blast. The Gaussian shapes are expected to be more realistic. With $W˙ini(θ)$ and $U˙ini(θ)$ as the dimensionless velocities imposed at t = 0, the initial kinetic energy associated with the blast is
$KEini2πD=12∫θ0π/2[W˙ini2+U˙ini2]cosθdθ$
(5)
The work done by the instantaneous impact is KEini. We have not investigated impacts applied over a short but finite time period Tini for which one would again want to compare the work done by the impact with the energy barrier. This would require a more detailed analysis of the behavior while the impact is acting. The dynamic responses for the instantaneous impacts should approximate the finite time impacts if Tini/T0 is sufficiently small, but the precise requirement has not been established here.

### 2.1 Blasts Targeted Toward the Saddle Post-Buckling State.

The targeted blast will illustrate that indeed it is possible to impact the shell and cause it to buckle with a disturbance energy that is just slightly above the energy barrier. The targeted blast will also allow us to judge how far from “optimal” the Gaussian-shaped blasts are.

Denote the displacement components in the stable and unstable (i.e., the node and the saddle) static equilibrium states at pressure p by (Wn, Un) and (Ws, Us), respectively. Then, define the targeted initial velocities and the initial conditions on the displacements by
$W˙ini(θ)=a[Ws(θ)−Wn(θ)],U˙ini(θ)=a[Us(θ)−Un(θ)]Wini(θ)=Wn(θ),Uini(θ)=Un(θ)$
(6)
where a is the amplitude of the initial velocity. This choice directs the initial velocity in proportion to the difference between the saddle and node states. Because the unstable buckled state is localized to the pole, Ws(θ) − Wn(θ) and Us(θ) − Un(θ) vanish outside the polar region, giving rise to an initial velocity distribution that is localized around the pole with the same width as the saddle buckling mode. It is unlikely that such targeted blasts will arise from an actual blast, in part because the initial radial velocity points outward for some angles near the pole (see the saddle mode shapes in Figs. 2 and 3) and because it contains a non-zero tangential component. The amplitude factor a is increased systematically until buckling is attained, as discussed below, and it is directly related to KEini by Eq. (5).

### 2.2 Gaussian-Shaped Blasts.

This set of blasts has an initial inward normal velocity component and no initial tangential component. Let β = π/2 − θ be the angular distance measured from the pole, and let
$βref=1.51−ν2R/h$
(7)
be a reference angle. For Gaussian-shaped geometric imperfections, this reference angle lies in the range of imperfection widths which are most deleterious for buckling, capturing the dependence on R/h and ν [11]. Take the initial velocity distribution to be
$W˙ini(θ)=−aexp(−(β/βG)2),U˙ini(θ)=0$
(8)
with the displacements at t = 0 given by Wini(θ) = Wn(θ) and Uini(θ) = Un(θ). In the simulations, the width βG of the blasts will be chosen to lie between 0.2βref and 1.2βref, and thus the blasts will be localized to regions of various widths centered at the pole. The velocity amplitude a is again tied to the kinetic energy imparted by the blast by Eq. (5).

### 2.3 Buckling Threshold.

Starting from the stable equilibrium position at fixed background pressure p and a given initial kinetic energy KEini, we numerically integrate the moderate rotation equations for some maximal time tE = 4T0 (recall that T0 is the period of sinusoidal volume oscillations of the spherically symmetric vibration mode of the perfect shell) but check if at any time before 4T0 the deflection at the pole Wpole(t) exceeds the deflection Ws(0) of the saddle by a factor of 1.5, in which case we classify the kinetic energy KEini as “leading to buckling.” Our simulations reveal that trajectories reaching this level of pole deflection will subsequently result in large buckling deformations characteristic of a snap buckling process. Otherwise, we classify the kinetic energy KEini as “not leading to buckling.” By systematically altering the initial velocity amplitude a, we perform a simple bracketing iteration until we find a sufficiently tight bracket [(KEini)low, (KEini)up] for the buckling threshold. Then, we consider (KEini)low the buckling threshold, labeling it KEC or the critical energy. The trajectory associated with this is called the critical trajectory or trajectory leading to buckling (short buckling trajectory).

### 2.4 Systematic Errors From Initial Numerical Dissipation in Buckling Thresholds.

Figure 3 illustrates the geometry of Gaussian blasts in the considered range of widths, compared with the deflection, w(θ)/h, of the static unbuckled (node) and buckled (saddle) states for p/pC = 0.6, shown as solid curves. The other two solid curves represent the “widest” Gaussian-shaped impact (1.2βref) and the most “narrow” one (0.2βref). For each of these Gaussian blasts, Fig. 3 also shows a dashed curve. This dashed curve was obtained by numerically time-stepping 5 time steps forward and 5 time steps backward from the initial blast. We make this distinction because our numerical time-stepping scheme causes a short sharp drop in energy during the first few time steps. This effect is caused by the finite step size Δτ = T0/32 because in the continuous-time system energy should be conserved. Accounting for this effect, we determine an effective initial velocity, which has a slightly lower kinetic energy and is shown as the dashed graph in Fig. 3. Any feasible time-stepping scheme will have similar behavior to ours in this respect. Consequently, some results will show two curves for buckling thresholds: one if the curve assumes that the full kinetic energy of the initial velocity shape is transferred. The other curve subtracts the amount of energy (SE + PE + KE) lost during five initial forward and backward time steps, resulting in an effective initial kinetic energy. See Figs. 4 and 5 below. The effect of damping on the threshold reported can be reduced by reducing the time step size (which in our study is T0/32 and T0/16). The thresholds in Figs. 46 are already close to the theoretical energy barrier, illustrating the main message of this paper, such that the effect of damping is already relatively small.

## 3 Buckling Thresholds for Targeted Blasts

The fact that the threshold targeted blast energy needed to buckle the shell, KEC, is only slightly above the energy barrier Eb is seen in Fig. 4. The buckling threshold energy is slightly larger than the theoretical energy barrier Eb, but dependent on p/pC. Since the theoretical energy barrier Eb is very small for p > 0.85pC, the ratio between KEC and Eb is a ratio between two small numbers. Thus, inaccuracies such as numerical dissipation and numerical tolerances create some uncertainty for p > 0.85pC. The ratio increases at lower pressures. Because the distance between the stable node and the unstable saddle equilibrium is larger at lower pressures, a blast in the linear direction from stable to unstable equilibrium may no longer be the blast leading to buckling with the minimal energy. Nevertheless, these results demonstrate that dynamic disturbances with energies only slightly above the energy barrier can induce buckling. In the discussion of the Gaussian blasts in Sec. 4, further insights into the dynamic responses including details of energy lost due to elastic waves spreading out from the impacted region will be provided. The results for the shell clamped at 45 deg will be discussed later.

## 4 Buckling Thresholds for Gaussian-Shaped Blasts

In this section, Gaussian-shaped impacts having a range of widths defined by Eqs. (7) and (8) will be imposed on the shell. We begin by reporting the results for the buckling thresholds for the perfect clamped hemispherical shell in Figs. 5 and 6 and follow up by presenting details of the dynamic response of the shell and a discussion of the level of error expected in computing the buckling thresholds.

The effect of the blast width on the threshold buckling energy for any fixed pressure in the range 0.25 ≤ p/pC ≤ 0.85 is displayed in the contour plots in Fig. 5. The plot on the left in Fig. 5 shows the contours of constant KEC/[(pCΔVC/2)Ch/R] with p/pC on the vertical axis and βG/βref on the horizontal axis, while the plot on the right has contours of constant KEC/Eb, with Eb given in Fig. 1. In the plot on the left, two sets of contours are shown: one which subtracts off the initial dissipation from the kinetic energy (solid line curves), as discussed earlier, and the other which does not (dashed line curves). The contours plotted on the right are for the lower estimate of the threshold kinetic energy computed by subtracting off the initial dissipation. The width producing the smallest blast threshold for buckling has βG ≅ 0.6βref over the entire range of pressures. The minimum value of KEC/Eb with respect to the width for fixed p/pC has been extracted from the contour plots and plotted in Fig. 6 for both sets of computations, accounting for initial dissipation and disregarding it. As the illustration in Fig. 2 had suggested, this initial dissipation is only strong for sharply localized blasts, e.g., βG/βref small. The blasts widths which generate the minimum buckling thresholds are not strongly affected by the numerical dissipation.

### 4.1 Approximate Dependence of Threshold Kinetic Energy on Blast Shape.

We observe that there is a distinct intermediate optimal blast shape for buckling, where the blast is neither too localized nor too wide. This optimum is shown as a black almost vertical curve in Fig. 5 at roughly 0.6βref. For larger pressures (p > 0.5pC), we can formulate a simple criterion for determining the optimal blast shape: what proportion of the blast is directed toward the saddle? This would correspond to measuring how similar the spatial shape of the initial velocity is to the difference between saddle and node shapes at the given pressure. A natural measure of similarity between two initial velocity shapes would be the scalar product based on the kinetic energy, scaled by the energy barrier
$⟨(W˙1,U˙1),(W˙2,U˙2)⟩KE=2πD2Eb∫θ0π/2[W˙1(θ)W˙2(θ)+U˙1(θ)U˙2(θ)]cosθdθ$
(9)
In this scalar product, the squared norm (“length”) of an initial velocity shape, i.e., $⟨(W˙ini,U˙ini),(W˙ini,U˙ini)⟩KE$, is its (scaled) kinetic energy, for example, KEini/Eb. The similarity of an initial velocity shape $(W˙ini,U˙ini)$ with a given amount of kinetic energy with the directions targeted toward the saddle is measured by the overlap
$OV:=⟨vini,vn↔s⟩KE,wherevini=(W˙ini,U˙ini),vs↔n=(Ws−Wn,Us−Un)/T0$
(10)
The dashed line curve in Fig. 5 shows where this overlap is maximal. We observe that this approximation is acceptably accurate (close to the true minimal curve dotted in Fig. 5) for p > 0.5pC. Conversely, the threshold initial velocity for a shape that has low overlap is expected to be associated with higher threshold KEC. In fact, it is reasonable to expect that for fixed background pressure the threshold kinetic energy and its associated initial velocity field, νini, satisfies approximately
$KECEb∼EbOV=Eb⟨vini,vn↔s⟩KE$
This means that the threshold kinetic energy for buckling (as measured by its “length”) of the initial blast is inversely proportional to its overlap with the optimal, targeted, blast. The proportionality constant will depend on the pressure. Figure 7 shows how close this proportionality is to a constant by plotting the graph of
$OVEbKECEb$
(11)
as a (nearly constant) function of the blast width βG and for several pressures p ≥ 0.4pc (observe the scale of the y-axis, which is a narrow range compared with changes by factor of up to 4 in Fig. 5). The largest deviation of the proportionality from a constant occurs for lower pressures. This deviation coincides with the observation that at lower pressures the blasts in the direction of the saddle $(vn↦s)$ can no longer cause buckling with an energy close to the energy barrier Eb. So, this deviation may be caused by the fact that the overlap here is computed with respect to a less optimal direction (with respect to buckling energy).

### 4.2 The Nature of Buckling Thresholds and Critical Trajectories.

In this subsection, we provide some details of the nature of the critical trajectory associated with the buckling threshold. Figure 8 shows two examples of critical blasts at the buckling threshold for background pressure p = 0.6pC. The examples have initial velocities with Gaussian shapes of different widths βG (βG = 0.2βref on the left, βG = 1.2βref on the right), resulting in different thresholds KEC = KE(0). These threshold values equal the value of KE at t = 0 in the upper panel: very slightly below 2Eb for βG = 0.2βref and at about 3Eb for βG = 1.2βref. The lower panels of Fig. 3 show the full spatio-temporal profile of the trajectory w(θ, t)/h as a contour plot.

In Fig. 8, the time profiles of (wpole(t) − wpole,s)/h of the buckling trajectory at the threshold perform small oscillations around the saddle equilibrium for times t/T0 > 1. Study [12], which applied spatially uniform pressure steps, discussed the typical spatio-temporal shape of critical trajectories in detail. A critical trajectory does not reach the saddle equilibrium (in the limit of zero damping) but rather its center-stable set (center-stable manifold) at an energy level that is slightly higher than that of the saddle equilibrium. This extra energy is visible in the form of small-amplitude waves across the surface in Fig. 8. In Fig. 8 for both examples, the trajectory approaches the neighborhood of the saddle within t < T0 and then performs small-amplitude oscillations around the saddle, before escaping toward the buckling regime (which is not shown in Fig. 8).

The time profiles of the strain energy SE and the potential energy PE (proportional to the deformation volume ΔV) both start at the values for the stable equilibrium. They show that the nature of the small-amplitude oscillations around the saddle is different for the very sharply localized blast with βG = 0.2βref (left) compared with the large-width blast with βG = 1.2βref (right). The left blast causes small-amplitude traveling waves across the surface (diagonal stripes in lower left panel) before developing standing wave patterns near the pole (spots at the top right of the bottom left panel). The traveling waves do not contain much strain energy SE or potential energy PE such that SE and PE are close to their saddle value at least for t from T0 to 3T0. (Figure 8 panels show (− PE + PEs)/Eb).) The kinetic energy KE is also close to zero for t from T0 to 4T0. Note that the y-axis of the upper left panel and the range of shadings in the lower left panel are smaller than in the upper and lower right panels.

Comparatively, the large-width threshold blast in the right panels shows volume oscillations which are more spatially uniform, visible in the oscillations of PE and SE in the top right and the vertical strips in the bottom right. This is similar to the spatial patterns observed in applying spatially uniform pressure as a sudden step at t = 0 which was analyzed and discussed in our earlier paper [12], We conclude that the Gaussian blast with βG = 1.2βref can indeed be considered as having a large width. Further increases of the width will lead to patterns of larger spatial scale with volume (and PE and SE) oscillations of increasing amplitude at the threshold, with correspondingly increasing necessary initial kinetic energy KE. The study [12] of step loading of uniform pressure revealed that the work done by the pressure greatly exceeded the energy barrier associated with the pressure, so much so that the energy barrier was irrelevant for this loading. Unlike the blast-like impacts considered here which are local and comparable in width to the saddle buckling mode, the step loaded pressure is applied over the entire shell.

The study [12] found that the threshold value reported in this way is not necessarily the lowest threshold (such that the attribute “first in time” threshold value was used in Ref. [12]). More precisely, trajectories that do not escape at this first threshold (because their initial kinetic energy KE is below KEc) may pass another (lower) threshold at a later time (delayed buckling). The number and values of these later (and lower) buckling thresholds depend strongly on the amount and precise nature of the small damping (in study [12] and here we only have artificial numerical damping). The number of further delayed lower thresholds may become arbitrarily large as the damping goes to zero.

For this study, we focus only on the first, “immediate” buckling threshold. This first threshold KEc is robust in the sense that it converges to a limit when the damping (which is introduced by the finite step size of the numerical scheme) goes to zero.

### 4.3 Spherical Shells Clamped at $θ0=45deg$.

There are a variety of static and dynamic buckling phenomena displayed by clamped shallow spherical shells, or caps as they are often called, but this paper will not venture into the domain of shallow shells. The kinetic energy thresholds for buckling by the targeted blasts for the perfect shell clamped at $θ0=45deg$ have been included in Fig. 4, illustrating that at least for these impacts the clamping location plays a secondary role. We have also analyzed the effect of Gaussian-shaped impacts, Eqs. (7) and (8), on the buckling of a perfect spherical shell clamped at $θ0=45deg$. Our main objective again is to demonstrate that the clamping location has a minor effect on the buckling threshold for shells that are not shallow. Detailed numerical simulations have been carried out for shells with R/h = 100 and ν = 0.3. With this radius to thickness ratio even for the widest impacts considered (with βG = 1.2βref), the impacted area does not extend more than to about halfway from the pole to the clamped boundary. The impacts considered are shallow, consistent with the nature of the saddle buckling mode.

Simulations computing the threshold values KEC for the full range of widths of Gaussian blasts considered for the hemispherical shell have been carried out for the shells clamped at $θ0=45deg$ for values of pressure in the range 0.25 ≤ p/pC ≤ 0.8. As before, KEC has been evaluated with no allowance for initial dissipation and with it subtracted off. The contour plots of KEC are similar to those for the hemispherical shell in Fig. 5. In particular, the minimum blast energy KEC with respect to the blast width at each p/pC occurs when βG ≅ 0.6βref, just as for the hemispherical shell. The plot of minimum KEC/Eb versus p/pC for these shells has been included in Fig. 6. Separate calculations for Eb as a function of p/pC have been made for the shell clamped at $θ0=45deg$ (c.f., Fig. 1), and these have been used in normalizing the curves for $θ0=45deg$ in Fig. 6. In the range of 0.25 ≤ p/pC < 0.65, there is essentially no difference between the predictions of KEC/Eb for clamping at θ0 = 0 and $θ0=45deg$, nor is there any noticeable difference between Eb for the two cases in Fig. 1. A dependence of both KEC/Eb and Eb on the two clamping locations begins to emerge for p/pC > 0.65. The role of the clamping location will be increasingly influential for more shallow spherical caps. For static buckling, Huang’s [7] analysis reveals significant interaction between clamping location and buckling for caps with $λ≡23(1−ν2)H/h<20$ where H is the height of the cap, and this is likely to apply to the dynamic case too.

## 5 Buckling Thresholds of Clamped Geometrically Imperfect Hemispherical Shells Subject to Gaussian-Shaped Blasts

The examples presented thus far are all for perfect spherical shells. A more complete understanding of the buckling energy barrier and buckling thresholds under dynamic impact requires consideration of the strong sensitivity of the spherical shell under external pressure to initial geometric imperfections [13]. Axisymmetric geometric imperfections have recently been characterized [11,12,14,15] for uniform thickness shells with a middle surface having an initial stress-free, Gaussian-shaped inward dimple at the pole in the form
$wimp=−δimpexp(−(β/βimp)2)$
(12)
The knockdown in the maximum uniform pressure the shell can support, pmax, is exemplified by imperfection widths βimp = βref, where βref is given in Eq. (7). Representative knockdowns for the choice βimp = βref are pmax/pC ≅ 0.74, 0.61, 0.39, 0.24 for δimp/h = 0.1, 0.2, 0.5, 1, respectively [11]. These knockdowns apply to a full spherical shell or to a hemispherical shell clamped at the equator and they are independent of R/h for thin shells.

To illustrate the role an imperfection plays in the relation between the energy barrier and buckling thresholds under local impact, we consider a clamped hemispherical shell with an initial imperfection (12) having βimp = βref. The shell is quasi-statically pre-loaded to a fixed uniform pressure p = 0.5pC and then subject to Gaussian-shaped impacts of the type considered in Secs. 3 and 4 with width equal to the “optimal” value βG = 0.6βref. The shell can only support the pressure load, p = 0.5pC, if δimp/h < 0.325 because, with δimp/h = 0.325, the quasi-static maximum support pressure is precisely pmax = 0.5pC. Thus, the relevant range of imperfection amplitude for the pre-load being considered is 0 ≤ δimp/h ≤ 0.325.

First, some details of the quasi-static pre-load solution relevant to the node, saddle, and energy barrier are presented. Figure 9 presents the deflection modes and the volume changes associated with the node and the saddle equilibrium states for the relevant range of imperfection amplitudes noted above, all with the pressure fixed at p = 0.5pC. As indicated in the figure caption, the nodal states are associated with wpole/h > −0.75 and ΔVVC,hemisphere < 0.495, while the saddle states are associated with wpole/h < −0.75 and ΔVVC,hemisphere > 0.495. The buckling process is localized in the vicinity of the pole (1 < θπ/2 for this radius to thickness ratio), with clamping influencing the deflections only in the range θ < 0.5, independent of the imperfection amplitude.

For imperfection amplitudes approaching the limit, δimp/h = 0.325, the nodal and saddle states coalesce such that the energy barrier vanishes. This is evident in the plot of the dimensionless energy barrier as a function of imperfection amplitude in the left panel of Fig. 10: Eb/[(pCΔVC/2)Ch/R] = 0 for δimp/h = 0.325, consistent with the fact that the shell is unstable at this level of imperfection when the pre-load is p/pC = 0.5. The other limit in Fig. 10 for the perfect shell (δimp/h = 0) is Eb/[(pCΔVC/2)Ch/R] = 0.062 in agreement with the result for p/pC = 0.5 in Fig. 1.

Threshold buckling calculations for the imperfect shells pre-loaded to p/pC = 0.5 have been carried out following the procedures described earlier. The results for the threshold kinetic energy (with and without subtracting off the initial dissipation) of the Gaussian blasts with βG = 0.6βref are included in the left panel of Fig. 10. These same results are replotted in the right panel of Fig. 10 normalized by the energy barrier at the corresponding values of δimp/h. We have not attempted to extrapolate these curves to the limit δimp/h = 0.325. The important feature that stands out is that the ratio of threshold energy to the energy barrier, KEC/Eb, has a relatively weak dependence on δimp/h over most of the relevant range, while the energy barrier itself is a strong function of the imperfection amplitude. Stated otherwise, the trends in KEC/Eb for the imperfect shell are similar to those for the perfect shell, but the normalizing factor in the denominator, Eb, depends rather strongly on p/pC and δimp/h as seen in Fig. 11.

## 6 Implications for Buckling of Imperfection-Sensitive Shell Structures

This study has shown that the energy barrier against buckling for both perfect and imperfect elastic spherical shells subject to a prescribed uniform pressure provides a relevant estimate of the buckling resistance of the shell to extraneous localized blast-like impacts. The kinetic energy imparted to the shell must exceed the energy barrier at the pre-loaded pressure. For impacts whose widths are comparable with the width of the dimple buckle associated with the saddle in the energy landscape, the threshold energy to buckle the shell imparted by the impact may be only slightly larger than the energy barrier. Shells pre-loaded to pressures close to the quasi-static buckling pressure associated with the level of imperfection are precarious because the energy barrier becomes very small, as seen in Fig. 11. The energy barrier plot in Fig. 11 provides a rationale for NASA’s [16] long standing recommendation that thin spherical shells not be loaded to external pressures above p/pC ≅ 0.2. If one can be certain, the imperfection level of a shell is not too large, e.g., δimp/h < 0.5 for Gaussian-type imperfections, then one can be sure that there is a significant energy barrier for buckling if p/pC ≤ 0.2. While all the results in this paper have assumed a pre-load pressure that is held fixed during the impact event, the energy barrier is almost the same for a spherical shell that is loaded to the same pre-load pressure and then buckled with no subsequent change in volume, except in the range of low pressures, as discussed in detail in Refs. [2,5]. We fully expect that the findings in this paper relating threshold buckling energy to the energy barrier will carry over to spherical shells impacted under conditions where the internal volume of the shell remains unchanged, the so-called rigid volume constraint. This expectation is tempered by the fact that we have assumed that inertia of the fluid medium inside and outside the shell can be ignored in modeling the shell dynamics. There is a regime in which the dynamic fluid–structure interaction must be considered and that has not been addressed here. Arguments for assigning buckling knockdown factors for spherical shells under external pressure based on the energy barrier concept have recently been promulgated [3].

While the analysis in this paper has been restricted to be axisymmetric with respect to the pole of the spherical shell, the results apply equally well to thin shells with dimple-like imperfections and localized impacts located anywhere away from the clamped boundary. The widths of the saddle dimple buckle, the most critical imperfections, and the blast impacts, each scale with $Rh$, and are small compared with the size of the shell. Thus, the present results should be applicable to a dimple imperfection at any location with impacts focused on it if they have locally axial symmetry and lie outside the boundary layer influenced by clamping. An experimental demonstration of the susceptibility of clamped hemispherical shells to local dimple buckling at any location away from the boundary is given in Ref. [17]. There is another important consequence of the locality of dimple buckling. Regardless of the presence of other imperfections, if the pre-loaded shell is impacted in a locality that is free of imperfections extending over a region large compared with $Rh$, then the energy barrier should be essentially that of a perfect shell. The simple measure OV, given in Eq. (10), for similarity of the blast with the buckling dimple can also be applied to blast shapes different from Gaussian. Shapes for which OV is maximal would be expected to have the lowest buckling threshold kinetic energy. The measure OV may also provide insights into impacts that do not possess local axial symmetry which have not been addressed in this paper. One might anticipate that such impacts will require more energy to trigger buckling. It is also worth noting that the importance of localized buckling of cylindrical shells under axial compression has also emerged in recent studies [1820], and we anticipate that the principles seen here for the spherical shells will also apply to axially compressed cylindrical shells subject to localized impacts.

This study and the earlier one on step pressure loading [12] have been informed by developments in nonlinear dynamics. There is a rich array of phenomena associated with the dynamic responses of the highly nonlinear equations governing spherical shell buckling which are only hinted at in this paper. The thresholds for dynamic buckling in this paper were computed by limiting the computed responses to times no greater than four times the period of the spherically symmetric vibration mode of the full sphere, i.e., t/T0 ≤ 4. As noted in the paper, had this time limit been set to a larger value, almost certainly somewhat lower thresholds would have been obtained. The earlier paper on step loading [12] explores the dependence of dynamic buckling thresholds on the time limit and the nature of the buckling trajectories in more detail. Nevertheless, we believe the values presented here of the threshold buckling energy to energy barrier, KEC/Eb, for the blast-like impacts to spherical shells are realistic and representative assuming the impacts are aligned with the local imperfections. A further refinement would require greater attention to misaligned impacts and imperfections as well as dissipation associated with both the numerical scheme and physical damping processes in shell.

## Funding Data

• J.S.’ research was supported by funding from the EPSRC via Grant Nos. EP/N023544/1 and EP/N014391/1 (Funder ID: 10.13039/501100000266).

## References

References
1.
Thompson
,
J. M. T.
, and
Sieber
,
J.
,
2016
, “
Shock-Sensitivity in Shell-Like Structures: With Simulations of Spherical Shell Buckling
,”
Int. J. Bifurcation Chaos
,
26
(
2
), p.
1630003
. 10.1142/S0218127416300032
2.
Hutchinson
,
J. W.
, and
Thompson
,
J. M. T.
,
2017
, “
Nonlinear Buckling Behaviour of Spherical Shells: Barriers and Symmetry-Breaking Dimples
,”
Philos. Trans. R. Soc. Lond., Ser. A
,
375
(
2093
), p.
20160154
. 10.1098/rsta.2016.0154
3.
Evkin
,
A. Y.
, and
Lykhachova
,
O. V.
,
2017
, “
Energy Barrier as a Criterion for Stability Estimation of Spherical Shell Under Uniform External Pressure
,”
Int. J. Solids Struct.
,
118–119
, pp.
14
23
. 10.1016/j.ijsolstr.2017.04.026
4.
Virot
,
E.
,
Kreilos
,
T.
,
Schneider
,
T. M.
, and
Rubinstein
,
S. M.
,
2017
, “
Stability Landscape of Shell Buckling
,”
Phys. Rev. Lett.
,
119
(
22
), p.
224101
. 10.1103/PhysRevLett.119.224101
5.
Hutchinson
,
J. W.
, and
Thompson
,
J. M. T.
,
2018
, “
Imperfections and Energy Barriers in Shell Buckling
,”
Int. J. Solids Struct.
,
148–149
, pp.
157
168
. 10.1016/j.ijsolstr.2018.01.030
6.
Marthelot
,
J.
,
López Jiménez
,
F.
,
Lee
,
A.
,
Hutchinson
,
J. W.
, and
Reis
,
P. M.
,
2017
, “
Buckling of Pressurized Hemispherical Shells Subject to Probing Forces
,”
ASME J. Appl. Mech.
,
84
(
12
), p.
121005
. 10.1115/1.4038063
7.
Huang
,
N.-C.
,
1964
, “
Unsymmetrical Buckling of Thin Shallow Spherical Shells
,”
ASME J. Appl. Mech.
,
31
(
3
), pp.
447
457
. 10.1115/1.3629662
8.
Dankowicz
,
H.
, and
Schilder
,
F.
,
2013
,
Recipes for Continuation
,
SIAM
,
. (URL for code: https://sourceforge.net/projects/cocotools/)
9.
Sanders
,
J. L.
,
1963
, “
Nonlinear Shell Theories for Thin Shells
,”
Quart. Appl. Math.
,
21
(
1
), pp.
21
36
. 10.1090/qam/147023
10.
Koiter
,
W. T.
,
1966
, “
On the Nonlinear Theory of Thin Elastic Shells
,”
Proc. Kon. Ned. Ak. Wet.
,
B69
, pp.
1
54
.
11.
Hutchinson
,
J. W.
,
2016
, “
Buckling of Spherical Shells Revisited
,”
Proc. R. Soc. A
,
472
(
2195
). 10.1098/rspa.2016.0577
12.
Sieber
,
J.
,
Hutchinson
,
J. W.
, and
Thompson
,
J. M. T.
,
2019
, “
Nonlinear Dynamics of Spherical Shells Buckling Under Step Pressure
,”
Proc. R. Soc. A
,
475
(
2223
), p.
20180884
. 10.1098/rspa.2018.0884
13.
Koiter
,
W. T.
,
1969
, “
The Nonlinear Buckling Behavior of a Complete Spherical Shell Under Uniform External Pressure, Parts I, II, III & IV
,”
Proc. Kon. Ned. Ak. Wet.
,
B72
, pp.
40
123
.
14.
Starlinger
,
A.
,
Rammerstorfer
,
F. G.
, and
Auli
,
W.
,
1988
, “
Beulen und Nachbeulverhalten von Duennen Verrippten und Unverrippten Kugelschalen Unter Aussendruck
,”
Z. Angew. Math. Mech.
,
68
(
4
), pp.
257
260
.
15.
Lee
,
A.
,
Marthelot
,
J.
,
López Jiménez
,
F.
,
Hutchinson
,
J. W.
, and
Reis
,
P. M.
,
2016
, “
The Geometric Role of Precisely Engineered Imperfections on the Critical Buckling Load of Spherical Elastic Shells
,”
ASME J. Appl. Mech.
,
83
(
11
), p.
111005
. 10.1115/1.4034431
16.
NASA
,
1969
, “
Buckling of Thin-Walled Doubly Curved Shells
,”
NASA Space Vehicle Design Criteria, National Aeronautics and Space Administration
,
Washington, DC
,
Technical Report No. NASA SP-8032
.
17.
Marthelot
,
J.
,
Brun
,
P.-T.
,
López Jiménez
,
F.
, and
Reis
,
P. M.
,
2017
, “
Reversible Patterning of Spherical Shells Through Constrained Buckling
,”
Phys. Rev. Mater.
,
1
(
2
), p.
025601
. 10.1103/PhysRevMaterials.1.025601
18.
Horák
,
J.
,
Lord
,
G. J.
, and
Peletier
,
M. A.
,
2006
, “
Cylinder Buckling: The Mountain Pass as an Organizing Centre
,”
SIAM J. Appl. Math.
,
66
(
5
), pp.
1793
1824
. 10.1137/050635778
19.
Kreilos
,
T.
, and
Schneider
,
T. M.
,
2017
, “
Fully Localized Post-Buckling States of Cylindrical Shells Under Axial Compression
,”
Proc. R. Soc. A
,
473
(
2205
), p.
20170177
. 10.1098/rspa.2017.0177
20.
Groh
,
R. M. J.
, and
Pirrera
,
A.
,
2019
, “
On the Role of Localizations in Buckling of Axially Compressed Cylinders
,”
Proc. R. Soc. A
,
475
(
2224
), p.
20190006
. 10.1098/rspa.2019.0006