Abstract

Throughout U.S. Department of Energy (DOE) complexes, safety engineers employ the five-factor formula to calculate the source term (ST) that includes parameters of airborne release fraction (ARF), respirable fraction (RF) and damage ratio (DR). Limited experimental data on fragmentation of solids, such as ceramic pellets (i.e., PuO2), and container breach due to mechanical insults (i.e., drop and forklift impact), can be supplemented by modeling and simulation using high fidelity computational tools to estimate these parameters. This paper presents the use of Sandia National Laboratories' SIERRA solid mechanics (SM) finite element code to investigate the behavior of the widely utilized waste container (such as 7A Drum) subject to a range of free fall impact and puncture scenarios. The resulting behavior of the container is assessed, and the estimates are presented for bounding DRs from calculated breach areas for the various accident conditions considered. This paper also describes a novel multiscale constitutive model recently implemented in SIERRA/SM that simulates the fracture of brittle materials such as PuO2 and determines ARF during the fracture process. Comparisons are made between model predictions and simple bench-top experiments.

Introduction

Safety analysts throughout the U.S. Department of Energy (DOE) complex rely heavily on the data provided in the DOE Handbook DOE-HDBK-3010 [1] to determine radionuclide source terms (STs) in support of safety and risk analyses. A five-factor equation of the multiplication of five parameters form the ST
(1)

where MAR is the material-at-risk, which represents the materials of concern, DR is the damage ratio, which is a fraction of MAR contributing to the accident, ARF is the airborne release fraction, which is a fraction of material effects that results in airborne aerosol, RF is the respirable fraction, which is a fraction of ARF particulates or aerosol which are less than or equal to 10 μm in size, and LPF is the leak path factor, which represents the fraction of the suspended material released from a facility. In calculating source terms, analysts tend to use the DOE-HDBK-3010's bounding values on ARF and RF for various categories of insults (representing potential accident release categories). This is typically due to both time constraints and the avoidance of regulatory critique.

Advances in computing capability at U.S. national laboratories have enabled analysts to use computer simulations to model structural and thermal/fluid dynamic phenomena. This gives more insight into physics related to potential accident scenarios. The SIERRA software suite (hereafter referred to as SIERRA), developed at Sandia National Laboratories (SNL), is a prominent tool within SNL for numerical simulation, and can perform source term calculations of the ARFs [2]. In recent years, this tool has been used to substantiate data for the DOE-HDBK-3010 [36]. SIERRA solid mechanics (SIERRA/SM) has enabled us to better understand the physics of mechanical insults related to solids and waste containers compared to the DOE-HDBK-3010, which contains limited experimental data on fragmentation of solids, such as ceramic pellets (i.e., PuO2), and container breach due to mechanical insults. Explosion-induced fragment impact and forklift impact experiments can be supplemented by modeling and simulation using SIERRA/SM [7].

This paper discusses the progress made in the study of waste container breach and fragmentation of composite solids caused by mechanical insults. Several accident simulations were conducted to determine the extent of damage in the container when a breach occurred. We then used the breach area to estimate DR [6]. Another issue the paper addresses is when the contents of a container are brittle (or consist of ceramic) materials. With sufficient impart energy, they will fragment into fine particles and may release through the breached opening of the container. The second half of this paper gives a fragmentation model study for estimating the airborne release of the micron sized particulates from a brittle or ceramic solid [6].

Container Breach Simulation

As previously described, there is a limited discussion in DOE-HDBK-3010 regarding waste container breach and the subsequent airborne release estimate. DOE standard DOE-STD-5506 [8] discusses the DR data on drop and puncture data related to the mechanical insults of waste containers. However, limited experimental data (i.e., drop height, drop orientation, content type/weight, and impact configuration) is provided in this standard. Our motivation behind the simulation study is to provide estimates of the DR resulting from a breach of the waste containers from free-fall and puncture due to a forklift impact using SIERRA/SM in a range of scenarios that can occur in a DOE waste processing facility. The determination of the scenario selection and the waste contents were based on the recent waste container inventory information, which is described later. Once this information was processed, the simulations were performed. The breach area was calculated instead of the DR because it is difficult to create a generalized waste configuration from the simulations. No single waste container was generated or packaged in an identical way. Additionally, polyethylene bags and sealed liners may be used in many cases to prevent releases if the bag or liner remains intact. After that, an attempt was made to estimate DR from the calculated breach area using the available experimental data from Appendix C of the DOE-STD-5506 [8].

The accident scenarios were chosen based on likelihood of events that could occur at DOE facilities such as at the Waste Isolation Pilot Plant (WIPP) in Carlsbad, NM, U.S. Free-fall heights of the waste containers were distributed across a range of heights based on the stacking drum height and bridge crane height in WIPP (2.5, 5, 10, and 20 m drop). The drop angles of the drums were: 0, 45, 90, 135, and 180 deg. The forklift impact was distributed across a range of operational speeds for the vehicle (0.5, 2.5, and 5 m/s). The locations of the drum subject to forklift impact were lid-center, side-lid, side-middle, and bottom-center. To produce the most conservative estimates of damage and increase the likelihood of achieving breach, the impact locations for the drop and forklift tests were chosen to be the weakest parts of the containers. Additionally, some of the key features of the analysis supporting conservatism and computational tractability include:

  • Rigid bodies were used for the ground and wall which the containers impacted. Because the containers have less than ¼ in. thickness of steel and different parts are connected by bolts, while the ground and wall were assumed to be concrete with 0.3048 m (1 ft) of thickness, the rigid body assumption is reasonable.

  • Due to the computational expense of performing all the desired accident scenarios with the highly refined meshes used, several choices were made to reduce the computational cost. Among the choices are: A half symmetry model in YZ plane was used with the corresponding symmetry boundary condition. This decision reduced the number of elements and nodes by half and applicable to the forklift impact cases chosen (see Ref. [6] for details).

  • It was assumed that the conservatism introduced in the analysis is significantly larger than the potential discrepancy of predicted breaches between a full model and a half-symmetry model. Future calculations will address this assumption. For completeness, a mesh convergence study will be performed for a single accident scenario to investigate the adequacy of the mesh refinement used in this work.

Waste Stream Data and Container Selection.

Transuranic waste stream data were obtained from Los Alamos National Laboratory (LANL) for this analysis [9]. The container database was sorted by container type to determine the most common container types: 55-gal (7A) drum (volume is about 0.208 m3) and standard waste box, which has a volume of 1.89 m3. These two containers were selected in our study. Note that this paper will not discuss the standard waste box further (see Ref. [6] for more details). The detailed analysis of this LANL data is documented in Ref. [6].

Based on the analysis of the LANL data, only three content types showed significant percentages for the 7A: (1) contaminated soil/debris (CS) consisting of soils and other sand like debris, (2) solidified inorganics (SI) consisting of large chunks of inorganics, and (3) heterogeneous debris (HD) consisting of laboratory trash, including gloves, wipes, etc. Thus, these three waste contents were chosen in our study.

Constitutive Modeling of Container Contents.

The SI consisted mainly of concrete debris which has been studied extensively and the SIERRA concrete model kc_concrete was used with parameters based on Refs. [1012]. For CS, the soil/foam model in SIERRA was used with the model parameters determined by Refs. [1316].2 HD consisted of plastic, rubber, and cellulose, with the density based on CS, and the Young's modulus as an average of the three materials. The HD used the orthotropic_crush model provided in SIERRA with the crush parameters based on common values for plastic and rubber. A broad summary of the three material types is as follows: SI is a brittle solid block, HD behaves like a single elastic block, and CS behaves like particles. Even though the SI has a different density than CS and HD, the volume of material was chosen so that the final mass of all three would be 453.6 kg.

7A Model Development.

The container models for this study were developed using the Cubit mesh generator [17]. The 7A model was based on the detailed SIERRA/SM model developed previously [5]. We only describe the development of these container models briefly in the Model Development section (see Ref. [6] for more details).

Model Development.

The 7A model is based on the Skolink drum specification [18] given in Fig. 1, showing various components modeled. Figure 1 shows the closure ring which holds the lid and the drum body together. The ring itself clamps together the lid and body though the use of a closure bolt, which during a preload step in the calculation, closes the gap between the two ends of the ring. The bung hole and vent use the spot weld feature in SIERRA to connect them with their caps. Spot welds allow for the cap to separate from the bung hole or vent if the interface force exceeds a user defined limit. Note that the mesh is entirely made of hexahedral elements and with single-point integration.

Fig. 1
7 A component model mesh
Fig. 1
7 A component model mesh
Close modal

The specifications of the 7A identify construction using different steels [19]. To simplify our model, the cylindrical body, lid, and bung cap are modeled using ASTM A1008 steel [20]. The closure ring is also assumed to be ASTM A1008 steel but with modifications on properties such as yield stress and parameters that define cracking and other failure modes to capture the specified steel type. The bolt is modeled with ASTM A307 steel [21]. The seal is modeled as a Mooney–Rivlin material.

As previously described, the material for the content models were from (1) HD, (2) SI waste [1012], and (3) CS [1316].2 Their densities are 1108.5, 2325.0, and 1755.5 kg/m3, respectively.

As previously described, two accident scenarios are modeled in this study: free-fall drop and impact from a forklift tine. Below is the brief description of these models. Figure 2 shows examples of the 7A simulation geometries and the details of element removal due to gross deformation.

Fig. 2
Before and after highly deformed (dead) elements are removed
Fig. 2
Before and after highly deformed (dead) elements are removed
Close modal

The drop model is based on the free-fall of the 7A onto a rigid surface with various drop heights and angles of impact. As shown in Fig. 3, the 7A is oriented at a variety of angles above the target surface. To expedite the simulation of the drops, an initial velocity equivalent to a drop at the height was applied so that simulations could begin at the point just before hitting the target surface.

Fig. 3
Examples of 7A simulation geometries: (a)–(e) for drop orientation and (f) tine puncture: (a)0-degree, (b) 45-degree, (c) 90-degree, (d) 135-degree, (e) 180-degree, and (f) forklift tine hitting 7A mid-Level against rigid wall
Fig. 3
Examples of 7A simulation geometries: (a)–(e) for drop orientation and (f) tine puncture: (a)0-degree, (b) 45-degree, (c) 90-degree, (d) 135-degree, (e) 180-degree, and (f) forklift tine hitting 7A mid-Level against rigid wall
Close modal

The forklift tine accident scenarios varied the orientations of the 7A, the location where the tine contacted the 7A, and whether a rigid wall is present (see Fig. 3(f)). To make the tine have an impact equivalent to one which is connected to a forklift, the density of the tine mesh was set so its total mass would be the mass of an average forklift (4082 kg based on Ref. [21]). Because the half-symmetry model is assumed, the tine volume is reduced by half (see Fig. 3(f)). The material model for the tine was chosen to be linear elastic steel. With the full forklift mass concentrated in the tine, and some points of impact chosen to be the most vulnerable parts of the 7A like the interface between the lid and the body, our results are conservative. To provide the model adequacy to use this half-symmetry and mesh size model, we developed the model to compare against the known rod-drop experiment (see Fig. 4). As shown Fig. 4(a), the experiment was conducted from a 1 m penetration bar dropped onto near the center of the lid. Figure 4(b) shows the result of the half-symmetry simulation of this experiment. A numerical simulation of a rod dropping from one meter onto the lid of the drum was done. Instead of placing the rod one meter above the lid surface of the drum, we started the simulation assuming an initial velocity of the rod at 4.29 m/s, corresponds to a 1 m high drop. The content in the drum was assumed as maximum mass of CS. The simulation result showed no puncture of the lid in this drop height, which is consistent to the finding in the experiment [22]. Although we did not provide a mesh convergence study for this drum model, the full symmetry drum model had been used previously to analyze a drum explosion accident occurring in a U.S. nuclear waste disposal [23]. Then a half-symmetry drum model was applied to simulate the contaminant release from a breached drum in a jet fuel fire experiment conducted at SNL's FLAME facility [5]. This half-symmetry drum model was used by several SIERRA codes, including SM, to model the pressure and temperature of the drum under this fire condition, the drum burst due to the thermal and pressure conditions, and the release of contaminant from the drum burst. Then this half-symmetry drum model was applied to this current research as described in this article.

Fig. 4
Simulation of a steel rod dropped one meter onto a 7A container: (a) experiment Setup [23] and (b) SIERRA/SM Simulation (von Mises Stress)
Fig. 4
Simulation of a steel rod dropped one meter onto a 7A container: (a) experiment Setup [23] and (b) SIERRA/SM Simulation (von Mises Stress)
Close modal

Simulation Discussion and Results on Damage Ratio Estimates.

This section describes the simulation results of the half-symmetry models for the 7A. Four 7A drop failure modes were considered (see Fig. 5). As shown in Fig. 5, the observed failure modes due to the drop accidents include: lid rupture, complete loss of lid, lid vent breach, and bottom failure of the container.

Fig. 5
Typical 7A drop failure modes: (a) Failure due to lid and sidewall separation, (b) Complete lid failure, (c) Lid vent failure, and (d) Bottom failure
Fig. 5
Typical 7A drop failure modes: (a) Failure due to lid and sidewall separation, (b) Complete lid failure, (c) Lid vent failure, and (d) Bottom failure
Close modal

After completing the simulations, possible sources of uncertainty in the results must be enumerated: (1) there could be drop cases where the 7A rebounds and hits the ground again which may cause further damage. (2) The state of container breach can be subjective for some cases. If the lid came off the container, it was considered fully open, but for many cases, the lid only deformed without necessarily creating an opening. Video and still images were postprocessed for the simulations and used to help determine by visual inspection whether a breach occurred. (3) The breach itself was calculated manually from the partial opening of the lid or holes created by “element death,” a feature in SIERRA where highly deformed elements are removed. The results of the simulations for the 7A are described below. Table 1 shows the different SIERRA/SM scenarios simulated for the free-fall drop and forklift tine puncture for the 7A.

Table 1

Simulation scenarios for 7A in free-fall drop and forklift tine puncture

DropPuncturea
ParametersCasesCases
Content type: HD, SI, Cs3Content type: HD, SI, Cs3
Drop height (m): 2.5, 5, 10, 204Constraint: free or fixed2
Orientation (deg): 0, 45, 90, 135, 1805Impact location: lid-center, side-lid, side-middle, bottom-center4
Weight (kg): 453.61Impact speed (m/s): 0.5, 2.5, 53
Total6072
DropPuncturea
ParametersCasesCases
Content type: HD, SI, Cs3Content type: HD, SI, Cs3
Drop height (m): 2.5, 5, 10, 204Constraint: free or fixed2
Orientation (deg): 0, 45, 90, 135, 1805Impact location: lid-center, side-lid, side-middle, bottom-center4
Weight (kg): 453.61Impact speed (m/s): 0.5, 2.5, 53
Total6072
a

Volume filled based on content density and mass of 453.6 kg.

Example Breach Area Calculation.

As previously described, Fig. 2 shows before and after images where dead elements (elements that have been severely distorted) have been removed. A surface was manually created from the nodes defining the boundary of the hole. The area of that surface represents the breach area. Though subjective, some calculations were repeated by different team members to ensure consistency.

For the case of breaches involving the lid, full breach is obtained if the lid came off entirely. If the lid is partially attached, the method described above was used to calculate the breach surface in the opening formed between the lid and the container (see Fig. 5 for lid failure examples). This is one method among others to calculate breach area.

Damage Ratio Estimate Methodology.

Damage ratio is not calculated, but instead is a quantity obtained qualitatively from experiments, the DOE Handbook, or standards such as DOE-STD-5506 [8]. In addition to the breach area, the DR is also dependent on the waste material type and the impact scenario. Because DOE-STD-5506 did not cover the range of accidents in our simulations, we had to extrapolate from the available data to estimate the DR for the accidents in this paper. The method to estimate the DR considers the calculated breach area, the existing DR data from DOE-STD-5506 and engineering judgment—either from free fall or puncture as briefly described in the 7A Model Development section. However, a detailed formulation and extensive assumptions are described in Ref. [6].

Simulation Results.

Among the 60 simulations for the 7A drop cases as shown in Table 1, the results show only 10% of HD content type, 30% of SI and 25% of CS yielded failure by visual inspection of the graphical results from the simulations. The 10 m and 20 m drop height cases yielded failure for all content types. Only the SI content type yielded failure at 5 m drop. In terms of drop orientations, the 45 deg, 180 deg, and 135 deg yielded failure for all content types. Table 2 shows the failure cases and their mode of failure for the 7A simulations. As shown in this table, only HD yielded a failure mode of a complete lid-sidewall separation. Because the orientation of the 7A has the lid facing downward when hitting the ground target, all content types yielded failure of the vent which is located on the lid.

Table 2

Simulation results for free-fall drop cases yielding 7A failure

Content typeDrop height (m)Orientation (deg)Breach area (10−3 m2)Failure modeaBreach area (10−3 m2)Damage ratio
HD201354.605148LSS9.210.04b
HD201800.743268LVF1.490.01c
SI51800.808858LVF1.620.01c
SI10135127.13272CLF254.300.5d
SI101800.788763LVF1.580.02c
SI20454.387404BT8.780.04c
SI20135124.95101CLF250.000.5d
SI201800.815839LVF1.630.04c
CS10455.464377BT10.930.01c
CS10135125.68497CLF251.401.0d
CS204586.786316BT173.601.0d
CS20135124.45447CLF248.901.0d
CS201800.795928LVF1.590.02c
Content typeDrop height (m)Orientation (deg)Breach area (10−3 m2)Failure modeaBreach area (10−3 m2)Damage ratio
HD201354.605148LSS9.210.04b
HD201800.743268LVF1.490.01c
SI51800.808858LVF1.620.01c
SI10135127.13272CLF254.300.5d
SI101800.788763LVF1.580.02c
SI20454.387404BT8.780.04c
SI20135124.95101CLF250.000.5d
SI201800.815839LVF1.630.04c
CS10455.464377BT10.930.01c
CS10135125.68497CLF251.401.0d
CS204586.786316BT173.601.0d
CS20135124.45447CLF248.901.0d
CS201800.795928LVF1.590.02c
a

LSS = lid, sidewall separation, CLF = complete lid failure, LVF = lid vent failure, BT = bottom failure.

b

Even though there is clear separation, but because the opening is small, it is assumed only very small release.

c

Small breach area, such as lid vent failure and minor bottom separation from sidewall of the container. It is assumed at DR of 0.01, then progress doubling if the kinetic energy is also doubled.

d

Complete lid failure or bottom failure with a very large breach area calculated. Assumed that sufficient kinetic energy remains to push majority of the materials out for CS type content. For both HD and SI, it is assumed only 50% material release, due to solids and large debris, which unlike sands/soils can easily flows out.

Among the 72 simulations for the 7A puncture by the forklift tine (impact) cases, as shown in Table 3, the results show only 25% of HD content type, 42% of SI, and 25% of CS yielded failure. When the 7A was held in place by a wall, all three content types yielded failure. For SI type, 30% of failure cases are from cases without constraint. As shown in Table 3, lid and sidewall separation always occur when the puncture is located at the “side-near lid” scenario when the drum is upright and fixed at the wall, and the tine hits the lid.

Table 3

Simulation results for forklift tine puncture cases yielding 7A failure

Content typeConstraint, impact location, impactor speedaFailure modebTotal breach area (10−3 m2)Damage ratio
HDC(Fi), IL(LC), IS (2.5)P5.280.1
HDC(Fi), IL(LC), IS (5)P11.540.2
HDC(Fi), IL(SL), IS (2.5)LSS3.710.01c
HDC(Fi), IL(SL), IS (5)P, LSS28.540.02c
HDC(Fi), IL(BC), IS (2.5)P3.370.1
HDC(Fi), IL(BC), IS (5)P8.140.2
SIC(Fi), IL(LC), IS (2.5)P16.250.1
SIC(Fi), IL(LC), IS (5)P13.740.2
SIC(Fi), IL(SL), IS (2.5)LSS4.520.01c
SIC(Fi), IL(SL), IS (5)LSS98.310.02c
SIC(Fi), IL(SM), IS (5)P8.630.2
SIC(Fi), IL(BC), IS (2.5)P5.580.1
SIC(Fi), IL(BC), IS (5)P6.450.2
SIC(Fr), IL(LC), IS (5)P9.130.2
SIC(Fr), IL(BC), IS (2.5)P0.530.1
SIC(Fr), IL(BC), IS (5)P4.970.2
CSC(Fi), IL(LC), IS (2.5)P5.440.5d
CSC(Fi), IL(LC), IS (5)P8.920.5d
CSC(Fi), IL(SL), IS (2.5)LSS94.490.01c
CSC(Fi), IL(SL), IS (5)LSS165.900.02c
CSC(Fi), IL(BC), IS (2.5)P3.910.5d
CSC(Fi), IL(BC), IS (5)P8.010.5d
Content typeConstraint, impact location, impactor speedaFailure modebTotal breach area (10−3 m2)Damage ratio
HDC(Fi), IL(LC), IS (2.5)P5.280.1
HDC(Fi), IL(LC), IS (5)P11.540.2
HDC(Fi), IL(SL), IS (2.5)LSS3.710.01c
HDC(Fi), IL(SL), IS (5)P, LSS28.540.02c
HDC(Fi), IL(BC), IS (2.5)P3.370.1
HDC(Fi), IL(BC), IS (5)P8.140.2
SIC(Fi), IL(LC), IS (2.5)P16.250.1
SIC(Fi), IL(LC), IS (5)P13.740.2
SIC(Fi), IL(SL), IS (2.5)LSS4.520.01c
SIC(Fi), IL(SL), IS (5)LSS98.310.02c
SIC(Fi), IL(SM), IS (5)P8.630.2
SIC(Fi), IL(BC), IS (2.5)P5.580.1
SIC(Fi), IL(BC), IS (5)P6.450.2
SIC(Fr), IL(LC), IS (5)P9.130.2
SIC(Fr), IL(BC), IS (2.5)P0.530.1
SIC(Fr), IL(BC), IS (5)P4.970.2
CSC(Fi), IL(LC), IS (2.5)P5.440.5d
CSC(Fi), IL(LC), IS (5)P8.920.5d
CSC(Fi), IL(SL), IS (2.5)LSS94.490.01c
CSC(Fi), IL(SL), IS (5)LSS165.900.02c
CSC(Fi), IL(BC), IS (2.5)P3.910.5d
CSC(Fi), IL(BC), IS (5)P8.010.5d
a

For constraint type: Fi = fixed, Fr = free; impact location type: LC = Lid-center, SL = side near lid, BC = bottom center, SM = side middle; impact speed type in m/s: 2.5 or 5.

b

P = puncture, LSS = lid, sidewall separation.

c

For the impact location at SL, it is assumed a small release.

d

It is assumed to be 50%, because the impact location is midpoint, regardless the extent of the opening.

The breach area was calculated based on the failure results from Section Example Breach Area Calculation. In addition, the DR was estimated based on the experimental data from DOE-STD-5066, breach area estimates and engineering judgment. In brief, DR represents the amount of MAR that are involved in the accident. In this case, the amount of MAR released from the container through the breach area. Readers are encouraged to consult [6] for details of derivations of the breach area and DRs obtained. As shown in Tables 2 and 3 for the DR values, depending on the size of the opening and the waste type, the DR is estimated. For example, it is assumed that DR is 0.5 for both HD and SI waste types because the waste is either solid or large debris, which may not readily release out of the breach. However, CS waste type is assumed to behave like sand; it can more readily be released out of a breach. The footnotes in the tables provide the assumptions and rationale for the DR values provided. For example, in Table 3, when the location of the impact is at the centerline of the container, the largest DR that can be assumed is 0.5 for the CS cases because of its sand-like composition.

Macro/Microscale Fragmentation

SIERRA/SM can be used to study fragmentation of composite solids, such as ceramics. The approach, demonstrated here, uses a concurrently coupled fragmentation model consisting of both a macroscale model and a microscale model, applied in SIERRA/SM. This hierarchical model enables the break-up of a composite solid into large fragments (called macrofragments, length scale ≥ 1 mm) and small fragments (called microfragments, length scale < 1 mm). In the Approach section, coupling of the microscale model to an existing macrofragmentation model is described. A validation study of the coupled model and a sensitivity analysis of the microscale model are also presented. This section concludes with a recommendation for future study of this approach and its application to the Handbook.

Approach.

In the fragmentation formulation, the brittle material is represented at the macroscale using the gradient damage explicit (GDE) model [24]. The GDE model calculates the macroscale fragmentation behavior of the material and removes elements from the model when they are damaged past a specified threshold. A one-dimensional (1D) microscale fragmentation model is tightly coupled to the macroscale model. This coupled model is implemented as a Library of Advanced Materials for Engineering material model in SIERRA/SM called the hierarchical gradient damage explicit model. Only brief descriptions of the models are provided in this paper. Detailed descriptions are found in Ref. [6].

Starting from an energetic point of view, the Helmholtz free-energy (ψ) depends on strain (ϵ) and phase field damage (ϕ). If ϕ = 0, the material is undamaged, whereas ϕ = 1 indicates the material is fully damaged. The Helmholtz free-energy is defined as
(2)

where g(ϕ)=(1ϕ)21+(m2)ϕ+(1+pm)ϕ2 defines stiffness as a function of damage and m=32EGcσc2L. The parameters in these equations are the positive and negative components of the strain energy (ψe+andψe), the fracture length scale (L), the fracture energy (Gc), Young's modulus (E), and critical stress value (σc). The parameter p is subject to the limitations p1 and pm2. It influences the shape of the stress decay following the initiation of damage accumulation.

For example, if a stress-based approach is implemented, damage can only accumulate if the maximum element principal stress is greater than the element critical stress. In the strain-energy-based approach, damage evolution can only occur when strain energy exceeds a critical strain energy. The microscale model is initiated when the critical threshold (either stress or strain-energy based) is exceeded and the element begins to accumulate damage. After incremental solution of the microscale model, the element fracture energy at the macroscale is updated to reflect the energy that has been dissipated at the microscale. This results in the macroscale model dissipating the same amount of energy as the microscale model, which leads to an energy consistent implementation of the coupling between the two length scales.

The microscale model takes a 1D approach to describing fragmentation because of elastic wave propagation coupled with a cohesive failure process. This approach was developed by Zhou et al. [25] and was used for the previously implemented sequential model [5]. In the microscale model, damage can only initiate when the local node stress (σi) exceeds the local critical stress (σc,i), σiσc,i. When this threshold is exceeded, crack opening distance and the local node stress is determined by the cohesive law and related by the following two equations:
(3)
(4)

As illustrated in Fig. 6, δcoh,δc, and δmax are the current crack opening distance for the cohesive, critical, and maximum, respectively. D is the damage number associated with the node, which is defined as D=min(δmaxδc,1). When a node is damaged, (D>0), the energy that has been dissipated at the node is equal to Di2Gc,i where Gc,i is the nodal fracture energy. The total energy dissipated by fracture at the microscale is defined for N grain boundaries as Energydissipated=i=1NDi2Gc,imax(Di2Gc,i). The energy dissipated at the microscale is used to update the fracture energy at the macroscale, which ensures that both length scales are consistent. With Energydissipated as the energy dissipated at the microscale for a given element, the updated macroscale fracture energy is then Gc,macro=Gc,macro,initial+Energydissipated.

Fig. 6
Representation of cohesive law for fragmentation
Fig. 6
Representation of cohesive law for fragmentation
Close modal

The projection of the stress on the grain boundary surface is illustrated in Fig. 7. In the case where the grain boundary orientation angle is zero, this means that the grain boundary is orientated in such a way that all applied force acts on the grain boundary surface. In the case where the grain boundary orientation is π/2 radians, the grain boundary surface is not aligned with the applied force direction, so no force acts on the grain boundary surface. A uniform distribution was assumed for the grain orientation angles.

Fig. 7
Projection of the stress onto grain boundary surface at angle θ. F is applied load, A is cross-sectional area, σ is normal stress, τ is shear stress, x and y are the plane axes.
Fig. 7
Projection of the stress onto grain boundary surface at angle θ. F is applied load, A is cross-sectional area, σ is normal stress, τ is shear stress, x and y are the plane axes.
Close modal

The implementation allows for the user to directly specify the length of the 1D bar, instead of length being based on the element size. This was implemented to allow the 1D bar to be the length of the phase field width. For the gradient damage explicit model, this is twice the fracture length scale. If a 1D bar length is specified, the damage accumulated on the 1D bar must be scaled such that the microscale model dissipates an amount of energy that is consistent with the energy dissipated for a 1D bar with L=vol3. This behavior was enforced by scaling the amount of dissipated energy by the ratio of the characteristic length of the element to the length of the user specified 1D bar.

The values of critical stress and fracture energy at the macroscale are determined by the microscale model parameter distributions. At the microscale model, each node has a critical stress, fracture energy, and grain boundary orientation parameter. Because grain boundary orientation is assumed, the effective critical stress value that is used to determine the crack opening displacement is a projection of the node critical stress value onto the grain boundary surface, at a grain boundary orientation angle, θ.

Validation.

Three different validation cases are included in this paper. The first two validation cases compare the unirradiated UO2 and Pyrex experiments from Ref. [26] and the Handbook to current simulation results. The final case compares current simulation results of the crushing of a brittle sphere during an impact test to experimental data (see Fig. 8 for the simulation mesh and model).

Fig. 8
Macroscale finite element model mesh and components for UO2 and Pyrex impact tests
Fig. 8
Macroscale finite element model mesh and components for UO2 and Pyrex impact tests
Close modal

The Handbook provides two tests for validation involving the impact of unirradiated UO2 and Pyrex glass. The information for the Pyrex glass was used to establish the Handbook Eq. (4-1) for calculating the fraction of particles that are below 10 μm. However, for a better comparison to experimental data, the UO2 and Pyrex impact tests conducted by Jardine et al. [26] are used for comparison and to perform validation of the hierarchical gradient damage explicit model.

Figures 9 and 10 show the UO2 and Pyrex simulations and comparison to the experimental data, respectively. As shown in these figures, both the UO2 and Pyrex models could capture fragmentation behavior at the macro- and microscales, but the predicted cumulative mass distribution was different from that observed for the experimental data. Figures 9 and 10 show that the current modeling is overpredicting the cumulative mass distribution at the microscale. A proper validation of the microscale model is needed before the coupled model can be properly validated. As the models are coupled, it is difficult to determine the source of discrepancies between experiment and modeling. The major discrepancies in Figs. 9 and 10 occur at the microscale length scale, which is only modeled by the microscale model; however, the microscale model utilizes coupled feedback from the macroscale model; discrepancies that appear in the microscale length scale could originate in either of the coupled models. A complete validation of the microscale model against a different data source would make it easier to identify where the modeling discrepancies occur when validating the coupled model. If this model is to be further developed, this additional validation should be performed.

Fig. 9
Simulation results of UO2 impact test: (a) macroscale stress degradation (legend) predicted by the finite-element model for UO2, 500 μs after impact for the undeformed (left) and deformed (right) mesh and (b) predicted cumulative mass distribution (CMD) (combined macroscale and microscale in black, microscale in blue, macroscale in red) compared to test data (dotted) [26]
Fig. 9
Simulation results of UO2 impact test: (a) macroscale stress degradation (legend) predicted by the finite-element model for UO2, 500 μs after impact for the undeformed (left) and deformed (right) mesh and (b) predicted cumulative mass distribution (CMD) (combined macroscale and microscale in black, microscale in blue, macroscale in red) compared to test data (dotted) [26]
Close modal
Fig. 10
Simulation results of Pyrex impact test: (a) macroscale stress degradation (legend) predicted by the finite-element model for Pyrex, 500 μs after impact for the undeformed (left) and deformed (right) mesh and (b) predicted cumulative mass distribution (combined macroscale and microscale in black, microscale in blue, macroscale in red) compared to test data (dotted) [26]
Fig. 10
Simulation results of Pyrex impact test: (a) macroscale stress degradation (legend) predicted by the finite-element model for Pyrex, 500 μs after impact for the undeformed (left) and deformed (right) mesh and (b) predicted cumulative mass distribution (combined macroscale and microscale in black, microscale in blue, macroscale in red) compared to test data (dotted) [26]
Close modal

The second validation source is from impact tests conducted on brittle spheres performed by Wu et al. [27]. Figure 11 shows the simulation of an impact test performed by Wu to test the macroscale behavior. Only the macroscale behavior was compared to experimental data as Wu did not provide any very small fragment data [27]. Wu observed failure modes of three, four, and five lune shapes slices for impact energies between 35.4 to 39.9 J. This simulation (impact energy 35.7 J) produced four lune shaped slices, which is a plausible result given the failure modes observed by Wu et al. [27]. Additional simulation is needed to identify if increasing the simulation impact energy increases the number of fragments, which was the trend observed by Wu.

Fig. 11
Brittle sphere impact in macroscale simulation (a) model mesh and (b) fracture behavior: (a) macroscale finite element model mesh and components and (b) macroscale fracture behavior, 1 millisecond after impact. Mesh is shown undeformed, with degraded elements removed from view.
Fig. 11
Brittle sphere impact in macroscale simulation (a) model mesh and (b) fracture behavior: (a) macroscale finite element model mesh and components and (b) macroscale fracture behavior, 1 millisecond after impact. Mesh is shown undeformed, with degraded elements removed from view.
Close modal

Microscale Sensitivity With DAKOTA.

A stand-alone version of the microscale model was supplied initial boundary conditions from the SIERRA/SM macroscale model to set the initial one-dimensional strain rate. A stand-alone model was required as sensitivity analysis was too computationally expensive for the fully coupled models. DAKOTA, an SNL-developed computational toolkit for iterative studies such as sensitivity analysis, calibration, and uncertainty quantification, was used to drive the sensitivity analysis. A description of the procedure for this analysis, as well as the results, can be found in Ref. [6].

The sensitivity analysis examined the response of the quantity of interest (percentage of mass under 10 μm for the 1D bar) and the microscale model parameters. The model parameters that were examined are the following, defined at the macroscale (and used to determine the microscale parameters): Young's modulus, density, critical stress, fracture energy grain size, and initial strain rate.

Performing a sensitivity analysis on the microscale model is difficult as some of the parameters (i.e., critical stress and fracture energy) are correlated. The change of the Young's modulus, density, and the grain size will affect the number of time-steps taken per microscale iteration. In the case of the current implementation of the sensitivity analysis, while the same amount of macroscale time elapses per evaluation of the stand-alone model, a different number of time-steps (and total model evaluations) will occur at the microscale, which will influence the results. There are likely other correlations and relationships between input parameters that will affect the model output.

Recommendation for Future Improvement.

Equation (4-1) in DOE-HDBK-3010 is ARF×RF=A·ρ·g·h where A is an empirical correlation, ρ is the specimen density, g is acceleration due to gravity, and h is the fall height. Equation (4-1) shows a direct relationship between the input energy and the ARF×RF. This is consistent with what was seen during the sensitivity analysis of the stand-alone 1D model when varying the initial strain-rate and with the study performed with the coupled models in SIERRA/SM with Pyrex while increasing energy density. The relationship between the energy density and the percentage of the mass under 10 μm was linear for the coupled model, which is consistent with the form of Eq. (4-1). However, the relationship between the percentage of the mass under 10 μm and the energy (in the form of increasing initial strain-rate) in the stand-alone model was not linear.

The sensitivity analysis performed on the stand-alone microscale model did not show any dependence of density on the percentage of the mass of the 1D bar under 10 μm (see Fig. 12). This is consistent with the implementation of Eq. (4-1) in DOE-HDBK-3010; however, density dependence was observed by Koch et al. [28]. As was noted previously, the stand-alone microscale sensitivity model was performed without coupled feedback from the macroscale model. The coupled model may show sensitivity to density and different behavior when compared to the stand-alone model. Follow-up work should be performed to determine if the microscale model shows sensitivity to density.

Fig. 12
Cumulative mass of 1D bar under 10 μm, varying density from 2000 to 11,500 kg/m3, using UO2 properties
Fig. 12
Cumulative mass of 1D bar under 10 μm, varying density from 2000 to 11,500 kg/m3, using UO2 properties
Close modal

Conclusions

This paper demonstrates that the use of SIERRA/SM within the SIERRA engineering mechanics simulation code suite developed at SNL can be used to predict the failure of the waste containers (e.g., 7A) in variety of configurations and conditions related to mechanical insults including, but not limited to, the conditions described in DOE-STD-5506; i.e., free-fall drop and puncture. This paper also demonstrates that extending the capability of SIERRA/SM can be used to predict the fragmentation of both brittle and ceramic type solids, which can be used to estimate the ARF and RF if fragmentation of the brittle/ceramic in a breached drum can occur. However, additional improvement and validation testing are needed to ensure the success of micro/macroscale fragmentation models for use in this type of analysis.

Acknowledgment

This paper describes the objective technical results and analysis. Any subjective views or opinions that might be expressed in the paper do not necessarily represent the views of the U.S. Department of Energy or the U.S. Government. The authors would like to thank the SNL staff who worked on this study, including Remi Dingreville, Natalie Gordon, Dominic G. Fascitelli, and former staff, John Bignell, including a special thanks to John Wilkes and Michael Starr for reviewing this paper.

Funding Data

  • SNL is managed and operated by National Technology and Engineering Solutions of Sandia, LLC for U.S. DOE/NNSA under Contract No. DE-NA0003525.

  • This work was supported by the DOE Environmental, Health, Safety and Security (EHSS) Nuclear Safety Research and Development Program (WA No. 2017-AU33-SNL-NSRD Rev 1).

Nomenclature

A =

cross-sectional area (m2)

A =

empirical coefficient

ARF =

airborne release fraction

c =

subscript indicating critical value

coh =

subscript indicating cohesive zone

CS =

contaminated soil/debris

D =

damage number associated with a node

DR =

damage ratio

DOE =

U.S. Department of Energy

E =

Young's modulus, Pa

Energydissipated =

the energy dissipated at the microscale for a given element, J/m2

F =

applied load, N

G =

fracture energy, J/m2

GDE =

gradient damage explicit

g =

gravitational acceleration, m/s2

g =

stiffness parameter in GDE model

h =

fall height, m

HD =

heterogeneous debris

i =

subscript indicating ith node

initial =

subscript indicating initial value

L =

characteristic length (1D bar simulation), m

LANL =

Los Alamos National Laboratory

m =

parameter in GDE model that is used to define stiffness

macro =

subscript indicating value at macroscale

max =

subscript indicating maximum value

MAR =

material-at-risk

micro =

subscript indicating value at microscale

p =

parameter in GDE model that is used to define stiffness

P =

pressure, Pa

RF =

respirable fraction

SI =

solidified inorganics

SM =

solid mechanics

SNL =

Sandia National Laboratories

ST =

source term

WIPP =

Waste Isolation Pilot Plant

vol =

volume of element, m3

δ =

crack opening distance, m

ϵ =

strain

θ =

grain boundary orientation angle, rad

ρ =

density, kg/m3

σ =

stress, Pa

τ =

shear stress, Pa

ϕ =

phase field damage

ψ =

Helmholtz free-energy, J/m3

ψe =

strain energy, J/m3

1D =

one-dimensional

7A =

7A (55-gal) drum

References

1.
DOE
,
1994
, “
Airborne Release Fractions/Rates and Respirable Fractions for Nonreactor Nuclear Facilities, Volume I—Analysis of Experimental Data
,”
U.S. Department of Energy
,
Washington, DC
, Standard No. DOE-HDBK-3010-94.
2.
SIERRA Code Suite, “SNL, Mechanics 4.48 User's Guide,” Sandia National Laboratories, Albuquerque, NM, accessed Apr. 2, 2018,
https://sandia.gov/ ASC/integrated_codes.html
3.
Louie
,
D. L. Y.
, and
Brown
,
A. L.
,
2015
, “
NSRD-6: Computational Capability to Substantiate DOE-HDBK-3010 Data
,”
Sandia National Laboratories
,
Albuquerque, NM
, Report No. SAND2015-10596.
4.
Louie
,
D. L. Y.
,
Brown
,
A. L.
,
Gelbard
,
F.
,
Bignell
,
J.
,
Pierce
,
F.
,
Voskuilen
,
T.
,
Rodriguez
,
S. B.
,
Dingreville
,
R.
,
Zepper
,
E. T.
,
Juan
,
P.
,
Le
,
S.
, and
Gilkey
,
L. N.
,
2016
, “
NSRD-11: Computational Capability to Substantiate DOE-HDBK-3010 Data
,”
Sandia National Laboratories
,
Albuquerque, NM
, Report No. SAND2016-12167.
5.
Louie
,
D. L. Y.
,
Bignell
,
J.
,
Dingreville
,
R.
,
Zepper
,
E. T.
,
O'Brien
,
C. J.
,
Busch
,
R. D.
, and
Skinner
,
C. M.
,
2018
, “
NSRD-15: Computational Capability to Substantiate DOE-HDBK-3010 Data
,”
Sandia National Laboratories
,
Albuquerque, NM
, Report No. SAND2018-0436.
6.
Louie
,
D. L. Y.
,
Bignell
,
J.
,
Le
,
S.
,
Dingreville
,
R.
,
Gilkey
,
L. N.
,
Gordon
,
N.
, and
Fascitelli
,
D. G.
,
2019
, “
NSRD-16: Computational Capability to Substantiate DOE-HDBK-3010 Data
,”
Sandia National Laboratories
,
Albuquerque, NM
, Report No. SAND2019-0290.
7.
SIERRA Solid Mechanics Team
,
2019
, “
SIERRA/SolidMechanics 4.54 Theory Manual
,”
Sandia National Laboratories
,
Albuquerque, NM
, Report No. SAND2019-11007.
8.
DOE
,
2007
, “
Preparation of Safety Basis Documents for Transuranic (TRU) Waste Facilities
,” U.S.
Department of Energy
,
Washington, DC
, Standard No. DOE-STD-5506-2207.
9.
Soest
,
G. V.
,
2017
, Email to David Louie, ‘Re-Digital Data of Site TRU Waste Inventory,
Los Alamos National Laboratory
,
Los Alamos, NM
.
10.
Collins
,
M. P.
, and
Mitchell
,
M. P.
,
1991
,
Prestressed Concrete Structures
,
Prentice Hall
,
Englewood Cliffs, NJ
.
11.
Brannon
,
R. M.
,
2009
, “
Survey of Four Damage Models for Concrete
,”
Sandia National Laboratories
,
Albuquerque, NM
, Report No. SAND2009-5544.
12.
EPA
,
1996
, “
Stabilization/Solidification Process for Mixed Waste
,”
U.S. Environmental Protection Agency
,
Washington, DC
, Report No. EPA 402-R-96-014.
13.
Holtz
,
R. D.
, and
Kovacs
,
W. D.
,
1981
,
An Introduction to Geotechnical Engineering
,
Prentice Hall
,
Englewood Cliffs, NJ
.
14.
Geotech Data
,
2014
, “Cohesion,” accessed Oct. 31, http://www.geotechdata.info/parameter/cohesion
15.
Geotechdata.info
,
2013
, “Soil Young's Modulus,” accessed Sept. 9, 2013, http://www.geotesting.info/parameter/soil-young's-modulus.html
16.
Renier Cloete
,
2018
, “Elastic Properties of Soils,” Prokon, Ireland, accessed Mar. 22, 2018, https://support.prokon.com/portal/kb/articles/elastic-properties-of-soils
17.
Blacker
,
T.
,
Owen
,
S. J.
,
Staten
,
M. L.
,
Quadros
,
R. W.
,
Hanks
,
B.
,
Clark
,
B.
,
Hensley
,
T.
,
Meyers
,
R. J.
,
Ernst
,
C.
,
Merkley
,
K.
,
Morris
,
R.
,
McBride
,
C.
,
Stimpson
,
C.
,
Plooster
,
M.
, and
Showman
,
S.
,
2017
, “
CUBIT Geometry and Mesh Generation Toolkit 15.3 User Documentation
,”
Sandia National Laboratories
,
Albuquerque, NM
, Report No. SAND2017-6895W.
18.
Skolnik Industries, Inc.
,
2006
, “
55 Gallon Open Head Drum, 1.5–1.5–1.5 CRCQ, 1A2/X430/S & 1A2/Y1.5/175
,” Drawing CQ5508, Rev A5,
Skolnik Industries
,
Chicago, IL
, accessed Jan., 2019, www.skolnik.com
19.
ASTM
,
2016
, “
Standard Specification for Steel, Sheet, Coled-Rolled, Carbon, Structural, High-Strength Low-Alloy, High-Strength Low-Alloy With Improved Formability, Solution Hardened, and Bake Hardenable
,”
ASTM
, West Conshohocken, PA, Standard No. A1008/A1008M-16.
20.
ASTM
,
2014
, “
Standard Specification for Carbon Steel Bolts, Studs, and Threaded Rod 60000 Psi Tensile Strength
,”
ASTM
, West Conshohocken, PA, Standard No. ASTM A307-14.
21.
Mangan
,
B. W.
,
2006
, “
How Much Do You Know About Forklift Safety?
,” Construction and Equipment, Seattle, WA, accessed Feb. 1, 2018, http://www.djc.com/news/co/11176945.html
22.
Skolnik Industries, Inc.
,
2009
, “
Exhibit 11.2A–Type A Test Report, Rev 03, ID 14-089, Part No. CQ5508
,”
Skolnik Industries
,
Chicago, IL
, accessed Jan., 2019, www.skolnik.com
23.
Hobbs
,
M. L.
,
2014
, “
Internal Memo: Cookoff Modeling of a WIPP Waste Drum (68660)
,”
Sandia National Laboratories
,
Albuquerque, NM
, Report No. SAND2014-20129O.
24.
Lorentz
,
E.
,
Cuvilliez
,
S.
, and
Kazymyrenko
,
K.
,
2011
, “
Convergence of a Gradient Damage Model Toward a Cohesive Zone Model
,”
C. R. Mec.
,
339
(
1
), pp.
20
26
.10.1016/j.crme.2010.10.010
25.
Zhou
,
F.
,
Molinari
,
J. F.
, and
Ramesh
,
K. T.
,
2005
, “
A Cohesive Model Based Fragmentation Analysis: Effects of Strain Rate and Initial Defects Distribution
,”
Int. J. Solids Struct.
,
42
(
18–19
), pp.
5181
5207
.10.1016/j.ijsolstr.2005.02.009
26.
Jardine
,
L. J.
,
Mecham
,
W. J.
,
Reedy
,
G. T.
, and
Steindler
,
M. J.
,
1982
, “
Final Report of Experimental Laboratory-Scale Brittle Fracture Studies of Glasses and Ceramics
,”
Argonne National Laboratory
,
Argonne, IL
, Report No. ANL-82-29.
27.
Wu
,
S. Z.
,
Chau
,
K. T.
, and
Yu
,
T. X.
,
2004
, “
Crushing and Fragmentation of Brittle Spheres Under Double Impact Test
,”
Powder Technol.
,
143–144
, pp.
41
55
.10.1016/j.powtec.2004.04.028
28.
Koch
,
W.
,
Lange
,
F.
,
Martens
,
R.
, and
Nolte
,
O.
,
2004
, “
Determination of Accident Related Release Data
,” 14th International Symposium on the Packaging and Transportation of Radioactive Materials (
PATRAM 2004
), Berlin, Germany, Sept. 20–24, p.
6
.https://inis.iaea.org/collection/NCLCollectionStore/_Public/37/088/37088755.pdf?r=1&r=1