Information Only 
2 

The simplified underground stratigraphic model used in geomechanical numerical analyses since 

the CCA (Butcher 1997, Stone 1997b, Park and Holland 2007) was modified to include the panel 

intake drift and the distance between adjacent panels (Figure 1). Other than those two changes, 

the mesh is changed as little as possible from the original mesh used by Stone (1997b) to 

minimize possible errors introduced by the changes. 

Figure 1. Simplified underground stratigraphic, half symmetry model and boundary conditions 

used in the JAS3D modeling of crushed salt panel closures. 

Information Only 
3 

In order to be used in JAS3D, the model had to be made 3D by adding a depth of 0.5 m. Even 

though it is a 3D model, a plane strain condition is assumed in the lengthwise direction of the 

panel closure. As shown in the verification and validation plan / validation document for JAS3D 

(WIPP PA 2009), modifying the simplified underground mesh to 3D has very little effect on 

identical underground analyses made for porosity surface calculations. The plane strain 

assumption for the crushed salt panel closure was shown to be valid based on a simplified 3D 

analysis of a model panel closure subjected to a constant stress loading condition (DOE 2012). 

A plane of symmetry through the center of the panel closure is used to limit the number of 

elements in the mesh by having to discretize only half of the model, a half symmetry model. This 

is accomplished by a no horizontal displacement boundary condition imposed on the left 

boundary of the mesh in Figure 1. The horizontal width of the model corresponds to half the 

distance between the drifts of adjacent panels, which is 230 ft (70.10 m). Use of another 

symmetry boundary condition on the right side of the mesh takes into account the influence an 

identical drift in an adjacent panel. On the top of the model is a traction boundary condition 

equivalent to the vertical in situ stress at that elevation. On the bottom of the model is a no 

displacement boundary condition in the vertical direction. 

The drift where the panel closure will be is initially completely filled with run-of-mine salt. This 

salt is modeled using the volumetric creep model for crushed salt. Contact surfaces are invoked 

between the crushed salt of the panel closure and the intact salt surrounding the drift so that the 

materials can deform independently and separate if necessary. In addition, the air pressure in the 

drift is assumed to be at one atmosphere and obey the ideal gas law. 

2 The Crushed Salt Constitutive Model 

The run-of-mine salt making up the panel closure is modeled using the volumetric creep model 

in JAS3D. The volumetric creep model is based on the work of Sjaardema and Krieg (1987). 

They developed their model based on the hydrostatic consolidation tests of crushed salt by 

Holcomb and Shields (1987). The following description of the model is based on those found in 

Stone (1997a) who describes how the model was implemented into various Sandia structural 

mechanics codes. The model was used because of its similarity to the crushed salt model in 

FLAC3D (Itasca, 2003), used by RockSol Consulting Group, Inc. for the panel closure design 

document. 

The total strain rate is decomposed into volumetric and deviatoric parts. Because intact salt 

creeps deviatorically when subjected to a deviatoric stress condition, crushed salt should 

logically be expected to also creep deviatorically. This expectation becomes more reasonable as 

the density increases and the crushed salt approaches an intact condition. The deviatoric salt 

creep model used in the volumetric creep model is based on the power-law creep model. The 

volumetric creep model is formulated such that the material response of the consolidated crushed 

salt becomes identical to that of the intact material as the density of the crushed salt approaches 

that of intact salt. The FLAC/FLAC3D crushed salt model used by RockSol differs from the 

Sandia model in that it appears to include primary creep as well as secondary creep in its 

implementation (Itasca, 2003). The solution technique is also different. 

Information Only 
4 

The bulk, K, and shear, G, elastic moduli were found from Holcomb and Shields tests to depend 

on the density, ρ, of the material through relationships of the form 

KK K  0 1 exp   (1) 

GG G  0 1 exp   (2) 

where K0, K1, G0, and G1 are material constants. 

The power-law creep model is cast as a conventional power-law secondary creep model of the 

form 

  exp

          m Q A

R (3) 

where  is the effective deviatoric stress, A and m are material constants,  is the absolute 

temperature, R is the universal gas constant, and Q is the activation energy. Using forward Euler 

integration and a creep strain rate for von Mises flow, the expression for the deviatoric stress rate 

is given by: 

3 1 2 exp 2

 

                

  n m n

ij ij ij

Q SG S

R (4) 

The volumetric creep model is developed by extension of a one dimensional problem into three 

dimensions. A uniaxially compressed porous crushed salt sample is envisioned as being 

composed of cylinders of salt, each of which has intact salt secondary creep behavior separated 

by areas of open space. The local stress acting on the salt cylinders is stated in terms of the 

average stress on the porous sample. The cross-sectional area of the porous sample is expressed 

in terms of the net cross-sectional area of the salt cylinders. The final resulting continuum model 

for the rate of the deviatoric stress of the crushed salt is given as: 

3 1 2 ( ) exp 2  



                      

 

m

n mn

ij ij ij

Q SG A S

R (5) 

where the constants A, Q, m, and ρ∞ refer to values for the intact salt, ρ is the density. 

The volumetric part of the model is divided into elastic and inelastic parts as 

  

   c

kk kk

p d d

K (6) 

Information Only 
5 

where dkk is the volumetric strain rate, 3

kk p



   is the rate of change of the pressure, c

kk d is the 

volumetric creep strain rate, and K(ρ) is the density-dependent bulk modulus. 

Laboratory consolidation tests have shown that the crushed salt volumetric creep rate can be fit 

by an expression of the form 

1

0

1 1 

       c A B p

kk d Be e (7) 

where B0, B1, and A are material constants. The density is computed from the relationship 

0

0   exp

         

t

kk

t

d dt (8) 

where 0 is the density at time t0. 

Equation (6) is solved for p and combined with the definition of the volumetric creep strain 

rate, Eq. (7), to yield 

    1

0

1 1  



       

 B p A

kk p K d Be e (9) 

which is the expression that is integrated. 

2.1 Input Parameters for the Crushed Salt Model 

The input parameters used in the volumetric creep model and their sources are given in Table 1. 

The initial densities, 0 , corresponding to fractional densities of 65, 67, 70, 75, 80, and 85%, are 

1404, 1440, 1512, 1620, 1728, and 1836 kg/m3. These represent, in order, initial porosities of 35, 

33, 30, 25, 20, and 15%. 

 

Information Only 
6 

Table 1. Input parameters used for the volumetric creep model panel closure. 

Parameter Value Units Source 

2G0 2.12e+04 Pa Sjaardema and Krieg (1987), p. 59 

K0 1.76e+04 Pa Sjaardema and Krieg (1987), p. 59 

D 5.79e-36 1/(Pan·s) Krieg (1984), p. 18 

n 4.9 — Krieg (1984), p. 18 

Q/RT 20.1422 — Krieg (1984), p. 18 

G1 6.53e-03 m3/kg Sjaardema and Krieg (1987), p. 59 

K1 6.53e-03 m3/kg Sjaardema and Krieg (1987), p. 59 

B0 1.30e+08 kg/(m3·s) Sjaardema and Krieg (1987), pp. 14-15 

B1 8.20e-07 1/Pa Sjaardema and Krieg (1987), pp. 14-15 

A -1.73e-02 m3/kg Sjaardema and Krieg (1987), p. 14 

ρ∞ 2160 kg/m3 Krieg (1984), p. 14; Brodsky (1994), p. 2 

ρ0 1404 * kg/m3 Initial fractional density of 65% of intact salt 

* The value ρ0 = 1404 kg/m3 is for an initial fractional density of 65% intact, which corresponds to a 

porosity of 35%. Initial density values at time t=0 vary depending on the predetermined porosity value 

being analyzed. The six initial porosities used were 15, 20, 25, 30, 33, and 35% (or fractional densities 

85, 80, 75, 70, 67, and 65%, respectively). 

3 Intact Salt Constitutive Model 

A combined transient-secondary creep constitutive model for rock salt known as the multimechanism deformation (M-D) model proposed by Munson and Dawson (1979, 1982, 1984) and 

extended by Munson et al. (1989), was used for the clean and argillaceous salt. The model is 

decomposed into an elastic volumetric part and a deviatoric part. 

The elastic volumetric part is defined by 

K

kk

kk 3



  (10) 

where, kk  = the total strain components 

 kk = the total stress components 

K = the elastic bulk modulus 

The deviatoric part defined by 

2

2 2

cos 2 3 sin 2 2

cos3 cos3 3 ij ij s ij ip pj ij s G e F s ss J

J J

   

 

                    

   (11) 

where 

3

kk

ijij s



  = the deviatoric stress 

Information Only 
7 

G = the elastic shear modulus 

3

kk

ijij e



  = the deviatoric strain 

1for

0 for ij

i j

i j

    

 

 = Kronecker delta 

J2 and θ are the second invariant of the deviator stress and the Lode angle, respectively, and will 

be defined later. 

The second term of Equation (11) represents the creep contribution. F is a multiplier on the 

steady-state creep rate to simulate the transient creep response according to the following 

equation 























]/1[ *

*

,]/1[ *

2*

2*

1

t

t

t

t

t

e

e

F











 (12) 

where  work-hardening parameter 

  recovery parameter 

*

t  transient strain limit 

 is an internal state variable whose rate of change is determined by the following evolutionary 

equation 

F s     )1( (13) 

In Equation (12), the work-hardening parameter  is defined as 

      G)/log( (14) 

where,  and  are constants. The variable  is the equivalent Tresca stress given by 

  J 2 cos2  (15) 

where 



 



  2/3

2

3

)(2

33 arcsin

3

1

J

J  is the Lode angle limited to the range ( 66

    ). 

qppq ssJ

2

1

2  = second invariant of the stress deviator 

rpqrpq sssJ

3

1

3  = third invariant of the stress deviator 

Information Only 
8 

The recovery parameter,  , is held constant. The transient strain limit is given by 

cT M

t GeK )/( 0

*    (16) 

where K0 , c, and M are constants. 

The steady-state, or secondary creep strain rate, s  , is given by 

1 2

12 1 2 / / // 0)

1 2 12

( sinh

n n

Q R Q R QR Q R

s

q Ae Ae H Be B e

GG G

   

                           

 (17) 

Where Ai s and Bi s are constants 

Qis are activation energies 

 = the absolute temperature 

R = the universal gas constant 

ni s are the stress exponents 

q = the stress constant 

 0 = the stress limit of the dislocation slip mechanism 

H = the Heaviside step function with the argument )(  0

3.1 Input Parameters for the Intact Salt Model 

The material constants used in the analyses for the clean and argillaceous salt are given in Table 

2 and Table 3. 

Table 2. Salt elastic properties (Butcher 1997) 

Parameter Value Units 

G 

E 



12,400 

31,000 

0.25 

MPa 

MPa 

– 

 

Information Only 
9 

Table 3: Salt creep properties (Munson et al. 1989) 

Parameter Clean Salt Argillaceous Salt Units 

A1 8.386e+22 1.407e+23 1/sec 

Q1 25,000 25,000 cal/mole 

n1 5.5 5.5 n1

B1 6.086e+06 8.998e+06 1/sec 

A2 9.672e+12 1.314e+13 1/sec 

Q2 10,000 10,000 cal/mole 

n2 5.0 5.0 –

B2 3.034e-02 4.289e-02 1/sec 

 0 20.57 20.57 MPa 

q 5335 5335 –

m 3.0 3.0 –

K0 6.275e+05 2.470e+06 –

c 9.198e-03 9.198e-03 1/T 

 -17.37 -14.96 – 

 -7.738 -7.738 – 

 0.58 0.58 – 

4 Run Control 

Herrick 2012 contains the run control information for the JAS3D analyses in which the emplaced 

panel closures use a crushed salt fill with initial porosities of 33, 25, 20, and 15% (or initial 

fractional densities of 67, 75, 80, and 85%, respectively) which were presented in the DOE’s 

response to EPA questions (DOE 2012). This section contains the run control information for the 

additional analyses in which the emplaced panel closures use crushed salt fill with initial 

porosities of 35 and 30% (or initial fractional densities of 65 and 70%, respectively). Of perhaps 

greatest importance is that the results module Jas3d/Output is the same for all analyses including 

those used in DOE (2012). 

For the additional two analyses, JAS3D (version 2.4.C-WIPP) and its associated utilities 

(Algebra2 (version 1.28), Aprepro (version 2.06), BlotII2_cps (version 1.59A), FASTQ_cps 

(version 3.18), and GEN3DII (version 1.37a)) were run on April 20, 2012 on seals.sandia.gov 

using the Python run control scripts PC3R_3.py and PC3R_4.py and the Python module rc.py. 

The Postscript files generated by BlotII2 (version 1.59A) were converted to PDF files using the 

UNIX utility ps2pdf. All input files for the runs were exported from the CVS repository 

$CVSLIB/Analyses/PC3R, module Jas3d/Input/Underground2. 

The run control scripts and an associated C shell scripts, RunJas3d_3.csh and RunJas3d_4.csh, 

were exported from the RunControl module of the same repository. JAS3D and associated utility 

codes and scripts were exported from code-specific CVS repositories (Table 4). The output files 

from JAS3D, Algebra, BlotII2, and ps2pdf were stored in the Output module of the repository 

$CVSLIB/Analyses/PC3R. Additional data files (.dat) produced by BlotII2 were added to the 

Jas3d/Output module of the repository after the PC3R.py script completed. The log files 

Information Only 
10 

PC3R_3.log and PC3R_4.log, and RTF versions of those files, are stored in the module 

RunControl/Log. All files were tagged with the CVS tag PC3R. 

Table 4. Run control code-specific CVS repositories. 

File Repository Comment

Builds/Linux64/algebra2 

(Ver:1.28) 

/nfs/data/CVSLIB/ACCESS_WIPP/ALGEBRA2 Manipulates EXODUS finite 

element output data by 

evaluating expressions 

Builds/Linux64/aprepro 

(Ver:2.06) 

/nfs/data/CVSLIB/ACCESS_WIPP/APREPRO An algebraic preprocessing 

program 

Builds/Linux64/blotII2_cps 

(Ver:1.59A) 

/nfs/data/CVSLIB/ACCESS_WIPP/BLOTII2 Postscript driver for blot 

Builds/Linux64/fastq_cps 

(Ver:3.18) 

/nfs/data/CVSLIB/ACCESS_WIPP/FASTQ Finite element 2-dimensional 

mesh generation program 

Builds/Linux64/gen3d2 

(Ver:1.37a) 

/nfs/data/CVSLIB/ACCESS_WIPP/GEN3DII Transforms a 2-dimensional 

GENESIS database into a 3-

dimensional GENESIS database 

Builds/Linux64/jas3d 

(Ver:2.4.C-WIPP) 

/nfs/data/CVSLIB/ACCESS_WIPP/JAS3D Finite element code for 3-

dimensional analyses 

RunControl/algebra2.sh /nfs/data/CVSLIB/ACCESS_WIPP/ALGEBRA2 Script for manipulating EXODUS 

finite element output data by 

evaluating expressions 

RunControl/blotII2.sh /nfs/data/CVSLIB/ACCESS_WIPP/BLOTII2 Script for running graphics code 

for examining meshes 

RunControl/fastq.sh /nfs/data/CVSLIB/ACCESS_WIPP/FASTQ Script for running finite element 

2-dimensional mesh generation 

program 

RunControl/gen3dII.sh /nfs/data/CVSLIB/ACCESS_WIPP/GEN3DII Script to transform a 2-

dimensional GENESIS database 

into a 3-dimensional GENESIS 

database 

RunControl/jas3d.sh /nfs/data/CVSLIB/ACCESS_WIPP/JAS3D Script for running jas3d.exe 

5 Results 

The following figures, from the PostScript files (Jas3d/Output/*.ps), and associated data files 

(Jas3d/Output/*.dat), where * indicates the base filename, were used in support of parameter 

value selections for PCS-2012 PA (Camphouse 2012). Base filenaming for the panel closure 

analyses takes one of two forms: “pclo_vcXX_por-void” or “pclo_vcXX_SIG-DY”. Here “plco” 

signifies that the analyses were of panel closures with the underground stratigraphy included, 

“vcXX” indicates that the model used for the run-of-mine salt was the volumetric creep model 

beginning with an initial fractional density of XX percent, “por-void” means that the analyses 

considered porosity and void volume properties of the crushed salt, and “SIG-DY” means that 

the analyses looked at stresses and vertical displacements of the back of the drift or top of the 

crushed salt. For the latest two analyses, the filename “vc65” indicates that the run-of-mine salt 

was modeled with the volumetric creep model having an initial fractional density of 65%, and 

filename “vc70” indicates that the run-of-mine salt was modeled as having an initial fractional 

density of 70%. The same filenaming convention is used for the previous set of analyses where 

the initial fractional densities 67, 75, 80, and 85% (Herrick 2012). 

Information Only 
11 

5.1 Porosity History Plots 

Starting with the smallest initial density (fractional density = 65%) and progressing to the highest 

(fractional density = 85%), the porosity histories through the first 500 years for all JAS3D 

underground analyses are given in Figure 2 through Figure 7. 

Figure 2. Porosity history through 500 years for a crushed salt panel closure in the intake drift 

having an initial porosity of 35%, which corresponds to an initial fractional density of 65%. 

Information Only 
12 

Figure 3. Porosity history through 500 years for a crushed salt panel closure in the intake drift 

having an initial porosity of 33%, which corresponds to an initial fractional density of 67%. 

Figure 4. Porosity history through 500 years for a crushed salt panel closure in the intake drift 

having an initial porosity of 30%, which corresponds to an initial fractional density of 70%. 

Information Only 
13 

Figure 5. Porosity history through 500 years for a crushed salt panel closure in the intake drift 

having an initial porosity of 25%, which corresponds to an initial fractional density of 75%. 

Figure 6. Porosity history through 500 years for a crushed salt panel closure in the intake drift 

having an initial porosity of 20%, which corresponds to an initial fractional density of 80%. 

Information Only 
14 

Figure 7. Porosity history through 500 years for a crushed salt panel closure in the intake drift 

having an initial porosity of 15%, which corresponds to an initial fractional density of 85%. 

As a check of the consistency of the results obtained from JAS3D by the method described 

herein, a plot of the time at which the crushed salt in the panel closure consolidates to a porosity 

of 5% for the six different initial porosities is shown in Figure 8. The high R2 value demonstrates 

the consistency of the analytical method. The relationship was fit by a linear trend line using 

Excel. 

Information Only 
15 

Figure 8. Plot of the time it takes for the crushed salt fill in the panel closures to reach 5% 

porosity starting from the six different initial porosities used in the JAS3D analyses. 

5.2 Backstress in the Back (Roof) of the Drift History Plots 

The backstresses generated by the crushed salt in the back (roof) of the intake drift resulting from 

consolidation of the crushed salt filling a panel closure was estimated graphically. The next 

series of plots, Figure 9 through Figure 14, are histories of the vertical stress in the back of the 

drift again through the first 500 years. The location of the element used to obtain the data is in 

the center of the drift. 

Information Only 
16 

Figure 9. Vertical stress history through 500 years in the back (roof) of the intake drift produced 

by a run-of-mine salt panel closure having an initial porosity of 35%. 

Figure 10. Vertical stress history through 500 years in the back (roof) of the intake drift produced 

by a run-of-mine salt panel closure having an initial porosity of 33%. 

Information Only 
17 

Figure 11. Vertical stress history through 500 years in the back (roof) of the intake drift produced 

by a run-of-mine salt panel closure having an initial porosity of 30%. 

Figure 12. Vertical stress history through 500 years in the back (roof) of the intake drift produced 

by a run-of-mine salt panel closure having an initial porosity of 25%. 

Information Only 
18 

Figure 13. Vertical stress history through 500 years in the back (roof) of the intake drift produced 

by a run-of-mine salt panel closure having an initial porosity of 20%. 

Figure 14. Vertical stress history through 500 years in the back (roof) of the intake drift produced 

by a run-of-mine salt panel closure having an initial porosity of 15%. 

Information Only 
19 

6 References: 

Brodsky, N.S. 1994. Determination of the Intact Density of WIPP Salt Specimens. Calculation 

No. 325/03/02. Respec Inc., Rapid City, SD 

Butcher, B.M. 1997. A Summary of the Sources of Input Parameter Values for the WIPP Final 

Porosity Surface Calculations. SAND97-0796. Sandia National Laboratories, Albuquerque, NM 

Camphouse, R.C., M. Gross, C.G. Herrick, D.C. Kicker, B. Thompson. 2012. Recommendations 

and Justifications of Parameter Values for the Run-of-Mine Salt Panel Closure System Design 

Modeled in the PCS-2012 PA. Sandia National Laboratories, Carlsbad, NM 

DOE (U.S. Department of Energy). 2012. Response to First Set of EPA Questions to DOE, 

December 22, 2011. Letter from Mr. Juan Franco, Manager, Carlsbad Area Field Office, U.S. 

Department of Energy, to Mr. Jonathan Edwards, Director, Radiation Protection Division, U.S. 

Environmental Protection Agency, dated April 17, 2012. 

Herrick, C.G. 2012. Calculations Performed in Support of Reconsolidation of Crushed Salt in 

Panel Closures. Memorandum to the WIPP Records Center dated March 29, 2012. ERMS 

557150. Sandia National Laboratories, Carlsbad, NM. 

Holcomb, D.J. and M.F. Shields. 1987. Hydrostatic Consolidation of Crushed Salt with Added 

Water. SAND87-1990., Sandia National Laboratories, Albuquerque, NM 

Itasca Consulting Group. 2003. FLAC3D: Fast Lagrangian Analysis of Continua in 3 

Dimensions User’s Guide. Version 2.1, Optional Features, Section 2 – Creep Materials Models. 

Minneapolis, MN 

Krieg, R. D. 1984. Reference Stratigraphy and Rock Properties for the Waste Isolation Pilot 

Plant (WIPP) Project. SAND83-1908. Sandia National Laboratories, Albuquerque, NM 

Munson, D.E. and P.R. Dawson. 1979. Constitutive Model for the Low Temperature Creep of 

Salt (with Application to WIPP). SAND79-1853. Sandia National Laboratories, Albuquerque, 

NM 

Munson, D.E. and P.R. Dawson. 1982, A Transient Creep Model for Salt during Stress Loading 

and Unloading. SAND82-0962. Sandia National Laboratories, Albuquerque, NM 

Munson, D.E. and P.R. Dawson. 1984. Salt Constitutive Modeling Using Mechanism Maps. The 

Mechanical Behavior of Salt, Proceedings of the First Conference. Pennsylvania State 

University, University Park, PA. November 9-11, 1981. H.R. Hardy, Jr. and M. Langer, eds. Karl 

Distributors, Rockport, MA. pp. 717-737. 

Munson, D.E., A.F. Fossum, and P.E. Senseny. 1989. Advances in Resolution of Discrepancies 

between Predicted and Measured in Situ Room Closures. SAND88-2948. Sandia National 

Laboratories, Albuquerque, NM 

Information Only 
20 

Park, B.Y. and J.F. Holland. 2007. Structural Evaluation of WIPP Disposal Room Raised to Clay 

Seam G. SAND2007-3334. Sandia National Laboratories, Albuquerque, NM 

Sjaardema, G.D. and R.D. Krieg. 1987. A Constitutive Model for the Consolidation of WIPP 

Crushed Salt and Its Use in Analyses of Backfilled Shafts and Drift Configurations. SAND87-

1977. Sandia National Laboratories, Albuquerque, NM 

Stone, C.M. 1997a. SANTOS – A Two-Dimensional Finite Element Program for the Quasistatic, 

Large Deformation, Inelastic Response of Solids. SAND90-0543. Sandia National Laboratories, 

Albuquerque, NM 

Stone, C.M. 1997b. Final Disposal Room Structural Response Calculations. SAND97-0795. 

Sandia National Laboratories, Albuquerque, NM 

WIPP PA. 2009. Verification and Validation Plan / Validation Document for JAS3D Version 

2.4.C-WIPP, document version 2.4.C-WIPP. ERMS 545606. Sandia National Laboratories, 

Carlsbad, NM. 

Information Only