April 2012 

NASA/CR-2012-217562

NIA Report No. 2012-01 

 

Development and Application of Benchmark 

Examples for Mixed-Mode I/II Quasi-Static 

Delamination Propagation Predictions

Ronald Krueger 

National Institute of Aerospace, Hampton, Virginia 
NASA STI Program . . . in Profile

Since its founding, NASA has been dedicated to 

the advancement of aeronautics and space science. 

The NASA scientific and technical information (STI) 

program plays a key part in helping NASA maintain 

this important role. 

 The NASA STI program operates under the 

auspices of the Agency Chief Information Officer. It 

collects, organizes, provides for archiving, and 

disseminates NASA’s STI. The NASA STI program 

provides access to the NASA Aeronautics and Space 

Database and its public interface, the NASA Technical 

Report Server, thus providing one of the largest 

collections of aeronautical and space science STI in 

the world. Results are published in both non-NASA 

channels and by NASA in the NASA STI Report 

Series, which includes the following report types: 

TECHNICAL PUBLICATION. Reports of 

completed research or a major significant phase 

of research that present the results of NASA 

programs and include extensive data or 

theoretical analysis. Includes compilations of 

significant scientific and technical data and 

information deemed to be of continuing 

reference value. NASA counterpart of peerreviewed formal professional papers, but having 

less stringent limitations on manuscript length 

and extent of graphic presentations. 

TECHNICAL MEMORANDUM. Scientific 

and technical findings that are preliminary or of 

specialized interest, e.g., quick release reports, 

working papers, and bibliographies that contain 

minimal annotation. Does not contain extensive 

analysis. 

CONTRACTOR REPORT. Scientific and 

technical findings by NASA-sponsored 

contractors and grantees. 

CONFERENCE PUBLICATION. Collected 

papers from scientific and technical 

conferences, symposia, seminars, or other 

meetings sponsored or co-sponsored by NASA. 

SPECIAL PUBLICATION. Scientific, 

technical, or historical information from NASA 

programs, projects, and missions, often 

concerned with subjects having substantial 

public interest. 

TECHNICAL TRANSLATION. Englishlanguage translations of foreign scientific and 

technical material pertinent to NASA’s mission. 

 Specialized services also include creating custom 

thesauri, building customized databases, and 

organizing and publishing research results. 

 For more information about the NASA STI 

program, see the following: 

Access the NASA STI program home page at 

http://www.sti.nasa.gov

E-mail your question via the Internet to 

help@sti.nasa.gov

Fax your question to the NASA STI Help Desk 

at 443-757-5803 

Phone the NASA STI Help Desk at 

443-757-5802 

Write to: 

 NASA STI Help Desk 

 NASA Center for AeroSpace Information 

 7115 Standard Drive 

 Hanover, MD 21076-1320 
National Aeronautics and 

Space Administration 

Langley Research Center Prepared for Langley Research Center 

Hampton, Virginia 23681-2199 Under Cooperative Agreement NNL09AA00A 

 

April 2012 

NASA/CR-2012-217562

NIA Report No. 2012-01 

Development and Application of Benchmark 

Examples for Mixed-Mode I/II Quasi-Static 

Delamination Propagation Predictions

Ronald Krueger 

National Institute of Aerospace, Hampton, Virginia 
Available from: 

NASA Center for AeroSpace Information 

7115 Standard Drive 

Hanover, MD 21076-1320 

443-757-5802 

The use of trademarks or names of manufacturers in this report is for accurate reporting and does not constitute an 

official endorsement, either expressed or implied, of such products or manufacturers by the National Aeronautics 

and Space Administration.
DEVELOPMENT AND APPLICATION OF BENCHMARK EXAMPLES FOR MIXEDMODE I/II QUASI-STATIC DELAMINATION PROPAGATION PREDICTIONS 

Ronald Krueger*

ABSTRACT

The development of benchmark examples for quasi-static delamination propagation 

prediction is presented and demonstrated for a commercial code. The examples are based on 

finite element models of the Mixed-Mode Bending (MMB) specimen. The examples are 

independent of the analysis software used and allow the assessment of the automated 

delamination propagation prediction capability in commercial finite element codes based on 

the virtual crack closure technique (VCCT). First, quasi-static benchmark examples were 

created for the specimen. Second, starting from an initially straight front, the delamination 

was allowed to propagate under quasi-static loading. Third, the load-displacement 

relationship from a propagation analysis and the benchmark results were compared, and 

good agreement could be achieved by selecting the appropriate input parameters. Good

agreement between the results obtained from the automated propagation analysis and the 

benchmark results could be achieved by selecting input parameters that had previously 

been determined during analyses of mode I Double Cantilever Beam and mode II End 

Notched Flexure specimens. The benchmarking procedure proved valuable by 

highlighting the issues associated with choosing the input parameters of the particular 

implementation. Overall, the results are encouraging, but further assessment for mixedmode delamination fatigue onset and growth is required. 

1. INTRODUCTION

Over the past two decades, the use of fracture mechanics has become common practice to 

characterize the onset and growth of delaminations. In order to predict delamination onset or 

growth, the calculated strain energy release rate components are compared to interlaminar fracture 

toughness properties measured over a range from pure mode I loading to pure mode II loading. 

The virtual crack closure technique (VCCT) is widely used for computing energy release rates 

based on results from continuum (2D) and solid (3D) finite element (FE) analyses and to supply the 

mode separation required when using the mixed-mode fracture criterion [1, 2]. The virtual crack 

closure technique was recently implemented into several commercial finite element codes. As new 

methods for analyzing composite delamination are incorporated into finite element codes, the need 

for comparison and benchmarking becomes important since each code requires specific input 

parameters unique to its implementation. These parameters are unique to the numerical approach 

chosen and do not reflect real physical differences in delamination behavior. 

An approach for assessing the mode I, mode II and mixed-mode I and II, delamination 

propagation capabilities in commercial finite element codes under static loading was recently 

presented and demonstrated for VCCT for ABAQUS®1 [3,4] as well as MD NastranTM and MarcTM2

[5]. First, benchmark results were created manually for finite element models of the mode I Double 

 

*

R. Krueger, National Institute of Aerospace, 100 Exploration Way, Hampton, VA, 23666, resident at Durability, 

Damage Tolerance and Reliability Branch, MS 188E, NASA Langley Research Center, Hampton, VA, 23681, USA. 

1

 ABAQUS® is a product of Dassault Systèmes Simulia Corp. (DSS), Providence, RI, USA 2

 MD NastranTM and MarcTM are manufactured by MSC.Software Corp., Santa Ana, CA, USA. NASTRAN® is a 

registered trademark of NASA. 

1
Cantilever Beam (DCB), the mode II End Notched Flexure (ENF) and the mixed-mode I/II Single 

Leg Bending (SLB) specimen. Second, starting from an initially straight front, the delamination was 

allowed to propagate using the automated procedure implemented in the finite element software. 

The approach was then extended to allow the assessment of the delamination fatigue growth 

prediction capabilities in commercial finite element codes [4,6]. As for the static case, benchmark 

results were created manually first for the mode I Double Cantilever Beam (DCB) and the mode II 

End Notched Flexure (ENF) specimen. Second, the delamination was allowed to grow under cyclic 

loading in a finite element model of a commercial code. In general, good agreement between the 

results obtained from the propagation and growth analysis and the benchmark results could be 

achieved by selecting the appropriate input parameters. Overall, the results were encouraging but 

showed that additional assessment for mixed-mode delamination cases is required.

The objective of the present study was to create additional benchmark examples, independent of 

the analysis software used, which allows the assessment of the quasi-static delamination 

propagation prediction capabilities in commercial finite element codes. For the simulation of mixedmode I/II fracture, the Mixed-Mode Bending (MMB) specimen was selected as shown in Figure 1. 

Dimensions, layup and material properties were taken directly from a related experimental study 

[7].

Static benchmark results were created for three mixed-mode ratios (GII/GT =0.2, 0,5 and 0.8) 

based on the approach developed earlier [3]. To create the benchmark results, two-dimensional 

finite element models were used for simulating the MMB specimens with different delamination 

lengths a0. For each delamination length modeled, the load, P, and the displacement, w, were 

monitored. The total strain energy release rate, GT, and mixed-mode ratio, GII/GT, were calculated 

for a fixed applied displacement. It is assumed that the delamination propagates when the computed 

energy release rate, GT, reaches the mixed-mode fracture toughness Gc. Thus, critical loads and 

critical displacements for delamination propagation were calculated for each delamination length 

modeled. From these critical load/displacement results, benchmark solutions were created. It is 

assumed that the load/displacement relationship computed during automatic propagation should 

closely match the benchmark cases.

After creating the benchmark cases, the approach was demonstrated for the commercial finite 

element code ABAQUS®. Starting from an initially straight front, the delamination was allowed to 

propagate under quasi-static loading based on the algorithms implemented into the software. Input 

control parameters were varied to study the effect on the computed delamination propagation. The 

benchmark enabled the selection of the appropriate input parameters that yielded good agreement 

between the results obtained from the propagation analysis and the benchmark results. Once the 

parameters have been identified, they may then be used with confidence to model delamination 

growth for more complex configurations. 

In this paper, the development of the benchmark cases for the assessment of the quasi-static 

delamination propagation is presented. Examples of automated propagation analyses are shown, and 

the selection of the required code specific input parameters are discussed. 

2. METHODOLOGY BASED ON FRACTURE MECHANICS 

For the current numerical investigation, the Mixed-Mode Bending (MMB) specimen, as shown 

in Figure 1, was chosen since it is simple, and exhibits the mixed-mode I/II opening fracture mode 

over a wide range of mixed-mode ratios GII/GT. For the current investigation, finite element models 

were created for three distinct mixed-mode ratios (GII/GT =0.2, 0,5 and 0.8). The methodology for 

2
delamination propagation, was applied to the MMB specimen to create the benchmark example [8, 

9]. For the current study, MMB specimens made of IM7/8552 graphite/epoxy with a unidirectional 

layup, [0]24, were modeled. The material, layup, overall specimen dimensions including initial crack 

length, a0, were identical to specimens used in related experimental studies [7,10,11]. The material 

properties are given in Tables I and II.

A quasi-static mixed-mode fracture criterion is determined by plotting the interlaminar fracture 

toughness, Gc, versus the mixed-mode ratio, GII/GT as shown in Figure 2. The fracture criterion is 

generated experimentally using pure Mode I (GII/GT=0) Double Cantilever Beam (DCB) tests [11], 

pure Mode II (GII/GT=1) End-Notched Flexure (ENF) tests [10], and Mixed Mode Bending (MMB) 

tests of varying ratios of GI and GII [7]. For the material used in this study, the experimental data 

(open symbols) and mean values (filled symbols) are shown in Figure 2. A 2D fracture criterion was 

suggested by Benzeggah and Kenane [12] using a simple mathematical relationship between Gc and 

GII/GT

Gc = GIc + ( ) GIIc GIc  GII

GT

 

   

 	 




 (1) 

In this expression, GIc and GIIc are the experimentally determined fracture toughness data for mode I 

and II as shown in Figure 2. The factor was determined by a curve fit using the LevenbergMarquardt algorithm in the KaleidaGraphTM graphing and data analysis software [13]. The 

parameters GIc, GIIc and are required input to perform a VCCT analysis in ABAQUS® as 

discussed in the appendix.

3. FINITE ELEMENT MODELING 

3.1 Model description 

An example of a two-dimensional finite element model of a Mixed-Mode Bending (MMB) 

specimen with boundary conditions is shown in Figure 3 for a mixed-mode ratio GII/GT =0.2. Based 

on previous experience [3,4,6], the specimen was modeled with solid plane strain elements (CPE4I) 

in ABAQUS® Standard 6.10. Along the length, all models were divided into different sections with 

different mesh refinement as shown in Figure 3a. The MMB specimen was modeled with six 

elements through the specimen thickness (2h) as shown in the detail of Figure 3b. The resulting 

element length at the delamination tip was a=0.5 mm. The load apparatus was modeled explicitly 

using rigid beam elements (R2D2) as shown in Figure 3a. Multi-point constraints were used to 

connect the rigid elements with the planar model of the specimen and enforce the appropriate 

boundary conditions as shown in Figures 3b and c. 

The plane of delamination was modeled as a discrete discontinuity in the center of the specimen. 

For the analysis with ABAQUS® 6.10, the models were created as separate meshes for the upper 

and lower part of the specimens with identical nodal point coordinates in the plane of delamination 

[14]. Two surfaces (top and bottom surface) were defined to identify the contact area in the plane of 

delamination as shown in Figures 3b and c. Additionally, a node set was created to define the intact 

(bonded nodes) region. Two models, used to simulate mixed-mode ratios GII/GT =0.5 and 0.8, are 

shown in Figure 4. The mesh of the specimen was kept the same for all three mode ratios. Only the 

lengths of the rigid elements used to simulate the load apparatus were changed as shown in 

Figure 4. 

3
Examples of three-dimensional finite element models of the MMB specimen are shown in 

Figures 5 and 6. Along the length, all models were divided into different sections with different 

mesh refinement. A refined mesh was used in the center of the MMB specimen as shown in the 

detail of Figure 5b. Across the width, a uniform mesh (25 elements) was used to avoid potential 

problems at the transition between a coarse and finer mesh [3-6]. Through the specimen thickness 

(2h), six elements were used as shown in the detail of Figure 5b. The resulting element length at the 

delamination tip was a=0.5 mm. The specimen was modeled with solid brick elements (C3D8I), 

which had yielded excellent results in previous studies [3,4,6]. The load apparatus was modeled 

explicitly using rigid plate elements (R3D4) as shown in Figure 5a. As before for the twodimensional model, multi-point constraints were used to connect the rigid elements with the solid 

model of the specimen and enforce the appropriate boundary conditions. Two models used to 

simulate mixed-mode ratios GII/GT =0.5 and 0.8 are shown in Figure 6. The mesh of the specimen 

was kept the same for all three mode ratios. Only the lengths of the rigid plate elements used to 

simulate the load apparatus were changed as shown in Figure 6. 

3.2 Static delamination propagation analysis 

For the automated delamination propagation analysis, the VCCT implementation in ABAQUS®

Standard 6.10 was used. The plane of delamination in three-dimensional analyses is modeled using 

the existing ABAQUS®/Standard crack propagation capability based on the contact pair capability 

[14]. Additional element definitions are not required, and the underlying finite element mesh and 

model does not have to be modified [14]. The implementation offers a crack and delamination 

propagation capability in ABAQUS®. It is implied that the energy release rate at the crack tip is 

calculated at the end of a converged increment. Once the energy release rate exceeds the critical 

strain energy release rate (including the user-specified mixed-mode criteria as shown in Figure 2), 

the node at the crack tip is released in the following increment, which allows the crack to propagate. 

To avoid sudden loss of stability when the crack tip is propagated, the force at the crack tip before 

advancement is released gradually during succeeding increments in such a way that the force is 

brought to zero no later than the time at which the next node along the crack path begins to open 

[14,15].

In addition to the mixed-mode fracture criterion, VCCT for ABAQUS® requires additional input 

for the propagation analysis. If a user specified release tolerance is exceeded in an increment 

(G-Gc)/Gc > release tolerance, a cutback operation is performed which reduces the time increment. 

In the new smaller increment, the strain energy release rates are recalculated and compared to the 

user specified release tolerance. The cutback reduces the degree of overshoot and improves the 

accuracy of the local solution [14]. A release tolerance of 0.2 is suggested in the handbook [14]. To 

help overcome convergence issues during the propagation analysis, ABAQUS® provides: 

• contact stabilization which is applied across only selected contact pairs and used to 

control the motion of two contact pairs while they approach each other in multi-body 

contact. The damping is applied when bonded contact pairs debond and move away 

from each other [14,15]

• automatic or static stabilization which is applied to the motion of the entire model and is 

commonly used in models that exhibit statically unstable behavior such as buckling 

[14,15]

• viscous regularization, which is applied only to nodes on contact pairs that have just 

debonded. The viscous regularization damping causes the tangent stiffness matrix of the 

4
softening material to be positive for sufficiently small time increments. Viscous 

regularization damping in VCCT for ABAQUS® is similar to the viscous regularization 

damping provided for cohesive elements and the concrete material model in 

ABAQUS®/Standard [14,15]. Further details about the required input parameters are 

discussed in the appendix. 

For automated propagation analysis, it was assumed that the computed behavior should closely 

match the benchmark results created below. For all analyses, the elastic constants (given in Table I) 

and the input to define the fracture criterion (given in Table II) were kept constant. Based on 

previous results [3], automatic or static stabilization was not used in this study. To limit the 

scope of the study, viscous regularization was also not considered. Instead, the study focused on 

using the default values and parameter combinations that had recently yielded good results [4]. 

This approach was used to gain confidence in the parameters that had been identified previously 

and thus determine if the selection of the parameters could be considered problem independent. 

Therefore, only the following items were varied to study the effect on the automated delamination 

propagation behavior during the analysis.

• The release tolerance (reltol) was varied. 

• Analyses were performed with and without contact stabilization. For analyses that 

included contact stabilization only, a single stabilization factor (cs=1x10-6) was used. 

• Two- and three-dimensional models with different types of elements were used. 

4. DEVELOPMENT OF THE STATIC BENCHMARK CASES 

4.1 Benchmark case for 20% mode II 

The static benchmark case was created based on the approach developed earlier [3]. Twodimensional finite element models simulating MMB specimens with 17 different delamination 

lengths a0 were created (25.4 mma070.6 mm). For each delamination length modeled, the load,

P, and displacement, w, were monitored as shown in Figure 7 (colored lines) for the case of 20% 

mode II (GII/GT =0.2). Using VCCT, the total strain energy release rate, GT, and the mixed-mode 

ratio GII/GT were computed at the end of the analysis as shown in Figure 7. The failure index GT/Gc

was calculated by correlating the computed total energy release rate, GT, with the mixed-mode 

fracture toughness, Gc, of the graphite/epoxy material. As discussed before (see Figure 2), the 

mixed-mode fracture toughness, Gc, is a function of the mixed-mode ratio GII/GT. Hence, the mixedmode fracture toughness, Gc for each computed mixed-mode ratio (GII/GT 0.2) was obtained from 

the curve fit of the material data (solid red curve) as illustrated in Figure 8 (blue arrows). It is 

assumed that the delamination propagates when the failure index GT/Gc reaches unity. Therefore, 

the critical load, Pcrit, can be calculated based on the relationship between load, P, and the energy 

release rate, G [16], 

G = P2

2 

CP

A (2) 

In equation (2), CP is the compliance of the specimen, and  A is the increase in surface area 

corresponding to an incremental increase in load or displacement at fracture. The critical load, Pcrit,

and critical displacement, wcrit, were calculated for each delamination length modeled

5
GT

Gc

= P2

Pcrit

2 Pcrit = P Gc

GT

 , wcrit = w

Gc

GT

 (3) 

and the results were included in the load/displacement plots as shown in Figure 9 (solid red circles). 

These critical load/displacement results indicated that, with increasing delamination length, less 

load is required to extend the delamination up to a0 48 mm. To further extend the delamination an 

increase in load is required past w=3.5 mm. 

As shown in Figure 10, these critical load/displacement results (solid black dots) can be used to 

create a benchmark solution (solid red line) for applied displacements, w. The benchmark result 

may also be visualized by plotting the prescribed displacements, w, at delamination growth onset 

versus the increase in delamination length, a*, as illustrated in Figure 11. It is assumed that the 

load/displacement, or applied displacement/delamination-length relationship computed during 

automatic propagation should closely match the benchmark cases. 

4.2 Benchmark case for 50% mode II

The benchmark case for 50% mode II (GII/GT =0.5) was created similarly. Two-dimensional 

finite element models simulating MMB specimens with 16 different delamination lengths a0 were 

created (25.4 mma070.6 mm) to study the case of GII/GT =0.5. For each delamination length 

modeled, the load, P, and displacement, w, were monitored as shown in Figure 12 (colored lines). 

Using VCCT, the total energy release rate, GT, and the mixed-mode ratio GII/GT were computed at 

the end of the analysis as shown in Figure 12. The mixed-mode fracture toughness, Gc for each 

computed mixed-mode ratio (GII/GT 0.5) was obtained from the curve fit of the material data (solid 

red curve) as illustrated previously in Figure 8 (green arrows). The critical load, Pcrit, and critical 

displacement, wcrit, were calculated for each delamination length modeled using equation (3), and 

the results were included in the load/displacement plots as shown in Figure 13 (solid red circles). 

These critical load/displacement results indicated that, with increasing delamination length, less 

load is required to extend the delamination up to a0 48 mm. To further extend the delamination an 

increase in load is required past w=2.0 mm. For the first three delamination lengths, a0, plotted in 

Figure 13, the values of the critical displacements also decreased at the same time. This means that 

the MMB specimen exhibits unstable delamination propagation under load control as well as 

displacement control in this region for GII/GT =0.5. The remaining critical load/displacement results 

indicated stable propagation. From these critical load/displacement results (solid grey circles and 

dashed line), two benchmark solutions can be created as shown in Figure 14. During the analysis, 

either prescribed displacements, w, or nodal point loads, P, are applied. For the case of prescribed 

displacements, w, (dashed blue line), the applied displacement must be held constant over several 

increments once the critical point (Pcrit, wcrit) is reached, and the delamination front is advanced 

during these increments. Once the critical path (dashed grey line) is reached, the applied 

displacement is increased again incrementally. For the case of applied nodal point loads (dashed red 

line), the applied load must be held constant while the delamination front is advanced during these 

increments. Once the critical path (dashed grey line) is reached, the applied load is increased again 

incrementally.

The benchmark result for prescribed displacements may also be visualized by plotting the 

prescribed displacements, w, at delamination growth onset versus the increase in delamination 

length, a*, as illustrated in Figure 15. For the case of applied nodal point loads, the benchmark 

result may be visualized by plotting the applied loads, P, at delamination growth onset versus the 

increase in delamination length, a*, as illustrated in Figure 16. It is assumed that the 

6
load/displacement, load/delamination-length or applied displacement/delamination-length 

relationship computed during automatic propagation should closely match the benchmark cases. 

4.3 Benchmark case for 80% mode II

The benchmark case for 80% mode II (GII/GT =0.8) was created based on the approach 

developed earlier [3] and discussed above for 20% and 50% mode II. Two-dimensional finite 

element models simulating MMB specimens with 15 different delamination lengths a0 were created 

(25.4 mma070.6 mm) to study the case of GII/GT =0.8. For each delamination length modeled, 

the load, P, and displacement, w, were monitored as shown in Figure 17 (colored lines). Using 

VCCT, the total energy release rate, GT, and the mixed-mode ratio GII/GT were computed at the end 

of the analysis as shown in Figure 17. The mixed-mode fracture toughness, Gc for each computed 

mixed-mode ratio (GII/GT 0.8) was obtained from the curve fit of the material data (solid red 

curve) as illustrated previously in Figure 8 (purple arrows). The critical load, Pcrit, and critical 

displacement, wcrit, were calculated for each delamination length modeled, using equation (3), and 

the results were included in the load/displacement plots as shown in Figure 18 (solid red circles). 

These critical load/displacement results indicated that, with increasing delamination length, less 

load is required to extend the delamination up to a0 48 mm. To further extend the delamination an 

increase in load is required past w=2.0 mm.For the first five delamination lengths, a0, plotted in 

Figure 18, the values of the critical displacements also decreased at the same time. This means that 

the MMB specimen exhibits unstable delamination propagation under load control as well as 

displacement control in this region for GII/GT =0.8. The remaining critical load/displacement results 

indicated stable propagation. From these critical load/displacement results (solid grey circles and 

dashed line), two benchmark solutions can be created as shown in Figure 19. During the analysis, 

either prescribed displacements, w, or nodal point loads, P, are applied. For the case of prescribed 

displacements, w, (dashed blue line), the applied displacement must be held constant over several 

increments once the critical point (Pcrit, wcrit) is reached, and the delamination front is advanced 

during these increments. Once the critical path (dashed grey line) is reached, the applied 

displacement is increased again incrementally. For the case of applied nodal point loads (dashed red 

line), the applied load must be held constant while the delamination front is advanced during these 

increments. Once the critical path (dashed grey line) is reached, the applied load is increased again 

incrementally.

The benchmark result for prescribed displacements may also be visualized by plotting the 

prescribed displacements, w, at delamination growth onset versus the increase in delamination 

length, a*, as illustrated in Figure 20. For the case of applied nodal point loads, the benchmark 

result may be visualized by plotting the applied loads, P, at delamination growth onset versus the 

increase in delamination length, a*, as illustrated in Figure 21. It is assumed that the 

load/displacement, load/delamination-length or applied displacement/delamination-length 

relationship computed during automatic propagation should closely match the benchmark cases. 

5. STATIC ANALYSIS BENCHMARKING 

5.1 Computed energy release rates and mixed-mode ratios 

First, static analyses were performed with propagation disabled to compute the total energy 

release rate, GT, the mode II component, GII, and the mixed-mode ratio GII/GT. For this static 

analysis, an applied maximum displacement of w=2.0 mm was chosen for all models. To ensure that 

the model geometry and material input data for all models produced consistent results, the 

7
computed results were evaluated. The results were obtained from models shown in Figures 3 

through 6.

5.1.1 Comparison of results obtained from two- and three-dimensional models 

For all mixed-mode ratios (results for GII/GT=0.2 in blue; 0.5 in green and 0.8 in black) the twodimensional models (results with open symbols and solid lines) yield total energy release rates, GT,

that are almost identical to the results obtained from three-dimensional models (solid symbols and 

dashed lines) as shown in Figure 22. Comparing the computed mode II strain energy release rates, 

GII (plotted in Figure 23), and the mixed-mode ratios, GII/GT, (plotted in Figure 24), confirms the 

good agreement between the different models. The computed mixed-mode ratios, GII/GT, are 

slightly lower than the targeted values for the test (dashed lines) [7] a trend that appears to increase 

with increasing applied displacement, w. Based on the current findings, it was assumed that results 

obtained from two- and three-dimensional models computed later during the assessment phase were 

easily comparable. 

For each of the two-dimensional finite element models representing different delamination 

lengths (25.4 mma076.2 mm), the mixed-mode ratios corresponding to an applied displacement 

w=2.0 mm, were computed. The results were plotted versus the increase in delamination length, a*,

as illustrated in Figure 25. For all three mixed-mode ratios studied (results for GII/GT=0.2 in blue; 

0.5 in green and 0.8 in black) the computed mixed-mode ratios increase slightly with increasing 

delamination length until a*20 mm when the mixed-mode ratios reach the targeted values for the 

test (dashed lines) [7]. The computed mixed-mode ratios increased further until the delamination tip 

is located under the loading point. For delaminations extending beyond the load point, the computed 

mixed-mode ratios increased for models with target values GII/GT=0.5 and 0.8 and decreased for 

models with target values GII/GT=0.2. It should be noted, that delaminations extending beyond the 

load point could be included in the current benchmarking process, which is entirely geared toward 

numerical analysis. For actual physical tests geared towards obtaining mixed-mode fracture 

toughness values however, delaminations extending beyond the load point become meaningless. 

5.1.2 Computed distributions across the specimen width 

For all the mixed-mode ratios studied (GII/GT=0.2, 0.5 and 0.8), the total energy release rate, GT,

the mode II component, GII, the mixed-mode ratio, GII/GT, and the failure index, GT/Gc, were also 

computed and plotted versus the normalized width, y/B, as shown in Figures 26 through 29. For a 

delamination length a0= 25.4 mm the computed total energy release rates rates, GT (results for 

GII/GT=0.2 blue circles; 0.5 green squares and 0.8 black diamonds), were plotted versus the 

normalized width, y/B, of the specimen as shown in Figure 26. The results obtained from twodimensional finite element models (open symbols) are plotted at the center of the specimen y/B =0. 

The values obtained from two-dimensional and three-dimensional analyses are in good agreement 

as discussed above (see results plotted in Figures 22 through 25). The results obtained from threedimensional models (solid symbols) indicate that qualitatively, the total strain energy release rate is 

fairly constant over almost the entire width of the specimen. Only near the edges of the specimen is 

a decrease in total energy release rate observed for the cases GII/GT=0.2 and 0.5 and an increase in 

total energy release rate is observed for the mode II dominated case GII/GT=0.8. The mode II 

components, GII, for each case were plotted in Figure 27. An increase in GII is visible near the edges 

for the mode II dominated case GII/GT=0.8 and also for the equal mode case GII/GT=0.5. The mixedmode ratios, GII/GT, plotted in Figure 28 indicate a relative dominance of GII near the edges of the 

specimens compared to the distribution in the center of the specimens. For all three-dimensional 

8
models, the failure indices, GT/Gc, plotted in Figure 29, were somewhat higher across the entire 

width than the results from two-dimensional planar finite element models. The failure indices 

dropped near the edges, which suggested that delamination propagation was expected to locally lag 

behind and turn the initially straight front into a thumbnail front for GII/GT =0.2 and 0.5 similar to 

fronts observed in a Double Cantilever Beam (DCB) specimen. For GII/GT =0.8, the drop of the 

failure index near the edges was less pronounced and a straight front across the entire width was 

expected. The computed front shapes will be discussed in detail in a later section. 

5.2 Results from automated delamination propagation analysis for 20% mode II 

The propagation analysis was performed in two steps using the models shown in Figures 3 and 5 

starting from an initial delamination length, a0=25.4 mm. In the first step, a displacement of w=1.5

mm was applied which nearly equaled the critical displacement, wcrit=1.64 mm, determined earlier. 

In the second step, the applied displacement was increased to w=8.0 mm. For this second step, 

automatic incrementation was used in ABAQUS® and a small increment size (0.5% of the total 

step) was chosen at the beginning of the step. To reduce the risk of numerical instability and early 

termination of the analysis, a minimum allowed increment size (10-18 of the total step) was also 

chosen. The analysis was limited to 5000 increments.

Initially, analyses were performed using two-dimensional planar models (shown in Figure 3) 

without stabilization or viscous regularization. Based on previous results [3,4], release tolerance 

values smaller than the default value (reltol=0.2) [14] were used since these values yielded better 

results. Using reltol=0.1, the load dropped at the critical point, however the computed path stayed 

above the benchmark as shown in Figure 30 (solid blue line) where the computed resultant force 

(load P) is plotted versus the applied displacement w. As the displacement continued to increase 

over 3.0 mm, the computed load/displacement path converged to the benchmark result (solid grey 

circles and solid grey line). To improve the results, the release tolerance was decreased to 

reltol=0.01. The results (solid red line) matched the benchmark case. Due to the fine mesh, only a 

small saw-tooth pattern was observed.

Based on problems identified during previous analyses, automated or static stabilization was not 

used in this study [3]. The results computed when contact stabilizationwas added are plotted in 

Figure 31. For a small stabilization factor (cs=1x10-6) and a release tolerance suggested in the 

handbook (reltol=0.2) [14], the load decreased, and delamination propagation started shortly after 

reaching the critical point of the benchmark solution (solid green line). The load/displacement path 

then ran parallel to the constant displacement branch of the benchmark result and later converged to 

the benchmark result. To reduce the observed overshoot, the release tolerance was reduced. For a 

stabilization factor of cs=1x10-6 and release tolerance reltol=0.1 (solid blue line), the overshoot was 

reduced, and the computed load/displacement path then ran closer to the constant displacement 

branch of the benchmark result. This result was practically identical to the results discussed above 

without stabilization. As above, further reducing the release tolerance (reltol=0.01) yielded results 

that were in excellent agreement with the benchmark (solid red line) and due to the fine mesh, only 

a small saw-tooth pattern was observed. 

An alternative way to plot the benchmark is shown in Figure 32 where the applied displacement 

w is plotted versus the increase in delamination length a*. This way of presenting the results is 

shown because it may be of advantage for large structures where local delamination propagation 

may have little effect on the global stiffness of the structure and may therefore not be visible in a 

global load/displacement plot. However, extracting the delamination length a from the finite 

element results required more manual, time consuming post-processing of the results compared to 

9
the relatively simple and readily available output of nodal displacements and forces. The results 

plotted in Figure 30 are the examples that were discussed above and were shown in the global 

load/displacement plots of Figure 31. The conclusions that can be drawn from this plot are identical 

to those discussed above.

Results obtained from three-dimensional models are shown in Figures 33 and 34. Based on the 

results from two-dimensional planar models shown previously, contact stabilization was added to 

the analyses to help overcome convergence issues. The results (solid blue line) computed for a 

release tolerance value (reltol=0.1) and a small stabilization factor (cs=1x10-6) are plotted in 

Figure 33. The load dropped and delamination propagation started upon reaching the critical point 

of the benchmark solution (solid circles and solid grey line). The load/displacement path then ran 

parallel to the benchmark result but was shifted slightly towards lower loads. Deviation from the 

benchmark may be explained by the fact that the benchmark results were created using twodimensional planar finite element models of the MMB specimen. As shown in Figure 29, the threedimensional models yield a failure index which is somewhat higher across the entire width of the 

specimen compared to the results from two-dimensional planar finite element models. Therefore, 

delamination propagation is expected to start prior to the benchmark results obtained from twodimensional planar finite element analysis, which ultimately leads to a shift of the entire results plot 

towards lower loads.

For three distinct locations across the width of the specimen (two edges and center as shown 

later in Figure 53) the results from above were plotted in Figure 34, where the applied displacement 

w is plotted versus the increase in delamination length a*. The results for all three locations (solid 

green line for edge 1, solid blue line for center location and solid red line for edge 2) are in excellent 

agreement with the benchmark case (solid circles and solid grey line). Any deviation from the 

benchmark results mentioned above is less visible compared to the load/displacement plot of 

Figure 33.

5.3 Results from automated delamination propagation analysis for 50% mode II

5.3.1 Computed delamination propagation for applied displacement 

As discussed above, the propagation analysis was performed in two steps using the models 

shown in Figures 4a and 6a starting from an initial delamination length, a0=25.4 mm. In the first 

step, a displacement of w=1.2 mm was applied which nearly equaled the critical displacement, 

wcrit=1.34 mm, determined earlier. In the second step, the applied displacement was increased to 

w=8.0 mm. For this second step, automatic incrementation was used in ABAQUS® and a small 

increment size (0.5% of the total step) was chosen at the beginning of the step. To reduce the risk of 

numerical instability and early termination of the analysis, a minimum allowed increment size (10-18

of the total step) was also chosen. The analysis was limited to 5000 increments.

Initially, analyses were performed using two-dimensional planar models (shown in Figure 4a) 

without stabilization or viscous regularization. Release tolerance values equal to the default value 

(reltol=0.2) [14] and smaller were used since these values had yielded better results in the past [3,4]. 

Using the default value (reltol=0.2), the computed load and displacement exceeded the critical point 

before the initial load drop occurred and delamination propagation started (solid green line) as 

shown in Figure 35 where the computed resultant force (load P) is plotted versus the applied 

displacement w. The computed path initially stayed above the benchmark, however, as the 

displacement continued to increase over 2.5 mm, the computed load/displacement path converged 

to the benchmark result (solid grey circles and dashed grey line). To reduce the observed overshoot, 

the release tolerance was decreased. For a release tolerance reltol=0.1, the overshoot was reduced 

10
and the results improved (solid blue line) and shifted towards the benchmark case. For a release 

tolerance reltol=0.01, the analysis terminated once the critical point was reached (solid red line) due 

to convergence problems. 

The results computed when contact stabilizationwas added are plotted in Figure 36. Using 

a stabilization factor of cs=1x10-6, and a release tolerance (reltol=0.01) yielded results that were in 

excellent agreement with the benchmark (solid red line). Due to the fine mesh, only a small sawtooth pattern was observed. The result confirms previous observations where small release tolerance 

values (reltol=0.01) in conjunction with a stabilization factor of cs=1x10-6 had yielded excellent 

agreement with the benchmark results [3,4]. 

An alternative way to plot the benchmark is shown in Figure 37 where the applied displacement 

w is plotted versus the increase in delamination length a*. The results plotted in Figure 37 are the 

examples that were discussed above and were shown in the global load/displacement plots of 

Figures 35 and 36. The conclusions that can be drawn from this plot are identical to those discussed 

above.

The result obtained from a three-dimensional analysis is shown in Figure 38. Based on the 

results from two-dimensional planar models shown above, contact stabilization was added to the 

analysis to help overcome convergence issues. Since the computer time increases with tighter 

release tolerances, a tolerance value (reltol=0.1) was selected. The result (solid blue line) computed 

for a release tolerance value reltol=0.1 and a small stabilization factor cs=1x10-6 are plotted in 

Figure 38. The load dropped and delamination propagation started upon reaching the critical point 

of the benchmark solution (solid circles and dashed grey line). The load/displacement path then ran 

parallel to the benchmark result but was shifted slightly towards lower loads. Deviation from the 

benchmark may be explained by the fact that the benchmark results were created using twodimensional planar finite element models of the MMB specimen. As shown in Figure 29, the threedimensional models yield a failure index which is somewhat higher across the entire width than the 

results from two-dimensional planar finite element models. Therefore, delamination propagation is 

expected to start prior to the benchmark results obtained from two-dimensional planar finite element 

analysis, which ultimately leads to a shift of the entire results plot towards lower loads.

5.3.2 Computed delamination propagation for applied quasi-static load

The propagation analysis was performed in two steps using the models shown in Figures 4a and 

6a starting from an initial delamination length, a0=25.4 mm. In the first step, a load P=360 N was 

applied which equaled nearly the critical load, Pcrit=385 N, determined earlier. In the second step, 

the total load was increased (P=600 N). Automatic incrementation was used with a small increment 

size at the beginning (0.5% of the total step) and a very small minimum allowed increment (10-18 of 

the total step) to reduce the risk of numerical instability and early termination of the analysis. The 

analysis was limited to 5000 increments.

The same steps discussed in the section on applied displacement were followed. Initially, 

analyses were performed using two-dimensional planar models without stabilization or viscous 

regularization. For the default value reltol=0.2, the load increase stopped after reaching the critical 

point, but the analysis terminated immediately (solid green line) due to convergence problems as 

shown in Figure 39. The release tolerance was not increased as suggested by the ABAQUS® error in 

the message (.msg) file. Previous analysis had shown that, by increasing the release tolerance 

(reltol>0.2), termination of the analysis could be avoided. However, the results had not been in good 

agreement with the benchmark [3,4]. Therefore, additional stabilization had to be introduced in 

order to obtain agreement with the benchmark case. 

11
The results computed when contact stabilizationwas added are plotted in Figure 40. A 

small stabilization factor (cs=1x10-6) was used for all cases since it had yielded good results in 

previous analyses [3, 4]. Initially, the release tolerance value was set at the default value (reltol=0.2)

(solid green line). For this parameter combination, the load increased up to the critical point, and 

delamination propagation started while the load remained constant (solid green line). Also for the 

stable propagation path, the result was in good agreement with the benchmark result (solid grey 

line). Then, the release tolerance was reduced to reltol=0.1 (solid blue line) and further to 

reltol=0.01 (solid red line). For all cases, the results were in good agreement with the benchmark 

results over the entire load/displacement range. 

An alternative way to plot the benchmark is shown in Figure 41 where the applied load P is 

plotted versus the increase in delamination length a*. This way of presenting the results may be of 

advantage for large structures where local delamination propagation may have little effect on the 

global stiffness of the structure and may therefore not be visible in a global load/displacement plot. 

However, extracting the delamination length a from the finite element results required more 

manual, time consuming, post-processing of the results compared to the relatively simple and 

readily available output of nodal displacements and forces. The results plotted in Figure 41 are the 

same as those that were discussed above and were shown as global load/displacement plots in 

Figure 40. The conclusions that can be drawn from the plots in Figure 41 are identical to those 

discussed in the above for Figure 40.

A result obtained from a three-dimensional model is shown in Figures 42. Based on the results 

from two-dimensional planar models shown above, contact stabilization was added to help 

overcome convergence issues. A release tolerance value (reltol=0.1) was used which had yielded 

good results previously. The load/displacement result computed when a small stabilization factor 

(cs=1x10-6) was added is plotted in Figure 42. The load increase stopped and delamination 

propagation started (solid blue line) shortly before reaching the critical point of the benchmark 

solution (solid grey line). As mentioned earlier, deviation from the benchmark may be explained by 

the fact that the benchmark results were created using two-dimensional planar finite element models 

of the MMB specimen. As shown in Figure 29, the three-dimensional models yield a failure index 

which is somewhat higher across the entire width than the results from two-dimensional planar 

finite element models. Therefore, delamination propagation is expected to start prior to the 

benchmark results obtained from two-dimensional planar finite element analysis, which ultimately 

leads to a shift of the entire results plot towards lower loads.

5.4 Results from automated delamination propagation analysis for 80% mode II 

5.4.1 Computed delamination propagation for applied displacement 

As discussed above, the propagation analysis was performed in two steps using the models 

shown in Figures 4b and 6b starting from an initial delamination length, a0=25.4 mm. In the first 

step, a displacement of w=1.5 mm was applied which nearly equaled the critical displacement, 

wcrit=1.65 mm, determined earlier. In the second step, the applied displacement was increased to 

w=8.0 mm. For this second step, automatic incrementation was used in ABAQUS® and a small 

increment size (0.5% of the total step) was chosen at the beginning of the step. To reduce the risk of 

numerical instability and early termination of the analysis, a minimum allowed increment size (10-18

of the total step) was also chosen. The analysis was limited to 5000 increments.

Initially, analyses were performed using two-dimensional planar models (shown in Figure 4a) 

without stabilization or viscous regularization. Using the default value reltol=0.2, the load and 

displacement exceeded the critical point and the analysis terminated (solid green line) as shown in 

12
Figure 43 where the computed resultant force (load P) is plotted versus the applied displacement w.

By increasing the release tolerance to reltol=0.5, – as suggested by the ABAQUS® error in the 

message (.msg) file – it was possible to extend the analysis without an error message until the load 

drop occurred and delamination propagation started as shown in Figure 43, (solid light blue line). 

However, the load drop occurred past the critical point and the analysis terminated soon afterwards 

due to convergence problems. The release tolerance was not increased any further as suggested in 

the ABAQUS® error in the message (.msg) file. Previous analysis had shown that by increasing the 

release tolerance (reltol>0.2), termination of the analysis could be avoided, however, the results had 

not been in good agreement with the benchmark [3,4]. Therefore, it was decided to introduce 

additional stabilization to obtain better agreement with the benchmark case. 

The results computed when contact stabilizationwas added are plotted in Figure 44. A 

small stabilization factor (cs=1x10-6) was used for all cases since it had yielded good results in 

previous analyses [3, 4]. Initially, the release tolerance value was set at the default value (reltol=0.2)

(solid green line). For this parameter combination, the load increased past the critical point and 

delamination propagation started while the load dropped (solid green line). The load/displacement 

path then ran parallel to the constant displacement branch of the benchmark result and later 

followed the stable propagation branch of the benchmark result. To reduce the observed overshoot, 

the release tolerance was reduced. For a stabilization factor of cs=1x10-6 and release tolerance 

reltol=0.1 (solid blue line), the overshoot was reduced, and the computed load/displacement path 

then ran closer to the constant displacement branch of the benchmark result. However, just before 

reaching the minimum load, the displacement increased slightly followed by a small load drop. The 

results then followed the stable propagation branch of the benchmark result. Further reducing the 

release tolerance (reltol=0.01) yielded results that were in excellent agreement with the benchmark 

(solid red line). Due to the fine mesh, only a small saw-tooth pattern was observed along the stable 

path.

An alternative way to plot the benchmark is shown in Figure 45 where the applied displacement 

w is plotted versus the increase in delamination length a*. The results plotted in Figure 45 are the 

examples that were discussed above and were shown in the global load/displacement plots of 

Figure 44. The conclusions that can be drawn from this plot are identical to those discussed above.

The result obtained from a three-dimensional analysis is shown in Figure 46. Based on the 

results from two-dimensional planar models shown above, contact stabilization was added to the 

analysis to help overcome convergence issues. Since the computer time increases with tighter 

release tolerances, i.e. smaller release tolerance values, first a release tolerance value (reltol=0.1)

was selected. The result (solid blue line) computed for a release tolerance value (reltol=0.1) and a 

small stabilization factor (cs=1x10-6) are plotted in Figure 46. The load dropped and delamination 

propagation started upon reaching the critical point of the benchmark solution (solid grey circles and 

solid grey line). The load/displacement path then ran parallel to the benchmark result but was 

shifted slightly towards lower loads. The analysis was stopped after 6.81 105 CPU seconds (7.9 

days). To improve the results, the release tolerance value was set to reltol=0.01. The results (solid 

red line) are practically identical to the results obtained for reltol=0.01, however the computation 

time increased and the analysis was stopped after 1.97 106

 CPU seconds (22.8 days). Deviation 

from the benchmark may be explained by the fact that the benchmark results were created using 

two-dimensional planar finite element models of the MMB specimen. As shown in Figure 29, the 

three-dimensional models yield a failure index, which is somewhat higher across the entire width 

than the results from two-dimensional planar finite element models. Therefore, delamination 

propagation is expected to start prior to the benchmark results obtained from two-dimensional 

13
planar finite element analysis, which ultimately leads to a shift of the entire results plot towards 

lower loads. The current results for 80% mode II resemble the results obtained from the analysis of 

the ENF specimen (pure mode II) [4]. 

5.4.2 Computed delamination propagation for applied quasi-static load

The propagation analysis was performed in two steps using the models shown in Figures 4b and 

6b starting from an initial delamination length, a0=25.4 mm. In the first step, a load P=730 N was 

applied which equaled nearly the critical load, Pcrit=751 N, determined earlier. In the second step, 

the total load was increased (P=1000 N). Automatic incrementation was used with a small 

increment size at the beginning (0.5% of the total step) and a very small minimum allowed 

increment (10-18 of the total step) to reduce the risk of numerical instability and early termination of 

the analysis. The analysis was limited to 5000 increments. 

The same steps discussed in the section on applied displacement were followed. Initially, 

analyses were performed using two-dimensional planar models without stabilization or viscous 

regularization. Using the default value (reltol=0.2), the load increase stopped after reaching the 

critical point, but the analysis terminated immediately (solid green line) due to convergence 

problems, as shown in Figure 47. The release tolerance was not increased as suggested by the 

ABAQUS® error in the message (.msg) file. Previous analysis had shown that by increasing the 

release tolerance (reltol>0.2), termination of the analysis could be avoided, however, the results had 

not been in good agreement with the benchmark [3,4]. Therefore, it was decided to introduce 

additional stabilization to obtain better agreement with the benchmark case. 

The results computed when contact stabilizationwas added are plotted in Figure 48. A 

small stabilization factor (cs=1x10-6) was used for all cases since it had yielded good results in 

previous analyses [3, 4]. Initially, the release tolerance value was set at the default value (reltol=0.2)

(solid green line). For this parameter combination, the load increased up to the critical point and 

delamination propagation started while the load remained constant (solid green line). Also for the 

stable propagation path, the result was in good agreement with the benchmark result (solid grey 

line). Then, the release tolerance was reduced to reltol=0.1 (solid blue line) and further to 

reltol=0.01 (solid red line). For all cases, the results were in good agreement with the benchmark 

results over the entire load/displacement range. 

An alternative way to plot the benchmark is shown in Figure 49 where the applied load P is 

plotted versus the increase in delamination length a*. This way of presenting the results may be of 

advantage for large structures where local delamination propagation may have little effect on the 

global stiffness of the structure and may therefore not be visible in a global load/displacement plot. 

The results plotted in Figure 49 are the same as those that were discussed above and were shown as 

global load/displacement plots in Figure 48. The conclusions that can be drawn from the plots in 

Figure 49 are identical to those discussed in the above for Figure 48.

A result obtained from a three-dimensional model is shown in Figures 50. Based on the results 

from two-dimensional planar models shown above, contact stabilization was added to help 

overcome convergence issues. Since the computer time increases with tighter release tolerances, i.e. 

smaller release tolerance values, first a release tolerance value reltol=0.1 was selected. The result 

(solid blue line) computed for a release tolerance value (reltol=0.1) and a small stabilization factor 

(cs=1x10-6) are plotted in Figure 50. The load increase stopped (solid blue line) and delamination 

propagation started shortly before reaching the critical point of the benchmark solution (solid grey 

line). The load/displacement path then ran parallel to the benchmark result but was shifted slightly 

towards lower loads. To improve the results, the release tolerance value was set to reltol=0.01. The 

14
results (solid red line) are practically identical to the results obtained for reltol=0.1, however the 

computational time increased. As discussed earlier, deviation from the benchmark may be explained 

by the fact that the benchmark results were created using two-dimensional planar finite element 

models of the MMB specimen. As shown in Figure 29, the three-dimensional models yield a failure 

index which is somewhat higher across the entire width than the results from two-dimensional 

planar finite element models. Therefore, delamination propagation is expected to start prior to the 

benchmark results obtained from two-dimensional planar finite element analysis, which ultimately 

leads to a shift of the entire results plot towards lower loads. The current results for 80% mode II 

resemble the results obtained from the analysis of the ENF specimen (pure mode II) [4].

The results from above were plotted for three distinct locations across the width of the specimen 

(two edges and center as shown later in Figure 56) as shown in Figure 51, where the applied load P

is plotted versus the increase in delamination length a*. This way of presenting the results is shown 

since it may be of advantage for large structures where local delamination propagation may have 

little effect on the global stiffness of the structure, and may therefore not be visible in a global 

load/displacement plot. However, extracting the delamination length a from finite element results 

required more manual, time consuming post-processing of the results - especially for the threedimensional models - compared to the relatively simple and readily available output of nodal 

displacements and forces. 

Initially during the unstable growth phase, the results for both simulations (reltol=0.1, cs=1x10-6

solid lines and reltol=0.01, cs=1x10-6 dashed lines) and all three locations (green line for edge 1, 

blue line for center location and red line for edge 2) were practically identical. All the result curves 

ran parallel to the benchmark result (solid grey line) but were shifted slightly towards lower loads. 

The unstable phase, however created longer delaminations compared to the benchmark result. For 

the stable path of delamination propagation, small differences can be observed between the two 

simulations (solid versus dashed lines). In this respect, more detailed information about the 

propagation can be obtained from the plot where the applied load is plotted versus the increase in 

delamination length than from the traditional load/displacement plot. The current results for 80% 

mode II resemble the results obtained from the analysis of the ENF specimen (pure mode II) [4]. 

5.5 Computed delamination front shape 

Besides matching the load displacement behavior of benchmark results, a delamination 

propagation analysis should also yield a delamination front shape that is representative of the actual 

failure. An example of delamination front shapes observed by opening tested MMB specimens are 

shown in Figure 52 [7]. From the initial straight delamination front, which is formed by the edge 

of the polytetrafluoroethylene (PTFE) insert, the delamination develops into a curved thumbnail 

shaped front (blue line) as shown in Figure 52a for 20% mode II. The front is similar to fronts 

observed in the mode I DCB specimen [3]. A pronounced thumbnail shaped front is observed for 

the pre-cracks (red line) shown in Figures 52 b and c for 50% and 80% mode II, respectively. The 

final front (blue line) exhibits less curvature as would be expected from the increased mode II 

contribution. The fronts are somewhat jagged, suggesting that growth happens in one location then 

stops and continues at another location across the width. On average, however, the delaminations 

appear to grow uniformly across the width. 

A straight front across almost the entire width of the specimen is to be expected when looking at 

the failure index distribution plotted in Figure 29. However, due to the drop of the failure index near 

the edges, delamination propagation is expected to locally lag behind and turn the initially straight 

front into a thumbnail front for GII/GT =0.2 and 0.5 similar to fronts observed in a Double Cantilever 

15
Beam (DCB) specimen. For GII/GT =0.8, the drop of the failure index near the edges is less 

pronounced, and a straight front across the entire width is expected.

Delamination propagation computed using the model with a uniform mesh across the width 

(Figure 5) is shown in Figure 53 for GII/GT =0.2. Plotted on the bottom surface (defined in Figure 5) 

are the contours of the bond state, where the delaminated section appears in blue and the intact 

(bonded) section in red. The transition between the colors (in green/yellow) indicates the location of 

the delamination front. The initial straight front is shown in Figure 53a and was added for 

clarification in Figures 53d and e. The first propagation was observed towards the edges of the 

specimen (Figure 53b) as expected from the distribution of the failure index (Figure 29) where two 

local maxima are observed at y/B ± 0.33. Then, as expected, the front advanced more in the center 

and lagged behind at the edges as shown in Figure 53c. Later in the analysis, the delamination 

continued to propagate through the entire specimen basically as a straight front as shown in 

Figures 53d and 53e. Locally however, always one element at each edge lagged behind as shown in 

Figures 53c and 53d. After the shape stabilized, new propagation usually started in the center first as 

shown in Figure 53e, the exact location being determined by the maximum of the failure index 

distribution across the width.

Delamination fronts computed using the model with a uniform mesh across the width (Figure 5) 

are shown in Figures 54 and 55 for GII/GT =0.5. For the case of applied displacement, w, the 

propagation was observed first in the center section of the specimen (Figure 54a) rather than at the 

location of the two local maxima in the failure index distribution (Figure 29) observed at y/B ± 0.33. 

Then, the front advanced across the entire width but lagged behind at the edges as shown in 

Figure 54b. Later, the delamination propagated across the width of the specimen as a straight front 

shape with only one element at each edge lagging behind as shown in Figure 54c. The final front as 

shown in Figure 54d, appeared jagged locally but on average moved forward as a straight front 

across the entire width, which may be explained by the increase in mode II for longer delamination 

lengths (increased GII/GT as shown in Figure 25). Delaminations under pure mode II generally 

exhibit straight fronts [4,10]. For the case of applied load, P, the propagation was observed first in 

the center section of the specimen (Figure 55a). Then, the front advanced across the entire width but 

lagged behind at the edges as shown in Figure 55b. Later, the delamination appeared jagged locally 

but on average moved forward as a straight front across the entire width as shown in Figures 55c 

and 55d, which may again be explained by the increase in mode II for longer delamination lengths 

(increased GII/GT as shown in Figure 25). 

Delamination propagation computed using the model with a uniform mesh across the width 

(Figure 5) is shown in Figure 56 for GII/GT =0.8. The initial straight front is shown in Figure 56a 

and 56b and was added for clarification in Figures 56d and 56e. The first propagation was observed 

towards the edges of the specimen (Figure 56c) as expected from the distribution of the failure 

index (Figure 29) where two local maxima are observed at y/B ± 0.33. Then, the delamination 

appeared jagged locally but on average moved forward as a straight front across the entire width as 

shown in Figures 56d, which may be explained by the already high mode II ratio (GII/GT =0.8) for 

this configuration. The front remained jagged until the end of the analysis as shown in Figure 56e. 

The current results for 80% mode II resemble the results obtained from the analysis of the ENF 

specimen (pure mode II) [4]. 

Computed delamination fronts obtained from three-dimensional solid models generally matched 

the experimentally observed shapes shown in Figure 52. Improved results may be obtained by 

refining the mesh across the width. 

16
6. SUMMARY AND CONCLUSIONS 

The development of benchmark examples, which allow the assessment of the static 

delamination propagation prediction capabilities was presented and demonstrated for ABAQUS®

Standard. The examples are based on finite element models of the mixed-mode I/II Mixed-Mode 

Bending (MMB) specimen. The models are independent of the analysis software used and allow the 

assessment of the automated delamination growth prediction capabilities in commercial finite 

element codes based on the virtual crack closure technique (VCCT).

First, quasi-static benchmark results were created based on the approach developed in reference 

3. Two-dimensional finite element models were used for simulating the MMB specimens with 

different mixed-mode ratio (20%, 50% and 80% mode II) and different initial delamination lengths, 

a0. For each delamination length modeled, the load and displacements were monitored. The mixedmode I/II strain energy release rate was calculated for a fixed applied displacement. It was assumed 

that the delamination propagated when the total strain energy release rate reached the fracture 

toughness value. Thus, critical loads and critical displacements for delamination propagation were 

calculated for each initial delamination length modeled. From these critical load/displacement 

results, benchmark solutions were created. It was assumed that the load/displacement relationship 

computed during automatic propagation should closely match the benchmark cases.

After creating the benchmark cases, the approach was demonstrated for the commercial finite 

element code ABAQUS®. Starting from an initially straight front, the delamination was allowed to 

propagate under quasi-static loading based on the algorithms implemented into the software. Input 

control parameters were varied to study the effect on the computed delamination propagation and 

growth.

The results showed the following: 

• The benchmarking procedure proved valuable by highlighting the issues associated 

with choosing the input parameters of the particular implementation. 

• In general, good agreement between the results obtained from the automated 

propagation analysis and the benchmark results could be achieved by selecting input

parameters that had previously been determined during analyses of mode I Double 

Cantilever Beam and mode II End Notched Flexure specimens. 

• In particular, the results for automated delamination propagation analysis under quasistatic loading showed the following: 

o Using the default release tolerance (reltol=0.2) as suggested in the ABAQUS®

handbook or increasing the value, as suggested in the user's manual, may help 

to overcome convergence problems, however, it leads to an undesired 

overshoot of the computed result compared to the benchmark. 

o A combination of release tolerance and contact stabilization is required to 

obtain more accurate results. 

o A gradual reduction of the release tolerance and contact stabilization over 

several analyses is suggested. 

o Good agreement between analysis results and the benchmarks could be achieved 

for release tolerance values (reltol<0.1) in combination with contact stabilization

(cs=1x10-6). 

o Computed delamination fronts obtained from three-dimensional solid models 

generally matched the experimentally observed shapes. Improved results may be 

obtained by refining the mesh across the width. 

17
o The analyses for three-dimensional models of a simple MMB specimen required 

several days of computation time. Improvements to the implementation are 

required to make analysis of larger structural components computationally 

affordable.

Overall, the results are promising and the current findings concur with previously published 

conclusions [3,4,5]. However, further assessment for mixed-mode delamination fatigue onset and 

growth is required. Additional studies should also include the assessment of the propagation 

capabilities in more complex mixed-mode specimens and on a structural level. 

Assessing the implementation in one particular finite element code illustrated the value of 

establishing benchmark solutions since each code requires specific input parameters unique to its 

implementation. Once the parameters have been identified, they may then be used as a starting point 

to model delamination growth for more complex configurations. 

ACKLOWLEDGEMENTS

This research was supported by the Subsonic Rotary Wing Project as part of NASA’s 

Fundamental Aeronautics Program. 

The analyses were performed at the Durability, Damage Tolerance and Reliability Branch at 

NASA Langley Research Center, Hampton, Virginia, USA. 

REFERENCES

1. E. F. Rybicki and M. F. Kanninen, "A Finite Element Calculation of Stress Intensity Factors by a Modified 

Crack Closure Integral," Eng. Fracture Mech., Vol. 9, pp. 931-938, 1977. 

2. R. Krueger, "Virtual Crack Closure Technique: History, Approach and Applications," Applied Mechanics 

Reviews, Vol. 57, pp. 109-143, 2004. 

3. R. Krueger, "An Approach to Assess Delamination Propagation Simulation Capabilities in Commercial Finite 

Element Codes," NASA/TM-2008-215123, 2008. 

4. R. Krueger, "Development and Application of Benchmark Examples for Mode II Static Delamination 

Propagation and Fatigue Growth Predictions," NASA/CR-2011-217305, NIA report no. 2011-02, 2011. 

5. A. C. Orifici and R. Krueger, "Assessment of Static Delamination Propagation Capabilities in Commercial 

Finite Element Codes Using Benchmark Analysis," NASA/CR-2010-216709, NIA report no. 2010-03, 2010. 

6. R. Krueger, “Development of a Benchmark Example for Delamination Fatigue Growth Prediction,” 

NASA/CR-2010-216723, NIA report no. 2010-04, 2010. 

7. J. G. Ratcliffe, personal communication. 

8. I. S. Raju and T. K. O'Brien, "Fracture Mechanics Concepts, Stress Fields, Strain Energy Release Rates, 

Delamination and Growth Criteria," in Delamination Behavior of Composites, S. Sridharan, Ed.: Woodhead 

Publishing in Materials, 2008. 

9. "Composite Fatigue Damage Onset," in Composite Materials Handbook CMH-17, Vol. 3, Section 12.6.4, 

draft of revision G, 2009. 

10. T. K. O'Brien, W. M. Johnston, and G. Toland, “Mode II Interlaminar Fracture Toughness and Fatigue 

Characterization of a Graphite Epoxy Composite Material,” NASA/TM-2010-216838, 2010. 

11. G. B. Murri, personal communication. 

12. M. L. Benzeggagh and M. Kenane, "Measurement of Mixed-Mode Delamination Fracture Toughness of 

Unidirectional Glass/Epoxy Composites with Mixed-Mode Bending Apparatus," Composites Science and 

Technology, Vol. 56, pp. 439-449, 1996. 

13. KaleidaGraph: Version 4.1, 2009.

14. Abaqus Analysis User’s Manual, ABAQUS® Standard, Version 6.10, DSS Simulia, 2010. 

15. Abaqus Theory Manual, ABAQUS® Standard, Version 6.10, DSS Simulia, 2010. 

16. D. Broek, The Practical Use of Fracture Mechanics, Kluwer Academic Publishers, 1991.

18
TABLE I. MATERIAL PROPERTIES [4]. 

Unidirectional Graphite/Epoxy Prepreg 

E11 = 161 GPa E22 = 11.38 GPa E33 = 11.38 GPa 

12 = 0.32 13 = 0.32 23 = 0.45 

G12 = 5.2 GPa G13 = 5.2 GPa G23 = 3.9 GPa 

The material properties are given with reference to the ply coordinate axes where index 11 denotes the ply principal 

axis that coincides with the direction of maximum in-plane Young’s modulus (fiber direction). Index 22 denotes the 

direction transverse to the fiber in the plane of the lamina and index 33 the direction perpendicular to the plane of the 

lamina. 

TABLE II. FRACTURE PARAMETERS. 

Fracture Toughness Data [7,10,11] – Figures 2, A1 

GIc = 0.212 kJ/m2 GIIc = GIIIc =0.774 kJ/m2 = 2.1 

19
APPENDIX

Delamination propagation analysis in ABAQUS®

Delamination propagation at the interfaces in laminated composites can be simulated in 

ABAQUS®[14]. The interface along which the delamination (or crack) propagates must be 

indicated in the model using a fracture criterion definition. The fracture energy release rates at 

the crack tips in the interface elements are calculated based on the virtual crack closure technique 

(VCCT) [14,15].

Required input for ABAQUS®

The input required to perform a delamination propagation analysis under quasi-static load in 

ABAQUS® Standard is discussed in the following paragraphs. It is assumed that the reader is 

familiar with ABAQUS® Standard and the syntax used in the input file (.inp). The focus is 

therefore on the specific input that relates to delamination propagation in ABAQUS®. An 

example input file is given at the end of this appendix to provide an overview of an entire 

analysis and assist the readers in creating their own analyses. 

Input for delamination propagation 

The interface along which the delamination (or crack) propagates must be indicated in the 

model using a fracture criterion definition:

*DEBOND,SLAVE=VCCT_TOP,MASTER=VCCT_BOT,FREQ=1 

*FRACTURE CRITERION,TYPE=VCCT,MIXED MODE BEHAVIOR=BK, TOLERANCE=<reltol> 

 <GIc>,<GIIc>,<GIIIc>,<eta> 

*CONTACT CONTROLS, STABILIZE, MASTER=VCCT_BOT, SLAVE=VCCT_TOP 

 1.E-06, 0, 0.1 

where VCCT_TOP and VCCT_BOT are interface surfaces as shown in Figure A1, and <reltol> is

the release tolerance within which the crack propagation criterion must be satisfied. The critical

energy release rates <GIc>,<GIIc>,<GIIIc> and the curve fit parameter <eta> are obtained from 

the mixed-mode failure criterion as shown in Figure A2. A quasi-static mixed-mode fracture 

criterion is determined for a material by plotting the interlaminar fracture toughness, Gc, versus the 

mixed-mode ratio, GII/GT as discussed earlier and shown in Figure A2. The input for contact 

stabilization is provided in the contact controls definition. 

Example input files

An input file is given to provide an overview of an entire analysis and assist the readers in 

creating their own analyses. The analysis was divided into two steps. In the first step, a 

displacement just below the critical value, wcrit, was applied to avoid delamination propagation in 

this step. Thus, large increments could be used to get up to the critical point. In the second step, 

the final desired displacement was defined. Automatic incrementation was used with a small 

increment size at the beginning (0.5% of the total step) and a very small minimum allowed 

increment (10-18 of the total step) to reduce the risk of numerical instability and early termination of 

the analysis. The analysis was limited to 5000 increments. 

For all analyses, the input to define the fracture criterion (<GIc>,<GIIc>,<GIIIc>,<eta>),

were kept constant. The ABAQUS® keywords shown in bold type were discussed in detail in the 

previous paragraphs. 

20
Input file for fatigue onset and growth analysis 

*HEADING

MMB-UD-IM7/8552, a=30 mm, ratio=0.2

units: mm, N, MPa

*** elements, nodes, material, etc 

*ELEMENT, TYPE=CPE4I , ELSET=EALL 

...

...

*NSET, NSET=BONDED, GENERATE 

 1069, 1533, 8 

**

  surface and contact definition ** ** 

*SURFACE, TYPE=ELEMENT, NAME=VCCT_BOT 

EL_BOT, S3 

*SURFACE, TYPE=ELEMENT, NAME=VCCT_TOP 

EL_TOP, S1 

*INITIAL CONDITION, TYPE=CONTACT 

VCCT_TOP, VCCT_BOT, BONDED 

*CONTACT PAIR, INTERACTION=VCCT, ADJUST=BONDED 

VCCT_TOP, VCCT_BOT 

**** VCCT -- STAGE I **** 

*SURFACE INTERACTION, NAME=VCCT 

 <width> 

**** rigid elements to simulate MMB fixture 

*** extra nodes for rigid elements

*NSET, NSET=MMBREF, UNSORTED 

**** MMB fixture loading 

 100006, 

**

*NSET, NSET=NPIN, UNSORTED 

**** MMB rollers 

 100001, 

**

*NODE, NSET=NFIXTURE 

100001, 2.25, 152.399994, 0. 

100002, 34., 152.399994, 0. 

100003, 2.25, 102.000168, 0. 

100004, 34., 102.000168, 0. 

** ratio=0.2 

100005, 34., 9.05, 0., 

100006, 3.25, 9.05, 0., 

** ratio=0.5 

**100005, 34., 60.72, 0., 

**100006, 3.25, 60.72, 0., 

** ratio=0.8 

**100005, 34., 74.75, 0., 

**100006, 3.25, 74.75, 0., 

*** rigid beam elements 

*ELEMENT, TYPE=R2D2, ELSET=MMB_RIG 

100001, 8, 100002 

100002, 100002, 100004 

**100003, 272, 100004 

100003, 100003, 100004 

100004, 100004, 100005 

100005, 100005, 100006 

21
*RIGID BODY, ELSET=MMB_RIG, REF NODE=MMBREF, PIN NSET=NPIN 

**

** mode II loading 

*EQUATION

2

616,1,1.0,100003,1,-1.0

**** VCCT input 

*parameter

** fracture toughness: 

 GIc = 0.212 

 GIIc = 0.774 

 GIIIc = 0.774 

** B-K parameter: 

 eta=2.1 

** Damage and tolerance parameters 

damv=0.0

reltol=0.01

** width in the plane stress/strain direction 

 width =25.4 

***

********************** STEP 1 ******************************** 

*STEP, NLGEOM, INC=500 

*STATIC

 0.1, 1., 0.000000001, 1. 

**** MMB loading 

*** applied displacement w for MMB, 2D-full model

*BOUNDARY

 100006, 1, 1, -1.5 

**

*CONTROLS,PARAMETERS=TIME INCREMENTATION 

 ,,,,,,,50 

*** fracture analysis using VCCT 

*DEBOND,SLAVE=VCCT_TOP,MASTER=VCCT_BOT,FREQ=1

*FRACTURE CRITERION,TYPE=VCCT,MIXED MODE BEHAVIOR=BK, TOLERANCE=<reltol> 

 <GIc>,<GIIc>,<GIIIc>,<eta> 

**

***CONTACT CONTROLS, STABILIZE, MASTER=VCCT_BOT, SLAVE=VCCT_TOP 

** 1.E-06, 0, 0.1 

**

*OUTPUT, FIELD, VARIABLE=PRESELECT, FREQ=1 

*CONTACT OUTPUT, MASTER=VCCT_BOT, SLAVE=VCCT_TOP 

 dbt,dbsf,dbs,openbc,crsts,enrrt,efenrrtr,bdstat 

**

*Output, history,VARIABLE=PRESELECT,freq=1 

*NODE output,NSET=MMBREF 

 RF1 

*NODE output,NSET=MMBREF 

 U1 

***NODE output,NSET=NSUP 

** RF1, RF2 

***NODE output,NSET=NSUP 

** U1, U2 

*CONTACT OUTPUT, MASTER=VCCT_BOT, SLAVE=VCCT_TOP, NSET=BONDED 

 enrrt 

*END STEP 

********************** STEP 2 ******************************** 

*STEP, NLGEOM, INC= 5000 

22
*STATIC

 0.00500000 , 1.0000 , 1.0000000E-18, 1.0

**** MMB loading 

*** applied displacement w for MMB, 2D-full model

*BOUNDARY

 100006, 1, 1, -8.0 

**

*CONTROLS,PARAMETERS=TIME INCREMENTATION 

 ,,,,,,,50 

*CONTACT PRINT 

*** fracture analysis using VCCT 

***DEBOND,SLAVE=VCCT_TOP,MASTER=VCCT_BOT,FREQ=1

*DEBOND,SLAVE=VCCT_TOP,MASTER=VCCT_BOT,FREQ=1

*FRACTURE CRITERION,TYPE=VCCT,MIXED MODE BEHAVIOR=BK, TOLERANCE=<reltol> 

 <GIc>,<GIIc>,<GIIIc>,<eta> 

*CONTACT CONTROLS, STABILIZE, MASTER=VCCT_BOT, SLAVE=VCCT_TOP 

 1.E-06, 0, 0.1 

**

*OUTPUT, FIELD, VARIABLE=PRESELECT, FREQ=1 

*CONTACT OUTPUT, MASTER=VCCT_BOT, SLAVE=VCCT_TOP 

 dbt,dbsf,dbs,openbc,crsts,enrrt,efenrrtr,bdstat 

*Output, history,VARIABLE=PRESELECT,freq=1 

*NODE output,NSET=MMBREF 

 RF1 

*NODE output,NSET=MMBREF 

 U1 

*CONTACT OUTPUT, MASTER=VCCT_BOT, SLAVE=VCCT_TOP, NSET=BONDED 

 enrrt 

*END STEP 

23
Figure 1. Mixed Mode Bending specimen (MMB) (dimensions from [7]).

layup: IM7/8552 [0]24

dimensions

b 25.4 mm

2h 4.5 mm

2L 100.8 mm

2L* 152.4 mm

a0 25.4 mm

P

2L

L

a0

w

c

base

loading 

arm

specimen 2h

width b

2L*

Figure 1. Mixed Mode Bending specimen (MMB) (dimensions from [7]).

layup: IM7/8552 [0]24

dimensions

b 25.4 mm

2h 4.5 mm

2L 100.8 mm

2L* 152.4 mm

a0 25.4 mm

P

2L

L

a0

w

c

base

loading

arm

specimen 2h

width b

2L*

0

0.2

0.4

0.6

0.8

1

1.2

1.4

1.6

0 0.2 0.4 0.6 0.8 1

experimental data [7,10,11]

averaged data and curve fit

DCB, Mode I MMB, Mode I and II

Gc ,

kJ/m2

mixed mode ratio GII/GT

ENF, Mode II

G

Ic

 =0.212

 GIIc = 0.774

curve fit: Gc = G

Ic

 + (GIIc-GIc)(GII/GT)

= 2.1

Figure 2. Mixed mode fracture criterion for IM7/8552.

24
Figure 3. Deformed 2D FE-model of a MMB specimen with mixed mode ratio GII /GT=0.2.

(c). Detail of load introduction zone.

rigid elements to 

simulate load 

apparatus

(a). View of full model - GII /GT=0.2, c=92.9 mm.

(b). Detail of specimen tip and crack tip zone. 

bonded nodes

initial delamination 

front

 bottom surface

top surface

w

2L

a0

w=0

2h

a

c

u,w=0

multi-point constraints 

to enforce same w 

displacement and allow u 

displacement (sliding)

x,u

z,w

x,u

z,w

rigid elements to 

simulate load 

apparatus

rigid elements to 

simulate load apparatus

rigid elements to 

simulate load apparatus

multi-point constraints 

to enforce same u and 

w displacement and 

allow free rotation

bonded nodes

a=0.5 mm

25
Figure 4. Deformed 2D FE models of MMB specimens with different mixed-mode ratios.

(a). Mixed mode ratio GII/GT=0.5, c=41.3 mm. 

(b). Mixed mode ratio GII/GT=0.8, c=27.3 mm.

rigid elements to 

simulate load 

apparatus

rigid elements to 

simulate load 

apparatus

w

c

2L

a0

rigid elements to 

simulate load 

apparatus

rigid elements to 

simulate load 

apparatus

c

w

2L

a0

26
rigid elements to 

simulate load 

apparatus

w

c

Figure 5. Deformed 3D FE-model of a MMB specimen with mixed mode ratio GII /GT=0.2.

(a). View of full model - GII /GT=0.2, c=92.9 mm.

(b). Detail of specimen tip and crack tip zone. 

2L

B

bonded

nodes

delamination front

x y

z a0

delamination front

bonded nodes

a=0.5 mm

a

 bottom surface

top surface

 bottom surface

top surface

27
Figure 6. Deformed 3D FE-models of a MMB specimen with different mixed mode ratios.

(a). Mixed mode ratio GII /GT=0.5, c=41.3 mm. 

(b). Mixed mode ratio GII /GT=0.8, c=27.3 mm.

rigid elements to 

simulate load 

apparatus

c

x

y

z

x

y

z

rigid elements to 

simulate load 

apparatus

c

delamination front

delamination front

w

w

28
0

0.2

0.4

0.6

0.8

1

0 0.2 0.4 0.6 0.8 1

Gc

 ,

kJ/m2

mixed mode ratio GII /GT

curve fit: Gc = G

Ic

 + (GIIc-GIc )(GII /GT ) 

= 2.1

GII /GT G =0.8 II /GT G =0.5 II /GT=0.2

Gc (GII /GT=0.8)

Gc (GII /GT=0.5)

Gc (GII /GT

=0.2)

Figure 8. Use of mixed-mode fracture criterion to obtain Gc for different mode ratios GII/G

T .

G

Ic

 =0.212 GIIc = 0.774

0.0 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0 9.0 10.0 11.0

0

20

40

60

80

100

120

140

160

180

25.4

25.9

26.9

27.9

28.9

30.4

35.4

40.4

45.4

46.9

48.4

49.9

51.4

55.4

60.5

65.5

70.6

load P, N

applied displacement w, mm

G

T

 in kJ/m2 delamination length a0

 in mm

0.34

0.32

0.44

0.41

0.32 0.37 0.38

0.31

Figure 7. Load-displacement behavior for MMB specimens (20% mode II)

with different delamination lengths a0.

29
0.0 2.0 4.0 6.0 8.0 10.0

0

20

40

60

80

100

120

140

160

critical, G

T

=Gc 

benchmark

load P, N

applied displacement w, mm

Figure 10. Calculated critical behavior and resulting benchmark case for 

applied displacement, w, for a MMB specimen (20% mode II).

wcrit=1.64

Pcrit

=128.5 

0.0 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0 9.0 10.0 11.0

0

20

40

60

80

100

120

140

160

180

critical, G

T

=Gc 

load P, N

applied displacement w, mm

G

T

 in kJ/m2

0.34

0.32

0.44

0.41

0.32 0.37 0.38

0.31

Figure 9. Calculated critical load-displacement behavior for a MMB specimen (20% mode II).

30
0.0 2.0 4.0 6.0 8.0 10.0

0

50

100

150

200

250

300

350

400

450

500

550

600

650

700

750

25.4

26.9

27.9

28.9

30.9

35.4

40.4

45.4

46.9

48.4

49.9

51.4

55.4

60.5

65.5

70.6

load P, N

applied displacement w, mm

G

T

 in kJ/m2

0.72

0.63

1.6

0.73

a

0

 in mm

0.87

0.7

1.0

Figure 12. Load-displacement behavior for MMB specimens (50% mode II)

with different delamination lengths a0

.

0 5 10 15 20 25 30 35 40 45 50 55

0

2

4

6

8

10

critical, G

T

=Gc 

benchmark

applied

displacement

w, mm

increase in delamination length a*, mm

G=Gc

wcrit=1.64

Figure 11. Benchmark case for applied displacement, w, plotted versus 

increase in delamination length, a*, for a MMB specimen (20% mode II).

31
0.0 2.0 4.0 6.0 8.0 10.0

0

100

200

300

400

500

600

critical, G

T

=Gc 

benchmark for applied load P

benchmark for applied displacement w

load P, N

displacement w, mm

Figure 14. Calculated critical behavior and resulting benchmark cases for 

applied load, P, and displacement, w, for a MMB specimen (50% mode II).

wcrit=1.34

Pcrit

=385

0.0 2.0 4.0 6.0 8.0 10.0

0

50

100

150

200

250

300

350

400

450

500

550

600

650

700

750

critical, G

T

=Gc 

load P, N

applied displacement w, mm

G

T

 in kJ/m2

0.72

0.63

1.6

0.73

0.87

0.7

1.0

Figure 13. Calculated critical load-displacement behavior for a MMB specimen (50% mode II).

32
0 5 10 15 20 25 30 35 40 45 50 55

0

100

200

300

400

500

600

critical, G

T

=Gc 

benchmark

applied

load 

P, N

increase in delamination length a*, mm

Figure 16. Benchmark case for applied load, P, plotted versus increase 

in delamination length, a*, for a MMB specimen (50% mode II).

Pcrit=385

0 5 10 15 20 25 30 35 40 45 50 55

0

1

2

3

4

5

6

7

8

critical, G

T

=Gc 

benchmark

applied

displacement

w, mm

increase in delamination length a*, mm

Figure 15. Benchmark case for applied displacement, w, plotted versus 

increase in delamination length, a*, for a MMB specimen (50% mode II).

wcrit=1.34

33
0.0 2.0 4.0 6.0 8.0 10.0

0

500

1000

1500

critical, G

T

=Gc 

load P, N

applied displacement w, mm

G

T

in kJ/m2

1.8

1.7

1.7

0.86

0.8

0.8

0.8

Figure 18. Calculated critical load-displacement behavior for a MMB specimen (80% mode II).

0.0 2.0 4.0 6.0 8.0 10.0

0

500

1000

1500

25.4

26.9

27.9

28.9

30.4

35.4

40.4

45.4

48.4

49.9

51.4

55.4

60.5

65.5

70.6

load P, N

applied displacement w, mm

a

0

in mm

G

T

in kJ/m2

1.8

1.7

1.7

0.86

0.8

0.8

0.8

Figure 17. Load-displacement behavior for MMB specimens (80% mode II)

with different delamination lengths, a0.

34
0 5 10 15 20 25 30 35 40 45 50 55

0

1

2

3

4

5

6

7

8

critical, G

T

=Gc 

benchmark

applied

displacement

w, mm

increase in delamination length a*, mm

wcrit=1.65

Figure 20. Benchmark case for applied displacement, w, plotted versus increase in 

delamination length, a*, for a MMB specimen (80% mode II).

0.0 2.0 4.0 6.0 8.0 10.0

0

200

400

600

800

1000

critical, G

T

=Gc 

benchmark for applied load P

benchmark for applied displacement w

load P, N

displacement w, mm

Figure 19. Calculated critical behavior and resulting benchmark cases for 

applied load, P, and displacement, w, for a MMB specimen (80% mode II).

wcrit=1.65

Pcrit

=751

35
0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 2.2

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

2D FE-results

2D FE-results

2D FE-results

3D FE-results

3D FE-results

3D FE-results

applied displacement w, mm

Figure 22. Computed total energy release rate.

GT,

 kJ/m2

GII /GT=0.2

GII /GT

=0.8

GII /GT

=0.5

0 5 10 15 20 25 30 35 40 45 50 55

0

200

400

600

800

1000

critical, G

T

=Gc 

benchmark

applied

load 

P, N

increase in delamination length a*, mm

Figure 21. Benchmark case for applied load, P, plotted versus increase in 

delamination length, a*, for a MMB specimen (80% mode II).

Pcrit=751

36
0.00 0.20 0.40 0.60 0.80 1.00 1.20 1.40 1.60 1.80 2.00 2.20

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

1.1

2D FE-results

2D FE-results

2D FE-results

3D FE-results

3D FE-results

3D FE-results

applied displacement w, mm 

Figure 24. Computed mixed-mode ratio.

GII /G

T

GII /GT=0.2

GII /GT

=0.8

GII /GT

=0.5

0.00 0.20 0.40 0.60 0.80 1.00 1.20 1.40 1.60 1.80 2.00 2.20

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

2D FE-results

2D FE-results

2D FE-results

3D FE-results

3D FE-results

3D FE-results

applied displacement w, mm

Figure 23. Computed mode II energy release rate component.

GII ,

 kJ/m2

GII /GT=0.2

GII /GT

=0.8

GII /GT

=0.5

37
0.0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1.0

1.1

-0.5 -0.4 -0.3 -0.2 -0.1 0 0.1 0.2 0.3 0.4 0.5

2D FE-results 3D FE-results

2D FE-results 3D FE-results

2D FE-results 3D FE-results

y/B

Figure 26. Computed total energy release rate distribution across the width of a MMB specimen.

GII /GT=0.2

GII /GT

=0.8

GII /GT

=0.5

GT,

 kJ/m2

0 5 10 15 20 25 30

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

1.1

2D FE-results

2D FE-results

2D FE-results

increase in delamination length a* ,mm

Figure 25. Dependence of computed mixed-mode ratio on delamination length.

GII /G

T

location of loading point

GII /GT=0.2

GII /GT

=0.8

GII /GT

=0.5

38
0.00

0.10

0.20

0.30

0.40

0.50

0.60

0.70

0.80

0.90

1.00

1.10

-0.5 -0.4 -0.3 -0.2 -0.1 0 0.1 0.2 0.3 0.4 0.5

2D FE-results 3D FE-results

2D FE-results 3D FE-results

2D FE-results 3D FE-results

y/B

Figure 28. Computed mixed-mode ratio distribution across the width of a MMB specimen.

GII /GT

GII /GT=0.2

GII /GT

=0.8

GII /GT

=0.5

0.0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1.0

1.1

-0.5 -0.4 -0.3 -0.2 -0.1 0 0.1 0.2 0.3 0.4 0.5

2D FE-results 3D FE-results

2D FE-results 3D FE-results

2D FE-results 3D FE-results

y/B

Figure 27. Computed mode II energy release rate distribution 

across the width of a MMB specimen.

GII /GT=0.2

GII /GT

=0.8

GII /GT

=0.5

GII

 kJ/m2

39
0.0 2.0 4.0 6.0 8.0 10.0

0

20

40

60

80

100

120

140

160

critical, benchmark

reltol=0.1

reltol=0.01

load P, N

applied displacement w, mm

wcrit=1.64

Pcrit

=128.5 

Figure 30. Computed critical load-displacement behavior for a MMB specimen (20% mode II)

obtained from two-dimensional planar models (CPE4I) with different release tolerance settings.

0.0

0.5

1.0

1.5

2.0

2.5

-0.5 -0.4 -0.3 -0.2 -0.1 0 0.1 0.2 0.3 0.4 0.5

2D FE-results 3D FE-results

2D FE-results 3D FE-results

2D FE-results 3D FE-results

y/B

Figure 29. Failure index distribution across the width of a MMB specimen.

GII /GT=0.2

GII /GT

=0.8

GII /GT

=0.5

GT /Gc

40
0 5 10 15 20 25 30 35 40 45 50 55

0

2

4

6

8

10

critical, benchmark

reltol=0.2, cs=1E-6

reltol=0.1, cs=1E-6

reltol=0.01; cs=1E-6

applied

displacement

w, mm

increase in delamination length a*, mm

Figure 32. Computed displacement-propagation behavior for a MMB specimen (20% mode II) 

obtained from two-dimensional planar models (CPE4I) with added contact stabilization.

0.0 2.0 4.0 6.0 8.0 10.0

0

20

40

60

80

100

120

140

160

critical, benchmark

reltol=0.2, cs=1E-6

reltol=0.1, cs=1E-6

reltol=0.01; cs=1E-6

load P, N

applied displacement w, mm

wcrit=1.64

Pcrit

=128.5 

Figure 31. Computed critical load-displacement behavior for a MMB specimen (20% mode II)

obtained from two-dimensional planar models (CPE4I) with added contact stabilization.

41
0 5 10 15 20 25 30 35 40 45 50 55

0

2

4

6

8

10

critical, benchmark 

from 2D FE analysis

edge 1

edge 2

center

applied

displacement

w, mm

increase in delamination length a*, mm

G=Gc

wcrit=1.64

Figure 34. Computed displacement-propagation behavior for a MMB specimen (20% mode II) 

obtained from three-dimensional solid models (C3D8I) with added contact stabilization.

0.0 2.0 4.0 6.0 8.0 10.0

0

20

40

60

80

100

120

140

160

critical, benchmark 

from 2D FE analysis

reltol=0.1, cs=1E-6

load P, N

applied displacement w, mm

wcrit=1.64

Pcrit

=128.5 

Figure 33. Computed critical load-displacement behavior forMMB specimen (20% mode II)

obtained from three-dimensional solid models (C3D8I) with added contact stabilization.

42
0.0 2.0 4.0 6.0 8.0 10.0

0

100

200

300

400

500

600

critical, benchmark

reltol=0.01; cs=1E-6

load P, N

applied displacement w, mm

wcrit=1.34

Pcrit

=385

Figure 36. Computed critical load-displacement behavior for MMB specimen (50% mode II)

obtained from two-dimensional planar models (CPE4I) with added contact stabilization.

0.0 2.0 4.0 6.0 8.0 10.0

0

100

200

300

400

500

600

critical, benchmark

reltol=0.2

reltol=0.1

reltol=0.01

load P, N

applied displacement w, mm

wcrit=1.34

Pcrit

=385

Figure 35. Computed critical load-displacement behavior for MMB specimen (50% mode II)

obtained from two-dimensional planar models (CPE4I) with different release tolerance settings.

43
0.0 2.0 4.0 6.0 8.0 10.0

0

100

200

300

400

500

600

critical, benchmark

reltol=0.1; cs=1E-6

load P, N

applied displacement w, mm

wcrit=1.34

Pcrit

=385

Figure 38. Computed critical load-displacement behavior for MMB specimen (50% mode II)

obtained from three-dimensional solid models (C3D8I) with added contact stabilization.

0.0 10.0 20.0 30.0 40.0 50.0

0

1

2

3

4

5

6

7

8

critical

benchmark

reltol=0.2

reltol=0.1

reltol=0.01; cs=1E-6

applied

displacement

w, mm

increase in delamination length a*, mm

Figure 37. Computed critical displacement-propagation behavior for MMB specimen (50% mode II)

obtained from two-dimensional planar models.

44
0.0 2.0 4.0 6.0 8.0 10.0

0

100

200

300

400

500

600

critical 

benchmark

reltol=0.2; cs=1E-6

reltol=0.1; cs=1E-6

reltol=0.01; cs=1E-6

applied

load 

P, N

displacement w, mm

wcrit=1.34

Pcrit

=385

Figure 40. Computed critical load-displacement behavior for MMB specimen (50% mode II)

obtained from two-dimensional planar models (CPE4I) with added contact stabilization.

0.0 2.0 4.0 6.0 8.0 10.0

0

100

200

300

400

500

600

critical 

benchmark

reltol=0.2

applied

load 

P, N

displacement w, mm

wcrit=1.34

Pcrit

=385

Figure 39. Computed critical load-displacement behavior for MMB specimen (50% mode II)

obtained from two-dimensional planar models (CPE4I) with different release tolerance settings.

45
0.0 2.0 4.0 6.0 8.0 10.0

0

100

200

300

400

500

600

critical 

benchmark

reltol=0.1; cs=1E-6

applied

load 

P, N

displacement w, mm

wcrit=1.34

Pcrit

=385

Figure 42. Computed critical load-displacement behavior for MMB specimen (50% mode II)

obtained from three-dimensional solid models (C3D8I) with added contact stabilization.

0 5 10 15 20 25 30 35 40 45 50 55

0

100

200

300

400

500

600

critical

benchmark

reltol=0.2; cs=1E-6

reltol=0.1; cs=1E-6

reltol=0.01; cs=1E-6

applied

load 

P, N

increase in delamination length a*, mm

Figure 41. Computed critical load-propagation behavior for MMB specimen (50% mode II) 

obtained from two-dimensional planar models.

46
0.0 2.0 4.0 6.0 8.0 10.0

0

200

400

600

800

1000

critical

benchmark

reltol=0.2; cs=1E-6

reltol=0.1; cs=1E-6

reltol=0.01; cs=1E-6

load P, N

applied displacement w, mm

Figure 44. Computed critical load-displacement behavior for MMB specimen (80% mode II)

obtained from two-dimensional planar models (CPE4I) with added contact stabilization.

wcrit=1.65

Pcrit

=751

0.0 2.0 4.0 6.0 8.0 10.0

0

200

400

600

800

1000

critical

reltol=0.5

reltol=0.2

benchmark

load P, N

applied displacement w, mm

Figure 43. Computed critical load-displacement behavior for MMB specimen (80% mode II)

obtained from two-dimensional planar models (CPE4I) with different release tolerance settings.

wcrit=1.65

Pcrit

=751

47
0.0 2.0 4.0 6.0 8.0 10.0

0

200

400

600

800

1000

critical

benchmark

reltol=0.1 cs=1E-6

reltol=0.01; cs=1E-6

load P, N

applied displacement w, mm

Figure 46. Computed critical load-displacement behavior for MMB specimen (80% mode II)

obtained from three-dimensional solid models (C3D8I) with added contact stabilization.

wcrit=1.65

Pcrit

=751

0 5 10 15 20 25 30 35 40 45 50 55

0

1

2

3

4

5

6

7

8

critical

benchmark

reltol=0.2; cs=1E-6

reltol=0.1; cs=1E-6

reltol=0.01; cs=1E-6

applied

displacement

w, mm

increase in delamination length a*, mm

Figure 45. Computed critical displacement-propagation behavior for MMB specimen (80% mode II)

obtained from two-dimensional planar models.

48
0.0 2.0 4.0 6.0 8.0 10.0

0

200

400

600

800

1000

critical

benchmark

reltol=0.2; cs=1E-6

reltol=0.1; cs=1E-6

reltol=0.01; cs=1E-6

applied 

load

P, N

displacement w, mm

Figure 48. Computed critical load-displacement behavior for MMB specimen (80% mode II)

obtained from two-dimensional planar models (CPE4I) with added contact stabilization.

wcrit=1.65

Pcrit

=751

0.0 2.0 4.0 6.0 8.0 10.0

0

200

400

600

800

1000

critical

reltol=0.2

benchmark

applied 

load 

P, N

displacement w, mm

Figure 47. Computed critical load-displacement behavior for MMB specimen (80% mode II)

obtained from two-dimensional planar models (CPE4I) with different release tolerance settings.

wcrit=1.65

Pcrit

=751

49
0.0 2.0 4.0 6.0 8.0 10.0

0

200

400

600

800

1000

critical from 2D FE analysis

benchmark from 2D FE analysis

reltol=0.1 cs=1E-6

reltol=0.01 cs=1E-6

applied

load 

P, N

displacement w, mm

Figure 50. Computed critical load-displacement behavior for MMB specimen (80% mode II)

obtained from three-dimensional solid models (C3D8I) with added contact stabilization.

wcrit=1.65

Pcrit

=751

0 5 10 15 20 25 30 35 40 45 50 55

0

200

400

600

800

1000

critical

benchmark

reltol=0.2; cs=1E-6

reltol=0.1; cs=1E-6

reltol=0.01; cs=1E-6

applied

load 

P, N

increase in delamination length a*, mm

Figure 49. Computed critical load-propagation behavior for MMB specimen (80% mode II)

obtained from two-dimensional planar models.

50
0 5 10 15 20 25 30 35 40 45 50 55

0

200

400

600

800

1000

critical from 2D FE analysis

benchmark from 2D FE analysis

reltol=0.1 cs=1E-6, edge 1

center

edge 2

reltol=0.01 cs=1E-6, edge 1

edge 2

center

applied

load 

P, N

increase in delamination length a*, mm

Figure 51. Computed critical load-propagation behavior for MMB specimen (80% mode II) 

obtained from three-dimensional solid models (C3D8I) with added contact stabilization.

51
PTFE film

Figure 52. Fracture surfaces of typical MMB specimens [7].

a. 20 % mode II.

B

B

PTFE film

PTFE film

B

initial front precrack final front

initial front precrack final front

initial front fiber bridging final front

b. 50% mode II.

c. 80% mode II.

10 mm

10 mm

10 mm

52
Figure 53. Computed delamination front shape for MMB specimen (20% mode II).

location of initial delamination front

delaminated

 section

intact section

b. Detail of rst delamination propagation. c. Detail of delamination propagation. 

edge 1

edge 2

delamination section not contributing to analysis

d. Detail of delamination propagation. e. Detail of delamination propagation. 

location of initial 

delamination front

location of initial 

delamination front

a. Initial delamination front shape (Bottom surface of FE model in Figure 5).

53
Figure 55. Computed delamination front for MMB specimen (50% GII - applied load).

a. Detail of rst delamination propagation. b. Detail of early propagation. 

c. Detail of advanced propagation. d. Detail of nal front.

location of initial 

delamination front

location of initial 

delamination front

Figure 54. Computed delamination front for MMB specimen (50% GII - applied displacement).

a. Detail of rst delamination propagation. b. Detail of early propagation. 

c. Detail of advanced propagation. d. Detail of nal front. 

location of initial 

delamination front

location of initial 

delamination front

initial front

initial front

54
Figure 56. Computed delamination front shape for MMB specimen (80% mode II).

a. Initial delamination front shape (Bottom surface of FE model in Figure 5).

location of initial delamination front

delaminated

 section

intact section

b. Detail of rst initial delamination front. c. Detail of rst delamination propagation. 

edge 1

edge 2

delamination section not contributing to analysis

d. Detail of delamination propagation. e. Detail of delamination propagation. 

location of initial 

delamination front

location of initial 

delamination front

55
initial delamination front

w

w=0

u,w=0

x,u

z,w

rigid elements to 

simulate load 

apparatus

(ELSET=MMB_RIG)

rigid elements

(R2D2)

*RIGID BODY, ELSET=MMB_RIG, REF NODE=MMBREF, PIN NSET=NPIN

 bottom surface

(VCCT_BOT)

top surface

(VCCT_TOP)

bonded nodes

(BONDED)

node set

(MMBREF)

rigid elements to 

simulate load 

apparatus

(ELSET=MMB_RIG)

(a). View of full model. 

(b). Detail of specimen tip and crack tip zone. 

Figure A1. Deformed 2D FE-model of a MMB specimen.

(c). Detail of load introduction zone.

x,u

z,w

bonded nodes

(BONDED)

rigid elements

(R2D2) *EQUATION

2

616,1,1.0,100003,1,-1.0

node set(NPIN)

56
0

0.2

0.4

0.6

0.8

1

1.2

1.4

1.6

0 0.2 0.4 0.6 0.8 1

experimental data [7,10,11]

averaged data and curve fit

DCB, Mode I MMB, Mode I and II

Gc ,

kJ/m2

mixed mode ratio GII/GT

ENF, Mode II

G

Ic

 =0.212

 GIIc = 0.774

Figure A2. Mixed mode fracture criterion for IM7/8552.

*FRACTURE CRITERION, TYPE=fatigue, MIXED MODE BEHAVIOR=BK, 

TOLERANCE=<reltol>

<c1>,<c2>,<c3>,<c4>,<r1>,<r2>,<GIc>,<GIIc>,

<GIIIc>,<eta>

curve fit: Gc = G

Ic

 + (GIIc-GIc)(GII/GT)

= 2.1

57
REPORT DOCUMENTATION PAGE Form Approved

OMB No. 0704-0188

2. REPORT TYPE 

Contractor Report

 4. TITLE AND SUBTITLE

Development and Application of Benchmark Examples for Mixed-Mode 

I/II Quasi-Static Delamination Propagation Predictions

5a. CONTRACT NUMBER

 6. AUTHOR(S)

Krueger, Ronald

 7. PERFORMING ORGANIZATION NAME(S) AND ADDRESS(ES)

NASA Langley Research Center

Hampton, VA 23681-2199

 9. SPONSORING/MONITORING AGENCY NAME(S) AND ADDRESS(ES)

National Aeronautics and Space Administration

Washington, DC 20546-0001

 8. PERFORMING ORGANIZATION

 REPORT NUMBER

NIA Report No. 2012-01

10. SPONSOR/MONITOR'S ACRONYM(S)

NASA

13. SUPPLEMENTARY NOTES

Langley Technical Monitor: Jonathan B. Ransom

12. DISTRIBUTION/AVAILABILITY STATEMENT

Unclassified - Unlimited

Subject Category 24

Availability: NASA CASI (443) 757-5802

19a. NAME OF RESPONSIBLE PERSON

STI Help Desk (email: help@sti.nasa.gov)

14. ABSTRACT

The development of benchmark examples for quasi-static delamination propagation prediction is presented and demonstrated for a commercial code. The 

examples are based on finite element models of the Mixed-Mode Bending (MMB) specimen. The examples are independent of the analysis software used and 

allow the assessment of the automated delamination propagation prediction capability in commercial finite element codes based on the virtual crack closure 

technique (VCCT). First, quasi-static benchmark examples were created for the specimen. Second, starting from an initially straight front, the delamination was 

allowed to propagate under quasi-static loading. Third, the load-displacement relationship from a propagation analysis and the benchmark results were 

compared, and good agreement could be achieved by selecting the appropriate input parameters. Good agreement between the results obtained from the 

automated propagation analysis and the benchmark results could be achieved by selecting input parameters that had previously been determined during 

analyses of mode I Double Cantilever Beam and mode II End Notched Flexure specimens. The benchmarking procedure proved valuable by highlighting the 

issues associated with choosing the input parameters of the particular implementation. Overall the results are encouraging, but further assessment for 

mixed-mode delamination fatigue onset and growth is required.

15. SUBJECT TERMS

benchmarking; composites; delamination; finite element analysis; fracture mechanics; virtual crack closure technique

18. NUMBER

 OF 

 PAGES

62

19b. TELEPHONE NUMBER (Include area code)

(443) 757-5802

a. REPORT

U

c. THIS PAGE

U

b. ABSTRACT

U

17. LIMITATION OF 

 ABSTRACT

UU

Prescribed by ANSI Std. Z39.18

Standard Form 298 (Rev. 8-98)

3. DATES COVERED (From - To)

5b. GRANT NUMBER

NNL09AA00A

5c. PROGRAM ELEMENT NUMBER

5d. PROJECT NUMBER

5e. TASK NUMBER

5f. WORK UNIT NUMBER

877868.02.07.07.05.01.01

11. SPONSOR/MONITOR'S REPORT

 NUMBER(S)

NASA/CR-2012-217562

16. SECURITY CLASSIFICATION OF:

The public reporting burden for this collection of information is estimated to average 1 hour per response, including the time for reviewing instructions, searching existing data sources, 

gathering and maintaining the data needed, and completing and reviewing the collection of information. Send comments regarding this burden estimate or any other aspect of this 

collection of information, including suggestions for reducing this burden, to Department of Defense, Washington Headquarters Services, Directorate for Information Operations and 

Reports (0704-0188), 1215 Jefferson Davis Highway, Suite 1204, Arlington, VA 22202-4302. Respondents should be aware that notwithstanding any other provision of law, no person 

shall be subject to any penalty for failing to comply with a collection of information if it does not display a currently valid OMB control number.

PLEASE DO NOT RETURN YOUR FORM TO THE ABOVE ADDRESS.

1. REPORT DATE (DD-MM-YYYY)

01- 04 - 2012