ML20056D845

Forwards CNWRA 93-005, Evaluation of Coupled Computer Codes for Compliance Determination

ML20056D845
Person / Time
Issue date:	08/06/1993
From:	Joseph Holonich NRC OFFICE OF NUCLEAR MATERIAL SAFETY & SAFEGUARDS (NMSS)
To:	Shelor D ENERGY, DEPT. OF
References
REF-WM-11 HLWR, NUDOCS 9308180144
Download: ML20056D845 (1)

Text

A9G 0 6 993 1

g 0010 Mr. Dwight E. Shelor, Associate Director I

for Systems and Co'mpliance Office of Civilian Radioactive Waste Management U. S. Department of Energy 1000 Independence, S. W.

Washington, D. C.

20585

Dear Mr. Shelor:

i Enclosed for your information is the Center for Nuclear Waste Regulatory Analyses (CNWRA) report entitled " Evaluation of Coupled Computer Codes for i

Compliance Determination." The enclosed report was prepared to document work i

performed by the CNWRA for the Nuclear Regulatory Commission. The Department of Energy may find the information in this report useful.

If you have any questions, please feel free to contact Charlotte Abrams of my staff at (301) 504-3403.

Sincgrely, Joseph J. Holonich, Director Repository Licensing and Quality Assurance Project Directorate Division of High-Level Waste Management Office of Nuclear Material Safety and Safeguards Enclosure.

As stated cc:

R. Loux, State of Nevada T. J. Hickey, Nevada Legislative Committee C. Gertz, DOE /NV I

M. Murphy, Nye County, NV M. Baughman, Lincoln County, NV D. Bechtel, Clark County, NV 3

D. Weigel, GA0 P. Niedzielski-Eichner, Nye County. NV B. Mettam, Inyo County, CA V. Poe, Mineral County, NV F. Sperry, White Pine County, NV R. Williams, Lander County, NV L. Fiorenzi, Eureka County, NV i

J. Hoffman, Esmeralda County, NV C. Schank, Churchill County, NV L. Bradshaw, Nye County, NV DISTRIBUTION CNWRA NMSS R/F HLPD R/F LSS LPDR ACNW PDR CENTRAL FILE BJYoungblood, HLWM JLinehan, HLWM RBallard, HLGE s.MFederline, HLHP b

HLP610 b

OFC HLPD b

HLPD 1/

NAME AGarci [/ @

CAhams Jhhch DATE 081 h /93 08/b/93 08/ 7/93 g

C = C0VER' E = COVER & ENCLOSURE N = NO COPY Sl_\\93-005 0FFICIAL RECORD COPY I

4%f

/

~

j$

r

\\h 93081B01'4d 330806 l]

PDR WASTE 7

WM-11 ppg

q

~

CNWR A 93-005 t

mvm mme

~ mywgmmmemmmergwmngemmwmmmmmm.x w a m a u a a;g s w. auw m sg.

'u 4s.

?h m a, a k s u:g h;i x agagugg;gav m

A A

h h

e e

e t

e e e-D A

8 3

r.

. msL he Mad &dadheML ~

~

m 1

l

'M9

1. $N f7$ e.~t',p ' & -'". JT%Y iE

[

D -

M DN N'.

,?

N I

't x

s e.u c n a a s a %2u.w a.ucamanua,cua;.a.nwa;ua..wa..wa u.2m u ua:;an w u sea m m2wasm 0

m

/*

a R

= reamed ic-Nuclear Regulatory Commission 3

2 Contract NRC-02-88-005 i8 e,

<,a a

c e

n4 B1 Center for Nuclear Waste Regulatory Analyses 3

San Antonio. Texas n

June 1993 l

y a

IP d

P y - =.,.-

. s. -;...

4 h t j,,

I CMWRA 934DS EVALUATION OF COUPLED COMPUTER CODES FOR COMPLlANCE DETERMINATION 4

l li

' I I

l Prepared for i

i Nuclear Regulatory Commission I

Contract NRC-02-88-005 I

I Prepared by Center for Nuclear Waste Regulatory Analyses l

San Antonio, Texas I

June 1993 g

g g

~ :::

CNWRA 93-005 I

I EVALUATION OF COUPLED COMPUTER CODES FOR COMPLIANCE DETERMINATION l

I Prepared for jl Nuclear Regulatory Commission i

Contract NRC-02-88-005 lI Prepared by Amitava Ghosh l

Sui-Min (Simon) Hsiung Mikko P. Ahola l

Asadul H. Chowdhury i

I Center for Nuclear Waste Regulatory Analyses San Antonio, Texas I

.g June 1993 g

I lI l

ABSTRACT The objective of this study is to select a computer code for simulating coupled processes that can be used by the U.S. Nuclear Regulatory Commission (NRC) and the Center for Nuclear Waste Regulatory Analyses (CNWRA), with necessary modifications, for determination of U.S. Department of Energy

{

(DOE) compliance with NRC regulations on thermal loads. This study will be conducted in two phases; only the first phase has been completed and is reported herein. In Phase I of the code selection study, fifteen computer codes were evaluated. Each of these codes has, at a minimum, the capability to simulate

[

thermomechanical effects. Three major evaluation criteria were developed as bases for the evaluation of L

the codes: number of coupled processes implemented, capability of including rock joints explicitly, and capability of simulating responses due to seismic loading. A subjective numerical scoring system was used p

for this evaluation. The " coupled processes" criterion has been broken down further, based on the nature L

of coupling (e.g., thermal to mechanical) included in the code. Two-way coupling between the mechanical processes and th fluid flow through joints has been deemed most important in the near-field environment p

of the propouxi nuclear waste repository at Yucca hiountain, Nevada. Consequently, it has received L

double weighting in the scoring and ranking analyses. Four codes - ROChiAS, ABAQUS, FEHhiS, and UDEC - received a score of 6 or more, out of a maximum possible score of 9. Neither FEHhiS nor ROCMAS is publicly available for further evaluation. The distinct element code UDEC and the finite

(

element code ABAQUS were selected for further evaluation in Phase II of the code selection study. One of these two codes will be selected for further development for use as a compliance determination code.

This code will be suitable for modeling thermal, mechanical, and hydrological processes. If necessary, the chemical processes will be incorporated at a later stage of code development.

J P

L t

{

[

CONTENTS Section Page TABLES iv ACKNOWLEDGEMENTS v

EXECUTIVE

SUMMARY

. vi 1

INTRODUCTION 1-1 lI 1.1 GENERAL

. 1-1 1.2 OBJECTIVES 1-2 1.3 EXISTING SURVEYS 1-2 l

2 GENERAL DESCRIPTION OF THE COUPLED PHENOMENA AT 2.1 GENE DESC Tibh OF THE SITE AT YUCC MOUNTAIN

2.2 DESCRIPTION

OF THE COUPLING PHENOMENA

. 2-1 2.2.1 Thermal to Mechanical Coupling

. 2-2 2.2.2 Mechanical to Thermal Coupling 2-3 2.2.3 Mechanical to Hydrological Coupling 2-3 2.2.4 Hydrological to Mechanical Coupling.

. 2-3 2.2.5 Thermal to Hydrological Coupling 2-4 2.2.6 Hydrological to Thermal Coupling 2-4 l

3 IDENTIFICATICN OF COMPUTER CODES 3-1 i

4 SELECTION CRITERIA.

. 4-1

]

4.1 COUPLED PROCESSES MODELED 4-1 4.2 JOINT CONSTITUTIVE LAWS 4-1 i

4.3 DYNAMIC AND SEISMIC CAPABILITIES

. 4-1 4.4 TWO-PHASE FLOW 4-2 y

4.5 SOURCE CODE 4-2 i"

5 CODE RANKING.

. 5-1 5.1 APPROACH OF CODE RANKING 5-1

~

5.2 RANKING OF COMPUTER CODES 5-1 6

DISCUSSION AND RECOMMENDATIONS 6-1 P

7 REFERENCES D

7-1 l

APPENDIX Additional Information on the Computer Codes

]

L E

111 s

N

L TABLES Table.

Page a

List of candidate coupled codes and their function capabilities.

4-3 4-1 5-1 Ranking analysis of the computer codes 5-2 5-2 Ranks of the computer codes.

5-3

[

r L

3 M

L f

)

d iV F

s s

I.

ACKNOWLEDGEMENTS I'

The authors thank Dr. Wesley C. Patrick for his technical review of this document. Thanks also go to Dr. Budhi Sagar for his programmatic review and Drs. Mysore Nataraja, Banad Jagannath, Shiann-Jang Chern, and William Boyle of the U.S. Nuclear Regulatory Commission (NRC) for their valuable I

suggestions during the evaluation of the codes and preparation of this document. The authors are also thankful to Rebecca A. Sanchez for skillful typing and formatting of the report, and to Mr. James W.

Pryor and Curtis Gray who provided a full range of expert editorial services in the preparation of the final document.

This report was prepared to document work performed by the Center for Nuclear Waste Regulatory Analyses (CNWRA) for the NRC under Contract NRC-02-88-005. The activities reported here were I

performed on behalf of the NRC Office of Nuclear Material Safety and Safeguards (NMSS), Division of High-Level Waste Management (DHLWM). This report is an independent product of the CNWRA and does not necessarily reflect the views or regulatory position of the NRC.

lI l

5 1

i lI lh J

e 1

L E

v y

M

I EXECUTIVE

SUMMARY

I The design and performance assessment of the proposed high-level nuclear waste (HLW) repository at Yucca Mountain, Nevada, likely will be based on numerical simulations of several coupled processes.

These simulations may consider thermal (T), mechanical (M), hydrological (H), and chemical (C)

~

processes and their interactions (TMHC) over a variety of spatial scales (e.g., individual excavation such as from the drifts to the whole repository) and temporal scales (e.g., O to 10,000 years). The U.S.

~

Nuclear Regulatory Commission (NRC) staff technical position on " Geologic Repository Operations Area Underground Facility Design - Thermal Loads" (Nataraja and Brandshaug,1992), outlined an acceptable L

method for U.S. Department of Energy (DOE) demonstration of compliance with the NRC regulations on thermal loads, as specified in 10 CFR 60.133(i). The objective of this study is to select a computer n

code that can be used by the NRC and the Center for Nuclear Waste Regulatory Analyses (CNWRA),

with necessary modifications, for determination of DOE compliance with NRC regulations on thermal l

loads. During Phase I of the code selection study reported herein, fifteen computer codes were evaluated.

Each of these codes has, at a minimum, the capability to simulate thermomechanical effects.

m 4

The fifteen codes were evaluated based on predetermined sets of criteria - number of coupled processes included, rock joint model, and capability of simulating responses due to dynamic and seismic loads. A rating system with these criteria was devised for comparing the capability of these codes in the context p

of reviewing design and performance of the repository. Among the coupled processes, two-way coupling i

between the mechanical process and the flow of fluid through rock joints (M*H ) has been considered j

very imponant in the near-field environment at the proposed nuclear waste repository at Yucca Mountain, Nevada. The symbol H; denotes joint flow. As a result, M*H received a double weight of 2. The j

L coupled model for simulating the coupling between the thermal (T) and matrix flow (H ) is under m

development by the Performance Assessment Element at the CNWRA. This model will be adopted and H

coded in the computer code selected. As a result, the coupled processes between thermal (T) and matrix l

flow (H ) were not included as criteria in this study. The M-T coupling was not included either since m

it was judged not important (Manteufel et al.,1993). Each of the remaining coupled processes under consideration in this study, namely T-M, T-H,, H -T, and M~H received an equal weight of 1. Of j

m l

these latter couplings T-H also been thought important with NRC's rule on thermal loads, due to the impact of repository-heatdriven hydrothermal flow on fracture-dominated flow (Buscheck and Nitao, 1993). Four codes - ROCMAS, ABAQUS, FEHMS, and UDEC - received a score of 6 or more, from l

a maximum possible score of 9. Although FEHMS and ROC.'1AS have been reported in the open literature, these codes are not in the public domain and are not available for further evaluation by the CNWRA. As a result, these tw o codes cannot be chosen at this time for further consideration in the Phase

[

11 study.

3 The distinct element code UDEC and the finite element code ABAQUS were selected for further evaluation during Phase II of the code selection study. UDEC has undergone some qualification study at the CNWRA (Brady et al.,1990). ABAQUS will go through a similar qualification program during the Phase Il study. One of these two codes will then be selected for further development for use as a F

compliance determination code. After appropriate modification, the selected code will be suitable for L

modeling thermal, mechanical, and hydrological processes to the extent needed to determine compliance with NRC regulations on thermal loads.

I L

b VI H

Section 1 of this report briefly discusses the code selection study and states the objectives. It also gives a list of publications from which the candidate codes were selected for evaluation in this study. These publications have extensive lists of codes capable of modeling structure in geological media. Section 2 describes the problems anticipated to be encountered in design and assessment of performance of the

(

proposed nuclear waste repository. It also briefly discusses the coupling among thermal, mechanical, and hydrological processes. The rationale for identification of the candidate codes is provided in Section 3.

The selection criteria are given in Section 4. Capabilities of each of the codes corresponding to the

[

individual selection criterion also are given. A ranking system and the rank of the codes are given in L

Section 5. Recommendations are given in Section 6. Additiont.1 information about each of the computer codes is given in the Appendix. Ahhough the information in the Appendix was not directly used in the

[

selection process, it makes the discussion of each code more complete.

t l

l

}

l I

t vii I

L

1 INTRODUCTION i

i 1.1 GENERAL I

The response of a rock mass in a high-level nuclear waste (HLW) geologic repository is a coupled phenomenon involving thermal (T), mechanical (M), hydrological (H), and chemical (C)

I processes. Coupled processes imply that one process affects another and that rock mass response in a I

repository environment cannot be predicted by considering each process independently. The importance of various processes will depend upon the thermal loading of the repository, the design of the engineered l

barriers, properties of the geologic medium, the time and spatial scales at which these processes are of interest, and the measure selected to characterize response.

l There is a general consensus in the literature concerning the importance of considering coupled processes at the proposed nuclear waste repository at Yucca Mountain (Manteufel et al.,1993). A recent I

U.S. Nuclear Regulatory Commission (NRC) staff technical position provided an acceptable methodology l

for systematically considering thermal loads and thermally induced mechanical, hydrological, and chemical processes (Nataraja and Brandshaug,1992). The primary purpose of the NRC technical position I

is to outline an acceptable method of comprehensively, systematically, and logically understanding and evaluating TMHC responses for the design and performance assessment of a geologic repository and for demonstration of compliance by the U.S. Department of Energy (DOE) with the NRC regulations on thermal loads (U.S. NRC,1992).

The spatial scale of interest in this code selection stucy is the near-field which includes both emplacement borehole and emplacement drift scales. The near-field complex environmental conditions at the repository horizon include mechanically disturbedjointed roc k, elevated temperature, and thermally induced mechanical, hydrological, and chemical processes - in.:luding phase changes of groundwater.

Ground motions due to earthquakes, underground weapons effect testing, etc., are superimposed on the in siru stresses, thermal loads, and thermally induced phenomena in a repository. It is necessary to have appropriate concep:ual models and associate computer codes describing the coupled phenomena to effectively determine compliance with various regulations.

The ultimate objective of this study is to select a computer code for NRC use in determining DOE compliance with NRC regulations on thermal loads. During Phase I of this code selection study, which is presented in this report, a number of coupled computer codes were evaluated based on reported l

information. Only two of these codes were selected for Phase II evaluation and one code eventually will be selected that could be used, possibly with some modifications, for determination of DOE compliance with NRC regulations on thermal loads.

l L

This study concentrates on the interactions among the thermal, mechanical, and hydrological l

processes relevant to HLW disposal at Yucca Mountain. The effect of chemical reactions on these l

processes and vice versa were not considered in this report, primarily because (i) no code was available that can deal with all 4 processes; and (ii) chemical processes are not considered to be important for design purposes. The final code which eventually will be selected in Phase II should have the capabilities, at a minimum, to simulate the thermomechanical responses of the fractured rock mass. Other capabilities that will be considered in code selection include the ability of the code to model stress-dependent fracture

{

and matrix flow and vice versa. In the near-field, modeling of flow through the fracture network is I

1-1

considered to be important (Buscheck and Nitao,1993; Norris,1989). However, calculation of flow through rock joints or fractures is complex and remains a subject of further research. Nevertheless, the code should have some capabilities to model the fracture flow even ifin a simplistic way. Matrix flow depends on the spatial distribution of porosity in the rock mass, on the fracture density and interconnectivity, and on fluid properties. Modeling matrix flow is a desirable feature in the selected code. 'Ite temperature distribution in the rock mass affects both the fluid flow and the mechanical stresses. Determination of thermal stress distribution in the medium and modeling of thermally induced

{

fluid flow are also desirable features in the selected code. The code should also be able to accept the actual time history of seismic load for calculating the response of the rock mass. Modeling of natural discontinuities is considered essential for the selected code. The code should be able to properly model the displacement discontinuity associated with the fracture and to redistribute the stresses accordingly in the surrounding media. It is understood that uncertainties in the above parameters will affect the ability of the selected computer code to predict long term repository performance under TMH processes. All of the codes considered in this study are purely deterministic in nature. Thus, uncertainty and sensitivity will likely be addressed through either parametric studies or possibly inclusion of Monte Carlo simulation capabilities within the code.

Fifteen computer codes were considered in this phase of the study. Each of these codes has, at a minimum, the capability to simulate thermomechanical processes. Most of the information given in this report comes from users' manuals of the codes and published reports. Wherever available, the personal experiences of the authors are also included in the evaluation of a particular code.

1.2 OBJECTIVES The objectives of this study are to:

(i)

Compile a list of computer codes that can simulate the coupled phenomena among TMH processes to be encountered at the proposed HLW disposal site at Yucca Mountain, Nevada. Each candidate computer code, at a minimum, should be able to simulate thermomechanical processes.

(ii)

Establish a set of criteria for code evaluation for simulating the coupled effects among TMH processes.

(iii) Produce a shon list of codes from the list developed in the first objective based on the evaluation criteria for a next level of detail evaluation (Phase II) to select a code for i

further development.

1.3 EXISTING SURVEYS L

As a starting point for this study, the following documents were reviewed to accomplish the first objective listed in Section 1.2:

(i)

Site Characterization Plan, Yucca Mountain Site, Nevada Research and Development Area, Volume VII, Part B, Chapter 8, Section 8.3.5.19.1 (Tables 8.3.5.19-1 and 8.3.5.19-2) (U.S. DOE.1988).

1-2

I (ii)

A Literature Review of Coupled Thermal-Hydrologic-Mechanical-Chemical Processes Peninent to the Proposed High-Level Nuclear Waste Repository at Yucca Mountain, NUREG/CR-6021 (Manteufel et al.,1993).

(iii) Rock Engineering Software, International Journal of Rock Mechanics and Mining Sciences & Geomedianics Abstracts, Vol. 25, No. 4. (International Society for Rock

.I Mechanics,1988).

(iv) Computer Aided Engineering Systems Handbook, Computational Mechanics Publications I

and Springer-Verlag, (Puig-Pey and Brebbia,1987).

(v)

Critical Assessment of Seismic and Geomechanics Literature Related to a High-Level Nuclear Waste Underground Repository, NUREG/CR-5440 (Kana et al.,1991)

In addition to the above publications, the personal knowledge of the authors was used in generating the list of candidate codes.

E s

5 I

E i.3 4

2 GENERAL DESCRIPTION OF THE COUPLED PHENOMENA AT YUCCA MOUNTAIN I

2.1 GENERAL DESCRIPTION OF TIIE SITE AT YUCCA MOUNTAIN A general description of the Yucca Mountain site for the proposed HLW repository is presented in the DOE Site Characterization Plan (U.S. DOE,1988). Yucca Mountain is located in southern Nevada about 160 km by road northwest of Las Vegas. The area is characterized by nonh to northwest-trending mountain ranges, composed of volcanic and volcan9 clastic strata that dip eastward. These strata are broken into en-echelon fault blocks. Arid climate prevails in the Yucca Mountain area with less than 25 l

em of rain per year. No perennial streams exist in the general vicinity.

The Paintbrush Tuff in the Yucca Mountain area is more than 460 m thick. It makes up nearly l

all of the strata exposed at Yucca Mountain (U.S. DOE,1988). The four members of the Paintbrush Tuff I

in ascending order are: Topopah Spring, Pah Canyon, Yucca Mountain, and Tiva Canyon. The Topopah Spring and Tiva Canyon members consist predominantly of devitrified and densely welded ash-flow tuffs.

They enclose between them nonwelded to partially welded and ash-flow tuff and tuffaceous sediments of the Pah Canyon and the Yucca Mountain members. The densely welded, devitrified part of the Topopah Spring tuff is currently considered as the likely emplacement horizon. The Topopah Spring tuff has lateral continuity and is above a major aquitard (Wahmonie formation) (U.S. DOE,1988). The fractures at the Yucca Mountain site are mostly stratabound. The preferential pathways for migration of radionuclides depend on complex intersections of fracture networks at strata boundaries as well as on some vertically continuous fractures. Cooling fractures of the Tiva Canyon member have formed in two distinct sets in 3-to 5-m wide swarms spaced at 150 to 200 m. Fractures with tectonic origin postdate the cooling fractures and have no well defined set. Fracture frequencies for short traverses in the densely welded middle part of the Tiva Canyon member are about 6 to 8 fractures per m and reduce to 2 to 4 fractures 3

3 per m for long traverses, which generally include less fractured parts of the exposures. The saturated liquid conductivity of fractured, densely welded tuffin the saturated zone is 3 to 8 orders of tragnitude larger than that of the rock matrix (U.S. DOE,1988).

m It has been assumed that, on a regional scale, two of the in siru principal stresses are e

approximately horizontal and the third is approximately vertical (U.S. DOE,1988). At a depth of 300 m, the venical stress varies from 5 to 10 MPa with a mean of 7 MPa. The ratio of minimum horizontal to vertical stress varies from 0.3 to 0.8 with a mean of 0.55. The maximum horizontal to vertical stress ratio has a mean value of 0.65 and varies from 0.3 to 1.0.

Data on temperature profile with depth show large variations in natural geothermal gradients near Yucca Mountain (U.S. DOE,1988). Upper and lower bounds are 37 *C/km and 20 *C/km, respectively. In situ temperature at the repository horizon (about 300 m beneath the mountain) is about 22 to 30 *C.

2.2 DESCRIPTION

OF THE COUPLING PHENOMENA In a previous study, coupling among the TMHC processes in a rock mass has been described by Manteufel et al. (1993). This study concentrates on the coupled interactions in the near-field environment among the TMH processes relevant to the proposed repository at Yucca Mountain. The effect of chemical reactions on these processes and vice versa were not included in this report, for reasons 2-1

=

I.

mentioned in the introduction section. In this report, only the important features of different aspects of coupling among these processes are described briefly. All the discussions given below are related to the selection of a code which will be used to simulate the rock mass responses caused by these coupled I

processes at Yucca Mountain due to the emplacement of nuclear waste. One-way coupling has been defined as the effect of one process on another. For example, one-way coupling between mechanical and hydrological processes means the effect of mechanical stresses in the rock mass only on the fluid flow.

,I it is represented by a one-way arrow from M to H: M-H. Similarly, other one-way coupled processes are: M-T, H-M, H-T, T-M, T-H. Two-way coupling between two processes means the effect of one on another and vice versa. This is represented by a two-way arrow: M*T, M~H, and H-T. Three-way I

TMH coupling represents interaction among the three processes.

2.2.1 Thermal to Mechanical Coupling I

Heat generated by emplacement of HLW in underground excavations will expand the surrounding rock mass. The restriction of the expansion by the surrounding rock will result in thermally induced mechanical stresses. This thermally induced stress field, in addition to the in situ stresses and the stresses induced due to excavation and repeated seismic effects, can induce normal and shear displacements of the rock joints. This increases the potential for rock mass failure resulting from excessive joint shear displacement. It will also induce microcracks in the rock which may lead to the formation or extension of a fracture network and, thereby, the formation of preferential flow paths.

I The thermal load is a result of heat generated due to decay of the nuclear waste. The magnitude of the thermal load depends on the age and/or the form of HLW, and number and configurations of spent fuel assembly for disposal. The effect of thermal load decreases rapidly with distance from the canister.

I High temperature gradients near the canisters are expected to cause high stress gradients. The thermal stresses are expected to persist during a significant part of the life of the repository (Manteufel et al.,

1993).

Results of preliminary investigations of the thermal propenies of the Topopah Spring member of the Paintbrush Tuff at Yucca Mountain show that the rock has low coefficients of thermal expansion and higher than expected thermal conductivity (Nimick,1990). As a result, emplacement of the waste I

will induce relatively low thermal stresses. It is likely that thermally induced stresses will be of the same order of magnitude as the stresses induced by the excavation at the proposed repository level assuming an initial power of 3.2 KW per canister (Christianson,1988) and seismic events. Further, the effect of I

thermal loads on the mechanical propenies of the rock in the repository horizon is not well understood and needs funher investigation.

Young's modules of a rock sample from the potential repository horizon decreased about 16 percent as temperature was increased from 22 to 150 *C under both uniaxial and biaxial compression (5 MPa) (Price et al.,1987). It is well known that the brittle-ductile transition pressure decreases with increasing temperature (Jaeger and Cook,1979). However, the temperature anticipated for HLW emplacement at Yucca Mountain (U.S. DOE,1988) is not large enough to cause a substantial decrease of the transition pressure.

Arulmoli and St. John (1987), Christianson and Brady (1989), and Bauer and Costin (1990) f have concluded that spalling and slip along discontinuities due to thermal stresses are expected to be minor for a temperature increase in the range of 200 to 240 *C. Results from the studies conducted by 2-2

E

.3 Kemeny and Cook (1990) contradict these conclusions. They concluded that, in a worst case scenario, about 35 percent of the boreholes in the proposed repoeiay a Yucca Mountain will experience some g

kind of rock failure. It should be indicated that this conclusion applies to a specific design, and also that 3

failure of borehole does not mean failure of waste package. Their study includes fracture formation by j

microcrack coalescence and degradation of material stiffness due to time-dependent thermal and mechanical loads.

In summary, the effects of thermal loads on mechanical processes may be significant during operation and containment periods of the proposed IILW repository at Yucca Mountain. It is also potentially important during the isolation period.

2.2.2 Mechanical to Therrnal Coupling

.I Mechanical processes can change the aperture of natural discontinuities of the rock mass which, in turn, can change the effective thermal conductivity. Mechanical energy also can transform into thermal l

enrgy through frictional dissipation, but it is expected to be a minor fraction of the dissipative heat load.

Therefore, mechanical to thermal coupling has not been considered important and is rarely mentioned in connection with the nuclear waste disposal at Yucca Mountain.

.I 2.2.3 Mechanical to Hydrological Coupling Mechanical processes can affect the flow of fluids in the rock mass by changing the joint aperture and the bulk porosity of the rock matrix. Changes in aperture, in turn, would change the permeability of the joints. The change in joint aperture may be due to both normal and shear displacements of the joints. Shear displacement causes dilation which increases the joint aperture. The j

apenure of a joint also increases with the decrease of normal stress acting on it.

h Creation of an opening in the rock mass redistributes the in siru stress field. Stress concentration

'W around the excavation changes the apertures of the existing joints. Seismic loading from earthquakes can also change the aperture of the joints. Change in apertures not only changes the hydraulic conductivities of the rock but may also change the preferential flow paths.

l Rock mass may also fail due to displacements along the joints. Depending on the strength of lg the rock, additional fractures can form which also change the hydraulic conductivities of the rock mass

'g surrounding the excavations. 'Ihe method of excavation can create additional fractures in the surrounding rock. These fractures reduce the load-bearing capacity and increase the bulk hydraulic conductivities of the rock mass. At the proposed repository horizon, which is approximately 300 m below the surface, the failed region around the excavations is expected to be small as the in situ stress field is relatively low and the Topopah Spring welded tuff is quite strong.

2.2.4 Hydrological to Mechanical Coupling l

The state of stress in a rock mass is coupled to the flow of groundwater. The proposed repository is located in an unsaturated zone. In an unsaturated condition, the fluid pressure is expected to be very small (close to atmospheric). From this sense, the H-+M coupling will be insignificant. Also, even if there is a perched water zone in the near vicinity of the repository, the water pressures will still 1

be close to atmospheric as a result of the perched water zone being unconfined.

2-3 I

I

Presence of fluid in the rock can change the mechanical proper Os. In unconfined compression tests (Olsson and Jones,1980), saturated samples of Grouse Canyon tun, a volcanic rock located at the Nevada Test Site, are 24 percent weaker than dry samples. Presence of water increases the coefficient of friction for some minerals and decreases it for other minerals (Jaeger and Cook,1979).

2.2.5 Thermal to Hydrological Coupling

{

Emplacement of HLW in an underground repository increases the temperatures and the thermal gradients near the waste canisters. This can create thermally-driven natural convection currents resulting in accelerated spread of radionuclides; however, this effect may be weaker in the unsaturated zone than-in the saturated zone. In unsaturated rock, the increase in temperature can dry the surrounding rock mass by vaporizing the water. The vapor can flow to the cooler regions of the rock mass and condense. If the rock matrix has a low permeability, condensate (water) may drain considerable distances along rock discontinuities before being totally imbibed by the matrix (hfanteufel et al.,1993). In field experiments (Zimmerman and Blanford,1986), downward liquid flow (presumably condensate flow) near heated zones has been observed. hianteufel et al. (1993) conclude that thermal to hydrological coupling is important at the proposed HLW repository at Yucca hiountain.

2.2.6 Hydrological to Thermal Coupling Hydrological processes may influence thermal processes through condensate dripping in fractures, buoyancy-driven natural convection, changing thermal properties of the medium, and heat pipe effects (hianteufel et al.,1993).

The condensate formed due to condensate drainage, discussed in Section 2.2.5, can flow back into a " boiling zone" influencing the thermal processes. A strong hydrological feedback to the temperature field may occur when condensate drains into a

• boiling zone," thereby maintaining the temperature near 100 *C. The heat will be transferred through conduction into the boiling region where it is removed by vaporization and vapor flow. Condensate dripping in fractures has been observed in field experiments (Zimmerman and Blanford,1986; Patrick et al.,1986; Buscheck et al.,1991).

Fluid flow can affect the temperature distribution in the rock mass by removing heat by convection. hicisture content of a porous medium can influence the thermal conductivity of rock.

wsmussen et al. (1990) report that the thermal conductivity for Apache Leap tuff decreases ap1,roximately 30 percent from fully saturated to oven-dried condition. The relationship between the j

thermal conductivity and the degree of saturation is nonlinear and becomes more important for drier samples.

Although hydrological processes can influence thermal processes in several ways, most do not appear to be significant at Yucca hiountain. The most prominent and important mechanism is condensate dripping in fractures.

2-4

I 3 IDENTIFICATION OF COMPUTER CODES I

Many public domain and commercial computer codes are available for numerical analysis of underground structures. These codes have been designed and developed to solve specific types of problems related to underground excavations. Most of these codes are general purpose codes and have capabilities that may I

not be directly relevant in simulating coupled phenomena in an HLW repository at Yucca hiountain. In order to reduce the number of computer codes for consideration in the selection process, a subjective screening criterion was developed and applied, based on Objective 1 given in Section 1.2 of this report.

I This subjective screening criterion is that a code qualified for consideration should be mechanically based (for stress analysis of underground structures) and coupled with at least the thermal processes. Based on this screening criterion, the following computer codes were identified as suitable for this phase of the study from the list of publications given in Section 1.3. The reported information and the authors'

'I experience with various computer codes were used in identifying these computer codes as candidates for further evaluation.

3DEC (version 1.2) (ITASCA Consulting Group, Inc.,1992c; 1992d)

ABAQUS (version 5.2-1) (Hibbit, Karlsson & Sorensen, Inc.,1992a; 1992b; 1992c; 1992d; 1992e; and 1992f)

ADINA (version 6.1) (ADINA R&D, Inc., 1987; 1992a; 1992b; 1992c; 1992d; 1992e; and 1992f)

ANSYS (version 5) (Swanson Analysis Systems, Inc.,1992a; 1992b; 1992c; and 1992d)

BEASY (version 4) (Computational Mechanics Publications, Inc.,1990)

BEST3D (version 3) (BEST Corp.,1989)

FEHMS (Kelkar and Zyvoloski,1991).

FLAC (version 3.0) (ITASCA Consulting Group, Inc.,1992e; 1992f)

GENASYS (Wijesinghe,1989).

MSC/NASTRAN (version 67) (MacNeal-Schwendler Corp., Inc.,1991a; 1991b; 1991c; 1991d; and 1991e)

ROCMAS (Noorishad et al., 1984; 1992; and Noorishad and Tsang,1989)

SANGRE (Anderson,1986)

STEALTH (version 4-1 A) (Hofmann,1981a; 1981b)

I

• THAMES (Ohnishi et al.,1990).

UDEC (version 1.8) GTASCA Consulting Group, Inc.,1992a; 1992b) 3-1 I

i i

l 4 SELECTION CRITERIA The computer codes identified in Section 3 vary considerably in their capabilities and in their states of development and verification. It is imponant to focus the code evaluation process so that the objectives set forth in Section I can be fulfilled. Toward this end, a list of code evaluation criteria, based on a philosophy that a code should have sufficient ability of modeling the imponant coupled phe:umena as identified in Section 2 and may require relatively less effort for further development, was established.

These selection criteria are listed in the first column of Table 4-1. It is recognized that the approach for

.E developir,g evaluation criteria is relatively subjective; nevenheless, it serves the purpose of this study.

!5 Based on these criteria, information regarding the features and capabilities of each code is given in Table 4-1.

A discussion of additional capabilities of each code, that may be ofinterest but not provided in Table 4-1, is provided in the appendix. Whenever possible, a discussion on the user friendliness of the program and readability of the manual (s) is given. It is recognized that this information is highly subjective, and may

!I.

depend on limited exposure to the code and a less than thorough familiarity with the manual; nevenheless, this information has been included here to provi@ a qualitativejudgement about the difficulties that may be experienced in using the code for modeling me proposed HLW repository at Yucca Mountain, Nevada.

4.1 COUPLED PROCESSES MODELED This is a critical criterion for selecting the code for modeling the interactions among the near-field coupled TMH processes for repository design and performance assessment. This criterion has been divided into several suberiteria following individual coupled processes: T-M, T-H, H -T, M*H),

j j j M*H, T-H, and H -T. A distinction has been made between fluid flow through the matrix, H,

m m

m m

proposed HLW repository (Manteufel et al.,j. As the M-T coupling is not considered and fluid flow through the rock joints, H 1993), it is not included here as a suberiterion. T-lj, H-T, j

T-H, and H *T are used as subtriteria. It was considered that two-way couplings between mechanical m

m and hydrological processes for both matrix and joint flow are important. Therefore, they are designated as two-way coupled processes.

I 4.2 JOINT CONSTITUTIVE LAWS This criterion checks the ability of the computer code to model joints in a jointed rock mass and the associated mechanical behavior. The model should be able to simulate the joint behavior with both shear and normal deformations. There are a number of rock joint models in the literature (Patton,1966;

'g Ladanyl and Archambault,1970; Jaeger,1971; Jaeger and Cook,1979; Goodman,1976; Bandis et al.,

W 1983; Cundall and Hart,1984; Barton et al.,1985; Bandis et al.,1985). Most of the codes include some form of mechanical Mohr-Coulomb or Coulomb rock joint model. Constitutive laws for the thermal and hydrological behavior along the joint are also imponant for the selected code.

4.3 DYNAMIC AND SEISMIC CAPABILITIES The selected code should be able to model the response from a given dynamic or time-dependent l

load. The ability of the code to incorporate the actual time history from an eanhquake loading is included here as a selection criterion.

I 1

1 I.

.g W

4.4 TWO-PHASE FLOW I

Two-phase flow has special sign;ficance in the near-field environment. The heat from the nuclear waste is expected to increase the temperature of surrounding rock mass, which will, in turn, l

increase the vapor pressure of in sim water. Due to increased temperature, part of the liquid phase (water) will vaporize and the pressure in the gas phase (vapor) will rise to achieve thermodynamic equilibrium.

I Significant drying of the rock will occur in the near-field. Water vapor will flow away from the waste package towards the cooler regions where it will condense. This condensate (water) may travel considerable distances along fractures (Manteufel et al.,1993). To simulate this process, the code should I

be able to model two-phase flow with the associated phase changes.

j 4.5 SOURCE CODE If the selected code does not have the capability to model the coupled TMH phenomena that are determined to be important, then further improvements of the code will be necessary. These I

improvements may include incorporation of the results of research activities, including those being carried out by various program elements at the CNWRA. In addition to integrating the results of the relevant research activities into the selected code, this code or abstracted version of it will also serve as a module of the Total System Performance Assessment Computer (TPA) code (Sagar and Janetzke,1991).

Depending on the extent and nature of required improvements, it may be necessary to have the source j

code available at the CNWRA. Alternatively, it may be possible to get the code enhancement work done by the supporting group for a commercial code, with verification and validation work being performed by the CNWRA.

j I

I I

I I

I E

II II I

Table 4-1. List of candidate coupled codes and their function capabilities Criterion UDEC (version 1.8) 3DEC (version 1.2)

Coupled Processes Modeled:

T-M Yes Yes T*H -

No No H -Y No No 3

M-H -

Yes No M*I/m No No l

l T-H No No m

H -T No No m

Joint Constitutive Laws Joint: Continuous-Yielding, Joint: Continuous-Yielding, j

Mohr-Coulomb, Mohr-Coulomb Barton-Bandis Dynamic and Seismic Capabilities Yes Yes Two-Phase Flow No No Source Code Available Not Available L

l ADINA MSC/NASTRAN L

Criterion ANSYS (version 5)

(version 6.1)

(version 67)

Coupled Processes Modeled:

F' T-M Yes Yes Yes T-H -

Yes Yes No H;-Y No Yes No M*H.

No No No L

M*IE No Yes No m

T-H No Yes No m

H -T No Yes No 7

m Joint Constitutive Laws Slide Line:

Slide Line:

Slide Line 2D and 3D Coulomb or l

interfaces / gaps frictionless Dynamic and Seismic Capabilities Yes Yes Yes i

l Two-Phase Flow No Partial (without No i

phase change) j

]

Source Code Not Available Not Available Not Available 4-3

I Table 4-1. List of candidate coupled codes and their function capabilities (cont'd)

Criterion ROCMAS GENASYS THAMES Coupled Processes Modeled:

I T-M Yes Yes Yes T*H.

Yes No No H -Y Yes No No 3

I M-H)

Yes No No M-H Yes Yes Yes m

T-H Yes Yes Yes m

H -T Yes Yes Yes m

Joint Constitutive Laws Goodman Joint model No No E

Dynamic and Seismic Capabilities No No No Two-Phase Flow No No No I

Source Code Not Available Not Available Not Available l

I I

STEALTH i

Criterion FLAC (version 3.0)

SANGRE (version 4-1 A)

Coupled Processes Modeled:

I T-M Yes Yes Yes i

T-H -

No Yes No H -Y No Yes No N-H No Yes No j

M*H Yes Yes Yes m

T-*H No Yes Yes m

H -T No Yes Yes m

Joint Constitutive Laws Interface / Slide Line:

Slide Line:

Slide Line:

Friction, no slip.

with or without Frictionless, tied I

Only a few simple frictior.

interfaces should be Dynamic and Seismic Capabilities Yes No Yes Two-Phase Flow No No No Source Code Not available Available Available g

I

I-Table 4-1. List of candidate coupled codc:; and their function capabilities (cont'd)

Criterion FEHMS BEST3D (version 3)

Coupled Processes Modeled:

I T-M Yes Yes T-H -

No No H -Y No No j

i M*H-Yes No M*H Yes No m

T-H Yes No m

H -T Yes No m

Joint Constitutive Laws Joint Slide Line: spring, friction Dynamic and Seismic Capabilities Yes Yes Two-Phase Flow No No N

Source Code Not Available Available d

I 1

Criterion BEASY (version 4)

ABAQUS (version 5.2-1)

Coupled Processes Modeled:

T-M Yes Yes T-H -

No No H;-Y No No M*H No Yes j

M*H No Yes m

T4H No No i

m H -T No No m

1 g

Joint Constitutive Laws Slide Line:

Interface, Slide Line:

B without friction friction Dynamic and Seismic Capabilities No Yes Two-Phase Flow No No l

Source Code Not Available Not Available B

l l

4-5

lI.

l 5 CODE RANKING g

5.1 APPROACH OF CODE RANKING The primary objective for this code selection is to have a code that can model near-fiwld rock

,I mass behavior to predict coupled TMH responses resulting from the disposal of HLW. The ability of the recommended codes to simulate rock joint behavior under both pseudostatic and dynamic loads and coupled TMH processes is considered to be of paramount importance. Therefore, the first three criteria ig of Table 4-1, that is, coupled processes modeled, joint constitutive laws, and dynamic and seismic lW capabilities, were used to calculate the score and to rank the codes. Codes with high scores will be examined in Phase Il to further determine their capabilities for final code selection. The result will be documented in a separate report. Two-phase flow was not considered as an evaluation criterion in this rating exercise as none of the candidate codes has the desired capability.

5.2 RANKING OF COMPUTER CODES A ranking system was developed based on the three criteria. Among the coupled processes, M-T is not considered significant in the Yucca Mountain environment (Manteufel et al.,1993).

. I Therefore, it was excluded from the ranking system. Also excluded from the ranking system are the coupled processes T-H, and H *T since the thermohydrology research under the Performance m

Assessment element at the CNWRA is developing and validating model(s) of them. The model(s) will be adopted and coded into the selected computer code.

The rock joint model and the two-way coupled process M*H criteria are each given a score j

of 2 while each of the remaining coupled processes and the dynamic capability criterion receives a score of 1. As discussed before, rock joint behavior is considered to be of primary importance in modeling the near-field response of excavations in the proposed repository. Therefore, this feature receives a weighted score of 2. In order to receive a score of 2, the rock joint model in a code should have the capability of modeling rock joint behavior with both shear and normal deformations of thejoints. If only a portion of the joint behavior is simulated, a score of 0 is assigned. M*H describes the coupled behavior of two j

'I individual processes in the near-field. As a result, it also receives twice the score. In the near-field environment around the waste package, the finid flow will be predominantly through the joints or fractures. As a result, the effects of mechanical processes on the fluid flow through the rock matrix and vice versa, M-H, have been considered somewhat less significant than M-H for this code selection m

j process. M*H receives the same score as one-way processes T-+M, T-+H), and Hj*T. Only the ADINA m

system of codes has the capability to model two-phase flow but they cannot simulate phase changes.

Therefore, this criterion is not shown in Table 5-1.

.I It is recognized that the scoring process described above is highly subjective. However, it is considered to be appropriate in the context of TMH modeling. Table 5-1 presents the score of the codes

.I listed in Section 3.

Based upon the scoring system, the candidate codes are ranked in Table 5-2 with I being the I

highest rank. More than one code received the same score and, therefore, the same rank. This subjective i

ranking system uses only seven features to differentiate among the codes. It is not sensitive enough to (I

break the tie among the codes receiving the same rank.

I 5-1 I

l l

I Table 5-1. Ranking analysis of the computer codes Criterion 3DEC ABAQUS ADINA ANSYS BEASY M *H-0 2

0 0

0 M-H 0

1 1

1 0

m T-+M 1

1 1

1 1

T-H; O

O O

O 0

H;-T 0

0 0

0 0

I Rock Joint 2

2 0

0 0

Model Dynamic 1

1 1

1 0

Composite Score 4

7 3

3 1

I Criterion BEST3D FElIMS FLAC GENASYS MSC/NASTRAN M*H; O

2 0

0 0

M*H 0

1 1

1 0

m T-+M 1

1 1

1 1

T-H, 0

1 0

1 0

I i

H;-T 0

1 0

1 0

Rock Joint 0

0 2*

0 2

I Model Dynamic 1

1 1

0 1

Composite Score 2

7 5

4 4

• Only a few interfaces should be used in a model I

I I

s2 I

Table 5-1. Ranking analysis of the computer codes (cont'd)

Criterion ROCMAS SANGRE STEALTH TIIAMES UDEC M -H-2 0

0 0

2 M*H 1

1 1

1 0

m T-M 1

1 1

1 1

T-H; 1

1 0

0 0

H;-T 1

1 0

0 0

Rock Joint Model 2

0 0

0 2

Dynamic 0

0 1

1 1

Composite Score 8

4 3

3 6

Table 5-2. Ranks of the computer codes Code Score Rank ROCMAS 8

1 ABAQUS 7

2 FEHMS 7

2 UDEC 6

3 SANGRE 4

4 3DEC 4

4 MSC/NASTRAN 4

4 ADINA 3

5 ANSYS 3

5 STEALTH 3

5 FLAC 3

5 THAMES 3

5 GENASYS 2

6 BEST3D 2

6 BEASY i

7 I

5-3 I

l.

I 6 DISCUSSION AND RECOMMENDATIONS I

Based on the criteria used in Table 5-1, the total available score for a code is 9. An arbi rary decision was made that if a code receives a total score of more than half of the total available score it should be qualified for the next level of code selection. This decision is imponant because the ranking of the codes is primarily based on the information provided by the corresponding users' manuals. It may be possible that some information in the manuals has been misinterpreted or misunderstood by the authors of this study due to lack of actual verification. Funher, the criteria in Table 5-1 are relatively general in nature.

For example, the comparison of a similar capability between two codes is not factored into the scoring process. Among the fifteen computer codes, four - ROCMAS, ABAQUS, FEHMS, and UDEC -

received a score greater than 5. Of these, only UDEC is a discrete element code; all the others are finite element codes. ROCMAS has been developed at the Lawrence Berkeley Laboratory, California.

I ABAQUS is a commercial code, with an executable version available on lease from Hibbit, Karlsson &

Sorensen,Inc., Rhode Island. FEHMS has been developed at the Los Alamos National Laboratory, New Mexico. UDEC is also a commercial code; both source and executable codes are available from Itasca Consulting Group, Inc., Minnesota. These codes need to be examined further to suppon selection and further development of a panicular code for compliance determination.

Of these codes, only UDEC has been subjected to a qualification study (Brady et al.,1990). Other codes should go through the same qualification process to be compared on an equal basis. There are some problems in subjecting ROCMAS and FEHMS to this qualification study. Although the capabilities of these codes to simulate coupled phenomena among TMH processes have been reported in the open literature (Noorishad et al., 1984, 1992; Noorishad and Tsang,1989; Kelkar and Zyvoloski,1991), the complete codes along with the manuals are not available in the public domain. Therefore, both ROCMAS and FEHMS cannot be examined nrther at this time.

1 l

ABAQUS will undergo a qualification study similar to that of UDEC (Brady et al.,1990). Details of this I

study will be described in a later repon. That report will also include appropriate comparisons between UDEC and ABAQUS and the selection of one of these two codes for further development.

I I

I I

6-1 I

7 REFERENCES I

ADINA R&D, Inc.1987. Automatic Dynamic Incremental Nonlinear Analysis: Theory and Modeling Guide. Watenown, MA: ADINA R&D, Inc.

I ADINA R&D, Inc.1992a. ADINA-INfor ADINA-F User's Manual. Watertown, MA: ADINA R&D, Inc.

ADINA R&D, Inc.1992b. ADINA-INfor ADINA User's Manual. Watertown, MA: ADINA R&D, Inc.

ADINA R&D, Inc.1992c. ADINA Veripcation Manual - Nonlinear Problems. Watertown, MA: ADINA R&D, Inc.

l l

ADINA R&D, Inc.1992d. ADINA Verification Manual-Linear Problems. Watertown, MA.

I ADINA R&D, Inc.

ADINA R&D, Inc.1992e. ADINA-T Ver#ication Manual. Watertown, MA: ADINA R&D, Inc.

ADINA R&D, Inc.1992f. ADINA-F Ver#ication Manual. Watertown, MA: ADINA R&D, Inc.

Anderson, C.A.1986. SANGRE: A Finite Element Codefor Fluid Migration, Heat Transport, and Faulting in Highly Deformable, Porous GeologicalMedia. LA-10666-MS. Los Alamos, NM: Los Alamos National Laboratory.

Arulmoli, K., and C.M. St. John.1987. Analysis ofHorizontal Waste Emplacement Boreholes of a Nuclear Waste Repository in Tuf. SAND 86-7133. Albuquerque, NM: Sandia National Laboratories.

Bandis, S.C., A.C. Lumsden, and N.R. Barton.1983. Fundamentals of rock joint deformation.

InternationalJournal ofRock Mechanics and Mineral Sciences & Geomechanical Abstracts 20(6):249-268.

Bandis, S.C., N.R. Barton, and M. Christianson.1985. Application of a new numerical model I

of joint behavior to rock mechanics problems. Fundamentals of rock joints. Proceedings of the International Symposium on Fundamentals of Rock Joints. Bjorkliden, Lulea, Sweden:

Centek Publishers: 345-356.

Barton, N.R., S.C. Bandis, and K. Bakhtar,1985. Strength, deformation and conductivity coupling of rock joints. International Journal of Rock Mechanics and Mineral Sciences &

I Geomechanical Abstracts 22(3):121-140.

Bauer, S.J., and L.S. Costin.1990. Thermal and Mechanical Codes First Benchmark Exercise Part1I Elastic Analysis. SAND 89-0757. Albuquerque, NM: Sandia National Laboratories.

BEST, Corp.1989. BEST3D User's Manual. Buffalo, NY: State University of New York at Buffalo and Best, Corp. United Technologies.

7-1 I

I Biot, M.A.1955. Theory of elasticity and consolidation for a porous anisotropic solid. Journal ofApplied Physics 26(2):182-185.

Brady, B.H.G., S.M. Hsiung, and A.H. Chowdhury.1990. Qualyication Studies on the Distinct Element Code UDEC Against Some Benchmark Analytical Problems. CNWRA 90-004, San Antonio, TX: Center for Nuclear Waste Regulatory Analyses.

Buscheck, T.A., R. Carlson, W. Daily, K. Lee, W. Lin, N. Mao, A. Ramirez, T. Ueng, H.

Wang, and D. Watwood.1991. Protonpe Engineered Barrier Systems Field Test (PEBSFT)

Final Report. A.L. Ramirez, ed. UCID-106159. Livermore, CA: Lawrence Livermore National Laboratory.

Buscheck, T.A., and J.J. Nitao.1993. The analysis of repository-heat-driven hydrothermal flow at Yucca Mountain. Proceedings of the 4th High-Level Radioactive Waste Management Conference. New York, NY: American Society of Civil Engineers.

I Christianson, M.C. I988. Sensitivity of the Stability of a Waste Emplacement Drift to Variation in Assumed Rock Joint Parameters in Welded Tuf. NUREGICR-5336. Minneapolis, MN:

ITASCA Consulting Group, Inc.

Christianson, M.C., and B.H.G. Brady.1989. Analysis of Alternative Waste isolation Concepts.

NUREG/CR-5389. Washington, DC: U.S Nuclear Regulatory Commission.

Computational Mechanics Publications, Inc.1990. BEASY 7he Boundary Element Analysis System User Guide. Boston, MA: Computational Mechanics Publications, Inc.

Cundall, P.A. and R.D. Hart. 1984. Analysis of Block Test No.1 Inelastic Rock Mass Behavior: Phase 2 - A Characterization of Joint Behavior (Final Report). Minneapolis, I

MN: ITASCA Consulting Group Inc.

Goodman, R.E.1976. Methods of Geological Engineering in Discontinuous Rock. St. Paul, E

MN: West Publishing Company.

Hibbit, Karlsson & Sorensen, Inc.1992a. ABAQUS Theory Manual Version 5.2. Pawtucket, I

RI: Hibbit, Karlsson & Sorensen, Inc.

Hibbit, Karlsson & Sorensen, Inc.1992b. ABAQUS Veryication Manual Version 5.2.

Pawtucket, RI: Hibbit, Karlsson & Sorensen, Inc.

Hibbit, Karlsson & Sorensen, Inc.1992c. ABAQUSiStandard User's I Manual Version 5.2.

Pawtucket, RI: Hibbit, Karlsson & Sorensen, Inc.

Hibbit, Natlsson & Sorensen, Inc.1992d. ABAQUS/ Standard User's 11 Manual Version 5.2.

Pawtucket, RI: Hibbit, Karlsson & Sorensen, Inc.

Hibbit, Karlsson & Sorensen, Inc.1992e. ABAQUSiExplicit User's Manual Version 5.2.

Pawtucket, RI: Hibbit, Karlsson & Sorensen, Inc.

7-2 I

i Hibbit, Karlsson & Sorensen, Inc.1992f. ABA QUS/ Explicit Example Problems Manual Version 5.2. Pawtucket, RI: Hibbit, Karlsson & Sorensen, Inc.

Hofmann, R.1981a. STEALTH - A Lagrange Explicit Finite Diference Code for Solids, Structural and 7hermohydraulic Analysis: Volume 1A: User's Manual - Theoretical Background and Numerical Equations. San Leandro, CA: Science Application, Inc.

Hofmann, R.1981b. STEAL 7H - A Lagrang# Explicit Finite Diference Code for Solids, Structural and Tnermohydraulic Analysis, Introduaion and Guide. EPRI-NP-2080-CCM-SY.

E Palo Alto, CA: Electric Power Research Institute.

International Society for Rock Mechanics. 1988. ISRM: Rock Engineering Software.

InternationalJournal ofRock Mechanics and Mineral Sciences & Geomechanical Abstracts 25(4).

ITASCA Consulting Group, Inc.1992a. UDEC UniversalDistinct Element Code Version 1.8 Volume I: User's Manual. Minneapolis, MN: ITASCA Consulting Group, Inc.

l ITASCA Consulting Group, Inc.1992b. UDEC Universal Distinct Element Code Version 1.8 Volume II: Venpcation and Example Problems Manual. Minneapolis, MN: ITASCA Consulting Group, Inc.

ITASCA Consulting Group, Inc.1992c. 3DEC Three Dimensional Distinct Element Code Version 1.2 Volume It User's Manual. Minneapolis, MN: ITASCA Consulting Group, Inc.

ITASCA Consulting Group, Inc.1992d 3DEC 1hree Dimensional Distinct Element Code Version 1.2 Volume II: Venpcation and Example Problems Manual. Minneapolis, MN:

ITASCA Consulting Group, Inc.

ITASCA Consulting Group, Inc.1992e. FIAC Fast IAgrangian Analysis Code Version 1.3 Volume It User's Manual. Minneapolis, MN: ITASCA Consulting Group, Inc.

ITASCA Consulting Group, Inc.1992f. FL4C Fast IAgrangian Analysis Code Version 1.3 Volume IIt Venpcation and Example Problems Manual. Minneapolis, MN: ITASCA E

Consulting Group, Inc.

Jaeger, J.C.1971. Friction of rocks and the stability of rock slopes. Geotechnique 21:97-134.

Jaeger, J.C. and N.G.W. Cook.1979. Fundamentals of Rock Mechanics. Third Edition.

London: Chapman and Hall.

!I

~

Kana, D.D., B.H.G. Brady, B.W. Vanzant, and P.K. Nair.1991. Critical Assessment of Seismic and Geomechanics Literature Related to a High-Level Nuclear Waste Underground j

Repository. NUREG/CR-5440. Washington, DC: U.S. Nuclear Regulatory Commission.

!I i

7-3

I Kelkar, S. and G.A. Zyvoloski.1991. An Eficient, Ihree Dimensional, Fully-Coupled Hydro-7hermo-Mechanical Simulator: FEHMS. LA-UR-90-3750. Los Alamos, NM: Los Alamos National Laboratory.

Kemeny, J.M., and N.G.W. Cook.1990. Demonstration of a Risk-Based Approach to High-Level Waste Repository Evaluation in Rock Mechanics and Crustal Stresses, R.K.

E McGuire, ed. EPRI-NP-7507. Palo Alto, CA: Electric Power Research Institute.

Ladanyi, B. and G. Archambault.1970. Simulation of the shear behavior of a jointed rock I

mass. Proceedings of the lith Symposium on Rock Mechanics. New York, NY: AIME:

105-125.

MacNeal-Schwendler Corp., Inc.1991a. NAS7RAN Theoretical Manual. Los Angeles, CA:

MacNeal-Schwendler Corp., Inc.

MacNeal-Schwendler Corp., Inc. I991b. NASTRAN Application Manual Volume I, Version 67.

Los Angeles, CA: MacNeal-Schwendler Corp., Inc.

MacNeal-Schwend1er Corp., Inc.1991c. NASTRAN Application Manual Volume 11, Version 67.

Los Angeles, CA: MacNeal-Schwendler Corp., Inc.

l MacNeal-Schwendler Corp., Inc.1991d. NASTRAN User's Manual Volume I, Version 67. Los Angeles, CA: MacNeal-Schwendler Corp., Inc.

l MacNeal-SchwendIet Corp., Inc.1991e. NASTRAN User's Manual Volume II, Version 67. Los Angeles, CA: MacNeal-Schwendler Corp., Inc.

Manteufel, R.D., M.P. Ahola, D.R. Turner, and A.H. Chowdhury.1993. A Literature Review I

of Coupled Thermal-Hydrologic-Mechanical-Chemical Processes Pertinent to the Proposed High-Level Nuclear Waste Repository at Yucca Mountain. NUREGICR-6021. Washington, DC: U.S. Nuclear Regulatory Commission.

Nataraja, M.S., and T. Brandshaug.1992. Staf Technical Position on Geologic Repository Operations Area Underground Facility Design - 7hermal Loads. NUREG-1466.

i Washington, DC: U.S. Nuclear Regulatory Commission.

Nimick, F.B.1990. The 7hermal Conductivity of Seven Thermal / Mechanical Units at Yucca E

Mountain, Nevada. SAND 88-1387. Albuquerque, NM: Sandia National Laboratories.

Noorishad,3. and C.F. Tsang.1989. Recent Enhancements of the Coupled Hydro-Mechanical Code: ROCMAS II. Technical Report 89:4. Stockholm, Sweden: Swedish Nuclear Power Inspectorate.

Noorishad, J., C.F. Tsang, and P.A. Witherspoon. 1984. Coupled thermal-hydraulic-mechanical phenomena in saturated fractured porous rocks: Numerical approach. Journal of Geophysical Research 89(B12):10365-10373.

g 74 I

Noorishad, J., C.F. Tsang, and P.A. Witherspoon.1992. Theoretical and field studies of coupled hydromechanical behavior of fra:tured rocks - 1. Development and verification of i

a numerical simulator. International Journal of Rock Mechanics and Mineral Science &

Geomechanical Abstracts 29(4):401-409.

Norris, A.E.1989. The use of chlorine isotope measurements to trace water movements at I

Yucca Mountain. Proceedings of the American Nuclear Society Topical Meeting on Nuclear Waste Isolation in the Unscrurated Zcme (Focus 89). Las Vegas, NV: American Nuclear Society.

Ohnishi, Y., M. Nishigaki, A. Kobayaski, and S. Akiyama.1990. Three dimensional coupled thermo-hydraulic mechanical analysis code with PCG method. Proceedings of the International Symposium GEOVAL-90. Sto:kholm, Sweden: Swedish Nuclear Power Inspectorate: 14-17.

Olsson, W. A., and A.K. Jones.1980. Rock Mechanics Properties of Volcanic Tuffsfrom the Nevada Test Site. S ANDS 0-1453. Albuquerque, NM: Sandia National Laboratories.

E Patri k, W.C. et al.1986. Spent Fuel Test - Climat: An Evaluation of the Technical Feasibility of Geologic Storage of Spent Nuclear Fuelin Granite. UCRL-53702. Livermore, CA: Lawrence Livermore National Laboratory.

Patton, F.D.1966. Multiple modes of shear failure in rock. Proceedings of fst Congress of international Societyfor Rock Mechanics. Lisbon: Imernational Society for Rock Mechanics:

1:509-513.

Price, J.G., S.T. Conlon, and C.D. Henry.1987. Tectonic controls on orientation and size of epithermal veins. North American Conference on Tectonic Control of Ore Deposits. Rolla, i-MO: University of Missouri: 36-46.

Puig-Pey, J., and C. A. Brebbia.1987. Computer Aided Engineering Systems Handbook. New E

York, NY: Computational Mechanics Publications and Springer-Verlag.

Rasmussen, T.C., D.D. Evans, P.J. Sheets, and J.H. Blanford.1990. Unsaturated Fractured

~I Rock Characterization Methods and Data Sets at the Apache Leap Tuff Site.

NUREG/CR-5596. Washington, DC: U.S. Nuclear Regulatory Commission.

Sagar B., and R.W. JanetzLe.1991. Total System Performance Assessment Computer Code:

Description of Executive Module. CNWRA 91-009. San Antonio, TX: Center for Nuclear Waste Regulatory Analyses.

.I Swanson Analysis Systems, Inc.19911. ANSYS User's Manualfor Revision 5.0 Volume 1 Procedures. Report D-R300:50-1. Houston, PA: Swanson Analysis Systems, Inc.

I Swanson Analysis Systems, Inc.1992b. ANSIS User's Manualfor Revision 5.0 Volume Il Commands. Report D-R300:50-2. Houston, PA: Swanson Analysis Systems, Inc.

E 7-5 I

I Swanson Analysis Systems, Inc.1992c. ANSIS User's Manualfor Revision 5.0 Volume 111 Elements. Repon D-R300:50-3. Houston, PA: Swanson Analysis Systems, Inc.

Swanson Analysis Systems, Inc.1992d. ANS15 User's Manualfor Rerision 5.0 Volwne IV 7heory. Report D-R300:50-4. Houston, PA: Swanson Analysis Systems, Inc.

U.S. Department of Energy.1988. Site Characteri:ction Plan: Yucca Mountain Site, Nevada Research and Development Area. DOE /RW-0199. Washington, DC: U.S. Department of Energy.

U.S. Nuclear Regulatory Commission.1992. Disposal of High-Level Radioactive Wastes in Geologic Repositories. Title 10, Energy, Part 60 (10 CFR Part 60). Washington, DC: Office of the Federal Register.

Wijesinghe, A.M.1989. Hydrothermomechanical simulator development task. Repository Technology Program Activities: FY 1988. 3. Yow, Jr. et al., eds. UClD-21600. Derkeley, CA: Lawrence Livermore National Laboratory: 3-14.

I Zimmerman, R.W., and M.K. Blanford.1986. Expected thermal and hydrothermal environments for waste emplacement holes based on G-Tunnel heater experiments. 27th U.S.

Symposium on Rock Mechanics. H. Hartman, ed. Littleton, CO: Society of Mining Engineers: 874-882.

I I

I I

I I

I i

i 7-6

.I I

E-APPENDIX Additional Information On The Computer Codes This appendix includes supplemental discussions on each code described in Table 4-1 of this report. Any additional information which is not induded in the table or is not adequately described in the I.

text is given here. As far as possible, the information given in Table 4-1 is not duplicated in this appendix, although some repetition has been necessary to make the discussion complete.

A.1 3DEC (VERSION 1.3) 3DEC (3 dimensional Distinct Element Code) is a three-dimensional distinct element code to simulate the behavior ef jointed rock masses under quasi-static and dynamic loads (ITASCA Consulting i

Group, Inc.,1992.c; 1992d). It is written mainly for 80386 based IBM-compatible microcomputers. 3DEC has been developed especially for geomechanical problems. Its features includ-I

• The rock mass modeled as a three-dimensional assemblage of rigid or deformable blocks

• The discontinuities regarded as distinct boundary interactions between these blocks I

An explicit solution algorithm that accommodates both large displacement and rotation j

Structural elements available to simulate rock reinforcement and internal supports One-way coupling from thermal to mechanical processes, assuming the material is thermally I

homogeneous and isotropic with thermally invariant material properties, and Point or line heat sources embedded in the material Material constitutive models include:

Isotropic linear elastic, and Elastic-plastic with Mohr-Coulomb failure criterion

,I 3DEC has the option of modeling nonreflecting or viscous boundaries in a problem. This is very j

useful in wave propagation problems to simulate the passage of the waves outside the modeled region.

A.2 ABAQUS SYSTEM (Version 5.2)

The ABAQUS system of programs indudes ABAQUS/ Standard, ABAQUS/ Explicit, and ABAQUS/ Post. ABAQUS/ Standard (Hibbit, Karlsson & Sorensen, Inc.,1992e; 1992d) is a general purpose finite element program. ABAQUS/ Explicit (Hibbit, Karlson & Sorensen, Inc.,1992e; 1992f) is an explicit dynamic finite element program fully vectorized for use on supercomputers. ABAQUS/ Standard can analyze both linear and nonlinear static and dynamic problems. Nonlinear stress analysis can have three sources for nonlinearity.

• Material nonlinearity: material models depend on the loading history which requires solution following the actual loading sequence.

Geometric nonlinearity: due to large displacement and rotations.

Boundary nonlinearity: due to contact or interface between two regions. In static analysis, there is no loss of energy when two nodes come in contact, although iteration is necessary to determine whether or not they are in contact. In dynamic applications, some energy is lost in I

generating stress waves when two nodes come in contact. ABAQUS automatically solves A-1 I

I impact equations at contact to provide new initial conditions for the continuation of the dynamic response.

ABAQUS has the rezoning capabilities to develop a new mesh for the problem when large deformation has distorted the original mesh considerably.

1 ABAQUS/ Standard has two-way coupling between:

Mechanical and thermal processes, and

,I Mechanical and hydrological processes 1

If the temperature distribution does not depend on the stress solution, then one-way coupled i

thermal stress analysis can be carried out. Heat transfer analysis includes:

• Temperature distribution, which can be time-dependent Material propenies. which can be temperature-dependent, and Conduction, forced convection caused by fluid flowing through the mesh, heat storage (specific heat and latent heat), and boundary radiation and convection Two-way coupled pore fluid diffusion and stress analysis can be carried out for problems involving partially and/or fully saturated fluid flow and include:

ilE Transient or steady state formulation Mechanical part of the model based on effective stress principle ig Continuity equation for the mass of the wetting fluid in a unit volume of the medium, and

'3 Partially saturated flow occurring when the wetting fluid is absorbed into or exsorbed from the medium by capillary action Features of interface elements available in ABAQUS include:

Interface elements for problems involving two deforming bodies that may undergo large I

relative motions Either temperature or pore pressure prescribed as an additional nodal degree of freedom to l

perform analysis coupled with the mechanical responses, and Slide line elements available to model the interaction between two deformable bodies along the slide line where separation and sliding of finite amplitude and arbitrary rotation of the surfaces may arise An equivalent continuum material model containing a high density of parallel joint :;urfaces in l

different orientations is available. The model provides for opening of the joints and frictional sliding in l

each of these systems. Failure and flow parameters are given as functions of temperature and other l

predefined variables. Under compressive stress, the joints can slide following Coulomb criterion. Bulk failure of the medium is based on the Drucker-Prager failure criterion.

Both ABAQUS/ Standard and ABAQUS/ Explicit can calculate the responses of the model from I

time-dependent applied loads such as an earthquake loading history (displacement, velocity, or IE acceleration). In explicit dynamic analysis, the equations of motion for the body are integrated using the j

jM explicit central difference method. The explicit dynamic analysis procedure is based upon the I

A-2 I

implementation of this explicit integration rule together with the diagonal or " lumped" element mass matrices. Infinite elements are used in boundary value problems defined in infinite or unbounded regions or problems in which the region ofinterest is small compared to the surrounding medium. Isotropic linear I

material properties are assumed in the infinite elements. Infinite elements provide a static force that is present at the start of the dynamic response analysis on this boundary. The infinite elements introduce additional normal and shear tractions at the boundary, with boundary damping constants chosen to minimize I

the reflection of dilational and shear wave energies back into the mesh of the finite element model.

A.3 ADINA, ADINA-T, AND ADINA-F (Version 6.1)

I ADINA, ADINA-T, and ADINA-F are general purpose finite element codes developed by ADINA R&D, Inc. (1987; 1992a; 1992b; 1992c; 1992d; 1992e; and 1992f). The programs are numerically coupled to simulate the coupled interactions among the thermal, mechanical, and hydrological processes.

ADINA is the mechanical code for calculating displacements, stresses, etc., due to loads applied on the model. ADINA-T calculates the distribution of temperature in the modeled region from given heat sources I

and heat fluxes. ADINA-F is designed for fluid flow and heat transfer analysis (conduction, convection, and radiation). ADINA uses the temperature distribution calculated by ADINA-T to determine the distribution of thermal strains and stresses over the modeled region. The ADINA system of codes does not I

solve the coupled equation for simultaneous stress and fluid flow analysis. Instead, ADINA and ADINA-F are coupled internally so that the results from fluid flow analysis (ADINA-F) are fed to ADINA for Stress calculations. Results from ADINA become input to ADINA-F. This process is repeated until the results converge.

ADINA has both static and dynamic analysis capabilities in two and three dimensions and can handle problems associated with material nonlinearity, large displacement and small strain (less than i

I 2 percent), and large displacement and large strain.

Material constitutive models include:

Linear elastic, linear orthotropic, thermo-isotropic and thermo-orthotropic elastic, rock model, cap model, von Mises isothermal plasticity, and curve description I

User-supplied material model in the form of FORTRAN code [ Code recompiling and relinidng are necessary]

Cracking model Isothermal plastic material models including the von Mises, Drucker-Prager, and Ilyushin, and Thermo-elasto-plastic and creep models including the effects of thermal strains, time-independent plastic strains, and time-dependent creep strains E

Other features of the ADINA system of codes include:

I Element birth and death options available for adding or taking out elements from the total system of elements.

• Two-and three-dimensional contact surfaces; both sticking and frictionless or frictional sliding can be modeled. Repeated contact and separation are permitted in any sequence.

ADINA-IN is the preprocessor and ADINA-PLOT is the postprocessor for the ADINA system of codes. Interface with PATRAN for preprocessing and postprocessing is available through another code I

TRANSOR.

A-3

,I

I The manuals available for the ADINA system of codes are far from user-friendly. They lack a clear discussion of the capabilities of each code. The theory manual is more than five years old and does not report the capabilities incorporated in the codes in the last five years. Most of the discussions, I

especially in user's manuals, are buried under individual keywords used in the program. Because of sketchy explanations of the capabilities and less than desirable discussions of the verification problems, it is hard for a new user to get a good grasp of the system in a relatively short time.

A.4 ANSYS (Version 5.0)

ANSYS is a general purpose finite element code developed by Swanson Analysis Systems,Inc.,

Houston, Pennsylvania (1992a; 1992b; 1992c; and 1992d). It is capable of analyzing two-and three-dimensional static and dynamic problems. Both linear and nonlinear problems can be analyzed by ANSYS. Three types of nonlinear problems can be analyzed:

Nonlinear behavior due to changing status - nonlinear elements (contact element), birth and death option Geometric nonlinearities - large deformation, and Material nonlinearities - nonlinear stress-strain relations Material constitutive models include:

I Elastic

• Plastic Viscoplastic Viscoelastic I

Creep Swelling, and

• Temperature-dependent properties ANSYS also can model coupling between thermal and mechanical processes in the structure. Fluid flow through a porous medium is analyzed using a separate program called FLOTRAN, a computational I

fluid dynamics program which has direct interface with ANSYS. Both batch and interactive modes are available for analyzing problems.

Both linear and nonlinear (material and geometric nonlinearity) analyses can be carried out for static problems. Dynamic analyses can be carried out with three different solution schemes. Five different damping schemes are available for analyzing transient dynamic problems. Contact elements are used to model interfaces between two regions. Contact elements fail in tension and follow Coulomb friction law in sliding. The material properties can be temperature dependent. Contact elements operate bilinearly only in a static or a nonlinear transient dynamic analysis. Infinite boundary elements are available to model infm' ite domain in thermal analysis.

ANSYS has the built-in capabilities of generating finite element models interactively.

I Postprocessing and graphical display of outputs are also possible in ANSYS, which interfaces with PATRAN for model generation and postprocessing of the results.

I I

^"

I J

A.5 BEASY t

l BEASY is a general purpose boundary element computer code (Computational Mechanics Publications, Inc.,1990) that allows continuum modeling of two-dimensional, axisynunetric, and three-dimensional problems for:

• Steady state potential flow (i.e., heat transfer and flow through porous media)

}

• Time-dependent (or transient) heat transfer Static, linear elastic stress analysis, and

• Thermal stress analysis (BEASY will perform the full heat transfer analysis, followed by the stress analysis)

A major application of BEASY is in the area of fracture mechanics.

f BEASY offers two different types of elements - continuous and discontinuous. Continuous elements have nodes on the edges and corners of the element. These nodes are shared with neighboring elements. Discontinuous elements do not share nodes with neighboring elements. Therefore, the nodes are I

not at the edges and corners of the element, but are displaced towards the element centroid. Thus, for discontinuous elements, continuity of variables is not enforced between neighboring elements. The use of discontinuous elements is necessary to obtain accurate results when there is an actual discontinuity in the L

problem, such as a crack tip.

Zoning within a BEASY model can be done to simulate, for instance, two or more distinct

[

materials or to create an interface within the interior of a model. If there is some design feature internal to the material, then it may be effective to split the model into zones to create an element in a particular location suitable for applying an interface condition to model the feature. The following three interface l

conditions are considered by BEASY:

No contact - the two surfaces are not in contact with each other.

l Sticking contact - the surfaces are in contact, but the tangential force is not sufficient to overcome the static frictional resistance.

Sliding contact - tne surfaces are in contact, and the tangential force has become large enough to overcome static frictional resistance. The surfaces slide tangentially and the dynamic frictional coefficient is used to calculate the frictional resistance.

Other conditions also can be imposed along interface elements, such as prescribed potential, added flux density, convection condition, thermal contact resistance (membrane), prescribed displacement, or u

added traction.

c' In simulating steady state or transient potential ficw problems with BEASY, the material properties may be nonlinear functions of the potential. For example, in heat transfer analysis, the conductivity, specific heat, density, and thermal diffusivity may be functions of temperature. However, for L

stress analysis, the code is limited to only static analysis with linear elastic material properties. Dynamic loading is considered in a pseudostatic way and is not adequate if the effects of true dynamic loading need to be simulated.

Preprocessing and postprocessing can be carried out with BEASY-1MS which has direct interface with BEASY.

A-5 I

A.6 BEST3D BEST3D (BEST Corp.,1989) is a three-dimensional boundary element method. The code is capable of elastic static and dynamic analysis, thermal analysis, as well as nonlinear stress analysis. For nonlinear analysis, material models, including a nonlinear strain-hardening Von Mises model, a multi-surface cyclic plasticity model, and a thermally sensitive viscoplastic model, are available. For stress I

analysis problems, material properties can be input as a function of temperature (e.g., elastic modulus, thermal expansion coefficient).

Interfaces also can be simulated easily with BEST3D. The following interface conditions are recognized in the code:

Fully bonded interface Sliding interface in which movement in the plane tangent to the interface is allowed - only normal displacement compatibility is imposed across the interface I

Spring resistance along an interfa:e between corresponding surfaces - the tractions across this interface are linearly related to the difference in displacements between the two surfaces, and

• Thermal resistance between the corresponding surfaces - the flux across this interface is linearly related to the temperature difference between the two surfaces Material constitutive models include:

lsetropic elastic Cross-anisotropic Isotropic plasticity with variable hardening I

• Kinematic plasticity with multiple yield surfaces, and

• Anisotropic viscoplastic with temperature-dependent material behavior The BEST3D code cannot simulate fluid flow. The only coupling mechanism is thermal to mechanical.

A.7 FEHMS I

FEHMS (Finite Element IIeat Mass Stress) (Kelkar and Zyvoloski,1991) is a finite element code capable of solving two-dimensional and three-dimensional coupled problems of TMH effects in fractured rocks. Its fully implicit coupled formulation allows for large time steps to be used in the analysis.

Fluid flow through the compressibn fractured rock is modeled using the mass balance equation I

combined with Darcy's law. 'Ihe permeability and porosity of the rock matrix can vary as a known function of the local fluid pressure. The permeability of the discontinuities is expressed as a power law function of the fracture aperture, calculated from the stress-displacement equations.

Assuming that the rock and fluid are in equilibrium everywhere in the modeled region, the energy of the system is balanced. Solid displacements are calculated using the equilibrium equation with Biot's I

poroplastic equations for small displacements (Biot,1955). Linear thermoelastic formulations are used to model the effects of thermal expansion or contraction. Points on the discontinuity surfaces also experience the effects of pore pressure present in the discontinuity. Using the theory of linear porcelasticity for small strain including thermal load, the strain tensor is calculated.

A-6 I

Thermal to hydrological coupling occurs through the pressure and temperature dependence of fluid density, viscosity, and enthalpy. The porosity and permeability of the rock vary with the change in pressure, tempe:ature, and displacement. As a result, tha mass balance and the conservation of energy equations are affected. Mechanical to hydrological coupling occurs due to the strong dependence of fracture permeability on the normal and shear displacements of the surfaces. Inclusion of the pore pressure term in Biot's equation accounts for the hydrological to mechanical coupling. Thermal to mechanical coupling I

is taken into account through the temperature term in the thermoelasticitv equation.

At the time of preparing this report, the relevant manual was not available. The above information I

was extracted from Kelkar and Zyvoloski(1991) l A.8 FLAC (Version 3.0) i I

FLAC (Fast Lagrangian Analysis of Continua) is a two-dimensional explicit finite difference program for modeling soil and rock behavior (Itasca Consulting group, Inc.,1992e; 1992f). FLAC is oriented specifically to geotechnical applications, with the following features:

Cable (grouted) elements and beam elements for modeling concrete, shoterete, reinforcement, and similar features i

Flexible and easily-controlled excavation sequence modeling Interfaces to represent joints or thin seams that are characterized by Coulomb sliding and/or I

tensile separation Automated mesh generation Groundwater flow and consolidation modeling Infinite elastic boundary conditions, and I

Built-in programming language (FLACish or FISH) which enables the user to define new variables and functions Material constitutive models include:

1 Null Isotropic elastic Transversely isotropic elastic Mohr-Coulomb plasticity Ubiquitous joint Strain-hardening / softening i

Double yield, and l

j Creep The standard version of FLAC provides approximately 2,000 elements of Mohr-Coulomb material with 2 Mbytes of RAM and approximately 10,000 elements with 4 Mbytes of RAM.

Specific features of FLAC include:

,g

'B Groundwater / Consolidation Model Upgrades I

Automatic calculation of phreatic surface Approximated unsaturated flow A-7 j

I I

Porosity and permeability dependent on volumetric strain Effective stress used on interfaces Fluid streamlines I

Fluid panicle tracking Infinite Elastic Boundary Condition Analytic solution used for infinite elastic boundary Experimental Double-Yield Constitutive Model to simulate, for example, volumetric yielding of backfill A.9 GENASYS GENASYS (Geotechnical Engineering Analysis System) is a two-and three-dimensional hybrid I

boundary element-finite element code for computing the coupled fracture flow, heat flow, and deformation response of fractured rock mass (Wijesinghe,1989). The discontinuities in the rock mass are modeled as one or more discrete curved surfaces. The fracture surfaces can have any arbitrary shape and orientation, I

and can intersect to form a fracture network. The far-field response is simulated using a linear elastic material model. The fractures show nonlinear behavior. The hybrid formulation keeps the computational problem at a manageable level.

At the time of preparing this report, the relevant manual was not available. The above information was extracted from Wijesinghe (1989).

A.10 MSC/NASTRAN MSC/NASTRAN is a general purpose finite element program, developed originally at the I

National Aeronautic and Space Administraticn and, at present, maintained and marketed by the MacNeal-Schwendler Corporation, Los Angeles, California. The code can model both linear and nonlinear problems (material and geometric nonlinearity) under static and dynamic loads, including piecewise linear analysis of nonlinear static response and transient analysis of linear dynamic problems. It has extensive error checking capabilities.

Equivalent loads due to thermal expansion are calculated by separate routines to analyze the heat transfer in the medium. Temperature at grid points, temperature gradients, and heat flow into the elements are calculated. Temperature-dependent thermal expansion coefficients and elastic moduli can be specified for each material type. Radiation from a distant source can be modeled by a prescribed flux into a surface element which depends upon the orientation of the radiation vector relative to the element.

Three types of dynamic analysis can be carried out in MSC/NASTRAN:

• Eigenvalue extraction I

• Frequency response analysis, and

• Transient response analysis I

A-8 I

4

I For each type of problem, both direct and modal formulations can be specified. Different load histories can be specified for each point in the structure. Loads can be either time-dependent or frequency-dependent.

Gap elements are used to model an interface with frictional properties. Two types of gap elements are available:

Adaptive, and Nonadaptive I

For adaptive elements, static and dynamic coefficients of friction are used. Adaptive gap elements have special features such as gap-induced stiffness update, gap-induced bisection, subineremental process, and adjustment of penalty values or axial stiffness. When the gap is open, a small stiffness is defined in the axial direction and there is no transverse stiffness. If the element is closed, that is, the relative displacement between the surface is less than the initial gap, the axial stiffness has a very large value I

relative to the adjacent structures. The gap has a transverse stiffness until the frictional force is exceeded and slippage starts to occur.

MSC/NASTRAN can simulate fully coupled fluid-structure interactions but the principal applications are in the areas of acoustics and noise-control. Development of the model and analysis of the results can be done using PATRAN.

I A.ll ROCMAS ROCMAS (ROCK Mass Analysis Scheme) is a finite element code developed at the Lawrence I

Berkeley Laboratories for solving two-and three-dimensional problems of coupled thermal, hydraulic, and mechanical processes in the geologic medium (Noorishad et al., 1984; 1992; Noorishad and Tsang,1989).

The discontinuities are represented explicitly as four-noded joint elements with strain-softening behavior I

for stress analysis and one-dimensional line elements for fluid flow in the discontinuities. The peak shear stress of the discontinuities is based on the criterion of Ladanyi and Archambault (1970). Normal stress and dilation behavior are modeled using Goodman's joint element (Goodman,1976). The discontinuities are modeled as parallel plates for fluid flow calculations.

At the time of preparing this repon, the relevant manual was not available. The above information was extracted from Noorishad et al. (1984,1992) and Noorishad and Tsang (1989).

A.12 SANGRE SANGRE is a two-dimensional finite element code for simulating fluid migration, heat transport, and faulting in highly deformable, porous geologic media (Anderson,1986). It can carry out coupled fluid I

flow and structural deformation, including large deformation and faulting, as well as heat transport by conduction and convection. It is claimed that the SANGRE code can be altered easily to include formation of fractures in the medium and subsequent changes of permeability.

SANGRE solves the consolidation equations of Biot (Biot,1955). These equations have been modified to account for inelastic creep of solid matrix. The pore pressure drives the fluid through a deforming solid. It also provides resistance to gravity load and traction applied to the boundary of the I

region being modeled.

A-9 I

SANGRE uses the

• leapfrog" time-stepping method to circumvent the riumerical instability problem caused by the large difference in rates for fluid flow and structural deformation. Although the accuracy of this method is less than that achieved by regular implicit methods that are unconditionally I

stable, it produces a 50 percent smaller system of equations for the global model. Consequently, about four times the savings in computer time can be achieved.

SANGRE models the interfaces (e.g., fault planes) as slide lines. These interfaces are characterized by reduced or no resistance to relative motion of the two surfaces along the tangential direction. The penalty function method is used to formulate the equations for dependent variables at nodes adjacent to the slide 1mes. No a priori assumption is made about the contact locations. Large relative motions of the regions can be modeled. The rezoning capability of SANGRE can be used in problems where severe mesh distortion occurs.

I Change of temperature in the modeled region occurs through heat changes in the solid as well as conduction of hea: in the solid and convection of heat by the fluid in the medium.

Material constitutive models include:

I Isotropic elastic

• Anisotropic elastic Mohr-Coulomb, and

=

Temperature-dependent material properties The available documentation does not give any information about the preprocessor and postprocessor to be used with SANGRE. The manual lacks sufficient details of m.deling techniques I.

required to use the program at its full potential. Only five verification problems are given in the manual.

{

A.13 STEALTH STEALTH (Hofmann,1981a; 1981b) is a general purpose explicit finite difference transient continuum computer code for solid, structural, and thermohydraulic analysis. It can carry out one, two,

and three-dimensional calculations. The STEALTH program is formulated based on the three conservation laws of physics and allows the simulation of "three-way* coupled thermal-hydrologic-mechanical behavior of a continuum. Since STEALTH is a continuum code, no fracture flow capability is provided. STEALTH has an automatic rezoning capability that gives it an arbitrary Lagrangian-Eulerian-like character that makes the code capable of performing large deformation calculations. Other code capabilities include:

1 A boundary interaction logic that helps defining internal sliding and debonding surfaces. The Coulomb friction law governs the sliding behavior. This boundary interaction logic does not include the ability of accepting normal and shear stiffnesses input. This logic gives a user the I

option to couple STEALTH to other types of codes.

Temperature and pressure dependency of viscosity coefficient and dilatational viscosity of fluid and gas.

Temperature and pressure dependency of shear and bulk modulus of solids.

Ability to restart with a change of material properties. One typical example of this capability is simulation of the potential effect of excavations and backfilling.

Temperature dependency of thermal conductivity and specific heat capacity.

A-10 I

Use of user provided material models, for instance, stress or pressure dependency of thermal conductivity.

Several limitations of the STEALTH code were also identified. The code is not appropriate for:

Finding " modes" in linear structural dynamics analyses I

Solving incompressible flow problems, and Efficiently solving steady-state boundary value problems and elastic static analyses The results can be analyzed with the postprocessor ADAPRO.

A.14 THAMES I

TH AMES (Thermal, IIydraulic, And Mechanical System analysis) is a finite element code for i

three-way coupled thermal, hydraulic, and mechanical processes in saturated and unsaturated geologic l

media (Ohnishi et al.,1990). The effect of discontinuities in rock mass is taken into account by using the concept of crack tensor. The coupled processes considered are based on the following assumptions:

i i

Medium is porous and elastic Darcy's law models the flow of water in saturated and unsaturated medium i

Fourier's law models the heat transfer in both solid and liquid phases with no consideration given to gaseous phase Phase change between water and vapor is not modeled, and Density of water is dependent on both temperature and pressure Heat flow causes thermal stress on the model producing deformation of the rock mass Mechanical work generates heat in the system. Heat flow causes buoyancy driven water flow in the model.

Flow of water causes heat flow by convection. Change in fluid pressure due to flow of water changes the I

effective stress in the system. Mechanical deformation changes the storativity of the rock mass.

At the time of preparing this report, the relevant manual was not available. The above information was extracted from Ohnishi et al. (1990).

A.15 UDEC (VERSION 1.8)

I The Universal Distinct Element Code (UDEC) is a two-dimensional distinct element code developed by Itasca Consulting Group, Inc. (1992a; 1992b). It is suecifically written for 80386 based

'I IBM-compatible computers, although an X-window version of the program is also available. UDEC simulates the response of discontinuous medium, such as jointed rock mass, subjected to either static or dynamic loads through an explicit solution algorithm, although a linited impF-it solution scheme is also available. UDEC is especially designed for geomechanics problems with the following features:

• The medium is repre. tented as an assemblage of either rigid or deformable discrete blocks.

The discentinuities are treated as the boundary conditions between the blocks: large displacements along the discontinuities and rotation of the blocks are permissible.

The user can develop the data structure interactively through the built-in commands or it can l

be read from an input file.

A boundary element model is available for simulating the boundary of an infinite elastic body.

l.

A-11 i

Structural elements are available to simulate rock reinforcement and interior support.

A visco-plastic model is available to simulate flow of cement grout in the discontinuity.

• The code simulates the transient flux of heat in materials and the induced thermal stresses I

(deformable blocks only). The flow of heat within the blocks is by conduction only. Heat transfer properties are independent of strain or pressure. Material properties and failure characteristics of blocks and joints are also temperature-independent.

I Heat sources can be added in real time and can decay exponentially with time.

Material constitutive models include:

I Null Isotropic elastic Drucker-Prager plasticity Ubiquitous joint Double-yield, and Strain-hardening / softening UDEC has an option for using non-reflecting boundaries in the model in dynamic analysis. This capability is very useful in wave propagation problems to simulate the passage of waves outside the modeled region. Loads, stresses, and velocities at the boundaries can be applied as a constant, linear, sinusoidal, or user-supplied function. Alternatively, loads, stresses, rnd velocities can be defined as a series of discrete poin's given in a tabular form or can be read from a file. This capability is very important in I

modeling the responses from actual time-history of applied load, stresses, or velocities.

1 1

I I

)

i I

I l

A-12 I