Stability of Cohesive Crack

Model: Part II--Eigenvalue

Z. P. Ba2ant
Walter Murphy Professor of Civil
Analysis of Size Effect on
Engineering and Material Science,
Fellow ASME, Strength and Ductility of
Yuan-Neng Li
Research Assistant Professor.
Structures
The preceding paper is extended to the analysis of size effect on strength and ductility
Department ol Civil Engineering, of structures. For the case of geometrically similar structures of different sizes, the
Northwestern University, criterion of stability limit is transformed to an eigenvalue problem for a homogeneous
Evanston, IL 60208
Fredholm integral equation, with the structure size as the eigenvalue. Under the
assumption of a linear softening stress-displacement relation for the cohesive crack,
the eigenvalue problem is linear. The maximum load o f structure under load control,
as well as the maximum deflection under displacement control (which characterizes
ductility of the structure), can be solved explicitly in terms of the eigenfunction of
the aforementioned integral equation.

1 Introduction influence functions are used to formulate the condition of stabil-

ity limit of a structure with a cohesive crack in the form of a
As explained in the preceding paper (Ba~ant and Li, 1995), homogenous Fredholm integral equation. The peak load, as well
the cohesive crack model is a nonlinear theory of fracture me- as the maximum displacement (which corresponds to snap-back
chanics in which the condition of stability limit is expressed in instability), is obtained. In addition, the cases of a structure
terms of the singularity condition of the second variation of
loaded through a spring coupled in series (i.e., the case o f a soft
the energy potential with respect to cohesive stresses or crack- loading device) and a structure restrained by a spring coupled in
opening displacements. Although the criterion of stability limit parallel are analyzed. Finally, some computational techniques
can also be formulated in terms of energy variation with respect are discussed and a numerical example of the size effect curves
to the crack length, the resulting equation is not very useful, for maximum deflection is given.
since the energy release rate in the cohesive crack model de-
pends on the cohesive stresses or crack-opening displacements.
For a given structure, the criterion of stability limit leads to 2 Dimensionless Process Zone Equations
a highly nonlinear equation for crack length. However, when a We consider a two-dimensional structure of a unit thickness
class of geometrically similar structures of different sizes is and introduce the following dimensionless variables:
considered and the relative crack length is given, the criterion
of stability limit can be treated as an equation for the structure a x o-
size at which the stability limit occurs at the given relative crack
length. In this manner, the criterion of the stability limit is
transformed into an eigenvalue problem, with the structure size w P D
w=--, P=-- O= (1)
as the eigenvalue. In the special case of linear softening, the wc Dft ' 2L0
eigenvalue problem is linear. It can be solved independently
of the solution of the cohesive crack model. Furthermore, the where L0 = E G j / f 2 = characteristic size of the process zone,
corresponding maximum value of the load or loading parameter ft = direct tensile strength of the material, and wc = threshold
can be expressed explicitly in terms of the eigenfunction. In value of the crack-opening displacement. All the notations from
this way, the size effect curve can be obtained readily, without the preceding paper (Ba~ant and Li, 1994) are retained. To
having to calculate the load-deflection curves for structures of simplify notations in the following text, we will drop the bars,
various sizes. with the understanding that all the variables are dimensionless
The eigenvalue problem of the cohesive crack model was unless specified otherwise.
studied by Li and Hong (1992), Li and Liang (1993) and Li For a generic elastic structure, the crack-opening displace-
and Ba~ant (1993). However, only the peak-load solution was ment w, the load-line displacement u, the load P, and the crack-
discussed in these previous papers. In the present paper, the bridging stress a must satisfy the compatibility equations:

w(() _ C~'°(~, ~')o-(~')d~' + C'~e(~)P (2)

Contributed by the Applied Mechanics Division of THE AMERICANSOCIETY OF D o
MECHANICAL ENGINEERS for presentation at the International Mechanical Engi-
neering Congress and Exhibition, the Winter Annual Meeting of the ASME, San
Francisco, CA, Nov. 12:17, 1995. u _ Ce°(~')cr(~')d~ ' + CPPP (3)
Discussion on this paper should be addressed to the Technical Editor, Professor D o
Lewis T. Wheeler, Department of Mechanical Engineering, University of Houston,
Houston, TX 77204-4792, and will be accepted until four months after final which represent the special case of Eqs. ( 11 ) and (10) or (14)
publication of the paper itself in the ASME JOURNALOF APPLIED MECHANICS. of the preceding paper for CI = 0; C°°(~, ~'), C°P(~), Ce°(~),
Manuscript received by the ASME Applied Mechanics Division, Mar. 21, 1995;
final revision, Apr. 29, 1995. Associate Technical Editor: G. J. Dvorak. C ee are dimensionless compliance influence functions (Green's
Paper No. 95-WA/APM-8. functions). The lower integration limit a0 is the relative length

Journal of Applied Mechanics DECEMBER 1995, Vol. 62 / 965

Copyright © 1995 by ASME
Downloaded From: http://pressurevesseltech.asmedigitalcollection.asme.org/ on 05/25/2015 Terms of Use: http://asme.org/terms
 of the initial traction-free crack (notch); a is the total relative v ( ( ) , with size D as the eigenvalue. The size D for which
crack length which includes both the process zone (crack-bridg- the given a corresponds to the maximum load is the largest
ing zone) and the stress-free crack. The problem can also be eigenvalue of ( 11 ). This approach, proposed by Li and Ba~ant
formulated as equilibrium conditions written in terms of stiff- (1994), makes it possible to avoid solving the load-deflection
ness influence functions: curves for various sizes D. It represents an efficient method of
calculating the size effect curve.
De(G) = RWW(~,~ ' ) w ( ~ ' ) d ~ ' + RW"(~)u (4) The dimensionless crack compatibility equation can be writ-
Ot0 ten as

DP = f2
0
RUW(~)w(~)d~ + R""u. (5) 1 - ~(~) = _
D f2o
ca~(6 Clo.(C)dC + c~P(()P. (12)

These equations represent the special case of Eqs. (26) and Multiplying this with the eigenfunction v(~) and then inte-
(25). Equation (4) for prescribed load P ensues by solving u grating with respect to ~, we obtain
from Eq. (5) and substituting it into Eq. (4). The dimensionless
stiffness functions are here defined with a unit value of Young's
modulus.
In the cohesive crack model, the cohesive stress o. is related
0
- C~P(~)P v(~)d~ =
y: f:[
0
o.(~)
1
0
8(~ - ~')

to the crack-opening displacement w by the stress-displacement C~a(~, ~ ' ) ] v ( ~ ' ) d ~ ' d ~ . (13)
relation, which can be described by either of the following
forms If the singularity condition is satisfied, then the applied load is
w = g(o'), o. = f ( w ) . (6a, b) at its maximum. This maximum value is found to be

Substituting (6a) into (2), we obtain what we call the crack

1 o v(~)d~
compatibility equation in terms of compliance functions:
P = -D- . (14)
c°e(~)v(~)d~
g[a(~)] = - C~(~, ~ ' ) a ( ~ ' ) d ~ ' + C~P(~)P (7) d~0
0

An equivalent expression for the peak load was obtained by Li

Substituting (6b) into (4), we obtain the crack equilibrium
and Hong (1992), and by Li and Ba~ant (1994). The eigen-
equation in terms of stiffness functions:
value problem (11 ) and the peak load solution (14) provide a
powerful set of equations for solving the size-effect curve of
Df[w(~)] = - RW~((, ~ ' ) w ( ( ' ) d ~ ' + RW"(()u. the cohesive crack model directly, without any need to solve
0
the load-deflection curve from the basic equations.
(8) The solution can also be generalized to include the case of
multiple (conservative) loads. They can vary arbitrarily but in
3 Peak-Load Solution by the Condition of Structural such a manner that there is no crack closure. Then the relation
among the load values at the stability limit of the structure is
Stability Limit linear. For instance, when a beam is subjected to combined
As established in Ba~ant and Li (1995), the singularity condi- action of lateral load P and axial load N, as shown in Fig. 1,
tion for the compliance formulation under load control can be the crack compatibility equation can be written as
expressed as the condition of finding a nonzero solution v(~)
of the following homogenous equation: 1 - o.(~)
D f:0
c~(~, C)~(C)d('

D
f:0
C~(C, Ov(C)dC = ag[~(~)-----!v ( O .
do"
(91 + C~P(()P + c ° N ( O N (15)

Since we are considering geometrically similar structures where the symbols are self-explanatory. Since the loading terms
only, (9) can be regarded as an eigenvalue problem if the rela- do not enter the criterion of stability limit, the equation for the
tive crack length ~ is given. The dimensionless quantity D structural stability limit remains the same. If the condition for
plays the role of an eigenvalue. In the actual calculation, the the stability limit is satisfied, the relation between these two
singularity condition should be solved simultaneously with the loads is found to be linear:
basic equations to obtain the nominal strength as the maximum P N
load parameter and the corresponding size for a given relative -- + - - = 1 (16)
crack length. Calculation of size effect curves in this manner P* N*
is very efficient. A discussion of the discrete form of the present where the denominators, defined as
formulation has been given by Li and Ba~ant (1994).
In the following, we restrict attention to the case of linear
softening, which is defined as
p*= 1 f; o
v(~)d~

w = g(o.) = 1 - or, a = f ( w ) = 1 - w. (10) D

flo C~P(~)v(~ld~
Since for linear softening dg/da = - 1 , the eigenvalue is no
longer coupled with the basic equations of the cohesive crack.
The eigenvalue problem can now be written as
N* = -- ~
1
,
f: o v(~)d~
(17)
D
f:
0
C'~a(~', ~ ) v ( ~ ' ) d ~ ' = v(~). (11) D f l C~N(~)v(~)d~
0

If the relative crack length is specified and geometrically similar represent the critical loads when P and N are applied to the
structures are considered, Eq. ( 11 ) represents a linear homoge- structure separately. Equation (16) is the general interaction
neous Fredholm integral equation (Tricomi, 1957) for function relation when the structure fails by tensile fracture and the

966 / Vol. 62, DECEMBER 1995 Transactions of the ASME

Downloaded From: http://pressurevesseltech.asmedigitalcollection.asme.org/ on 05/25/2015 Terms of Use: http://asme.org/terms

 softening stress-displacement law is linear. A relation of this 5 Stability Limit of Structure Loaded Through a
type was also reported by Li, Mtfller, and W6rner (1994) in a Spring
discrete (matrix) form. Generalization to an arbitrary number
of applied loads is self-evident. If the device that controls loading (e.g., the testing machine)
When the stress-displacement relation for a cohesive crack has finite compliance C/, the device can be represented as a
is nonlinear, one can use an iterative succession of linear ap- spring connected to the structure in series. In such a connection,
proximations representing tangents of the stress-displacement the device and the structure share the same force. Denote as u
curve according to the preceding approximation (this approach the total deflection that is controlled, which is the sum of the
was formulated for the maximum load in Li and Ba~ant, 1994). deflection us of the structure and the deflection of the device u
- u~. Using (3), we can solve load P in terms of u as
4 Solution of Maximum Deflection P = (C ~e + C:) -1 - CP°(()a(~)d( (25)
If the structure is loaded by controlled displacement (i.e., 0

with a rigid grip), the stability limit is reached when there is a In the dimensionless form, the process zone equation is
snap back in the diagram of load P versus load-line displace-
ment u. The crack equilibrium equation for this case is Eq. (30)
of the preceding paper which, in the case of linear softening,
1 - a(() = -D L 0
d~((, (')cr((')d('
yields
+ C ° e ( ( ) ( C ee + C:)-Ju (26)
[1 - w ( ~ ) l D = -
L 0
R~W(~, ~ ' ) w ( ~ ' ) d ~ ' + RWU(~)u. (18) where
d~°((, C)
The dimensionless condition of stability limit may now be writ-
ten as = C~°(~, ~') - C " e ( ~ ) ( C eP + C : ) - ' C e ' ( ~ ' ) . (27)

1
D
:i o
RWW(~, ~ ' ) v ( ~ ) d ~ = v ( ~ ' ) . (19)
The form of the eigenvalue problem is the same as (23) except
that the modified compliance function is defined by (27). The
maximum deflection is found to be
Since a is constant for geometrically similar structures, (19) is
a linear homogeneous Fredholm integral equation for the un- o
known cohesive stress v (~) in the process zone. This represents u = ( cp~ + G ) [,. (28)
an eigenvalue problem with 1/D as the eigenvalue. Only the C"ev(~)d~
smallest eigenvalue 1/D represents a stability limit. The maxi- L 0

mum deflection, characterizing snap back, is found to be This formula reduces to (24) when compliance C: approaches
zero.
D v(~)d~ On the other hand, if the spring is connected to the structure
u = ~, 0
(20) in parallel, it shares the same deflection with the structure.
SWUv(~)d~ Denote by P the total load applied to the structure-spring sys-
0
tem, which is the sum of the load PI which acts on the structure
and Szu where S: = 1/C:. Using (5) we can express u in terms
However, the maximum deflection can also be solved in terms of P as
of the compliance functions. To this end, we eliminate the load
parameter P from (3) and (7) and obtain the following crack u = R"W(~)w(~)d~ + DP ( R " " + R : ) - ' . (29)
compatibility equation under displacement control: 0

Substituting (29) into (8), one obtains the following crack

1 - cr(() = - D f' d ~ ( ~ , ~ ' ) a ( ~ ' ) d ~ ' + - C~P(~)
-~---u (21) compatibility equation:

fo £,~w(~,C)w(~')d~'
0

where
1 - w ( ~ ) = - ~ . ~10 (30)

C~°(~, ~') = C ~ ( ~ , ~') - C°P(~)CP'~(~ ') ~ . (22)

+ R w " ( ~ ) P ( R "" + Ry)-'
The corresponding eigenvalue problem now becomes where the modified stiffness function is defined as
/~w~(~,C)
v(~) = D da~(~ ', ~ ) v ( ~ ' ) d ~ ' . (23)
0 = Rw~(~, ~ ' ) - R ~ " ( ~ ) R " ~ ( ~ ' ) ( R "" + R f ) - ' . (31)
This is equivalent to the eigenvalue problem (19) of stiffness The eigenvalue problem is to find a nonzero eigenfunction
formulation, because the modified compliance function is the satisfying
inverse of the stiffness function R ~w. The maximum deflection
can be expressed as v(~) = ~1 f'o/~ww(~, C)v(C)d~'. (32)
The maximum total applied load can be calculated from the
u = C Pe o (24) following equation:
f ' c~'%(~)d~
0 f' v(~ld~
0
P=(R"+R:) ~ (33)
The compliance formulation is of course equivalent to the stiff-
ness formulation. In a similar way, we can also express the
0
maximum load in terms of the stiffness influence functions. The
details will not be given because the derivation is analogous. Of course when the spring constant of the connected spring

Journal of Applied Mechanics DECEMBER 1995, VoL 62 / 967

Downloaded From: http://pressurevesseltech.asmedigitalcollection.asme.org/ on 05/25/2015 Terms of Use: http://asme.org/terms

 approaches zero, (33) becomes the peak-load solution in the (ft/GfD)A (ft/G~D)A
stiffness formulation without a spring.

6 Numerical Implementation
As a numerical example, a three-point bent fracture specimen
(Fig. 1 ) is analyzed. The finite element method is used to obtain 15 /

the compliance functions in a discretized form (although other 10

10;
methods, such as the boundary element method, might also be 0
suitable). The four-node finite element, which is the simplest, 5' 5 42
is chosen to discretize one half of the beam. To determine
the nodal compliance matrix, the displacement solutions are ol
~ ~ 0 0.6 04 0.2 D/Lo
obtained for one unit load applied successively at each node
L ° a / l ~ a/D
along the potential crack path or at the load point.
Each column of matrix C ~" represents nodal displacements Fig. 2 Dimensionless deflection for (a) d0 = 0 and (b) do = 0.2
on the crack line when a unit load is applied to one node in the
process zone, C ~P represents the nodal displacements in the
process zone when a unit load is applied at the load point, and
C ee represents the load-line displacement when a unit load is mathematically similar, there exists one important difference.
applied at the load point. During the calculation, the total rela- For three-point-bent beams, the maximum load always exists,
tive crack length ce is first taken to correspond to the node that no matter how large the relative process zone length a - a0 is,
is farthest from the crack mouth as allowed by the compliance or how small the dimensionless size D / L o is. However, for
matrix, and then cracks reaching successively to nodes closer maximum deflections, the situation is different. As shown in
and closer to the crack mouth are considered. In each case Fig. 2, there is no maximum deflection if the relative length a
the nodal displacements that lie in the uncracked ligament are of the cohesive crack is large enough. The smallest dimen-
eliminated. In this way, the dependence of the compliance func- sionless size D below which there is no snap back will be called
tion on the crack length is reflected by the sizes of the compli- the critical size of the structure. The critical size is a function
ance matrices. of relative notch depth so as well as the slenderness ratio (span-
Starting with Hillerborg (1976), the zero-K condition has to-depth ratio of the beam).
been approximated by the condition that the elastic stress ahead The dependence of the critical size on the relative notch depth
of the cohesive crack tip be equal to the tensile strength. So in can also be seen in Fig. 2. Fig. 2(a) gives the deflection for beams
our dimensionless definition, o-tip = 1. In the space of continuous without a notch (so = 0), and Fig. 2(b) for beams with relative
functions, this condition is mathematically equivalent to the notch depth ce0 = 0.2. For s0 = 0, the critical size is found to be
condition that the stress intensity factor K at the crack tip be approximately 0.43 and for a0 = 0.2 approximately 1.4.
zero (Barenblatt, 1962). After discretization, however, these According to Eqs. (23) and (24), we can obtain the size
two conditions are equivalent only approximately. Thus the effect curves for maximum deflection for any given relative
use of the condition ~rtiv = 1 inevitably introduces additional length a of the cohesive crack. Figure 3 shows the size effect
numerical error into the discrete solution. But this small price curves for different initial notch ratios. Note that, paradoxically,
is quite justifiable, because we do not need the corresponding the curves extend even to the left of the critical sizes (dashed
stress intensity factors, which are not easy to calculate anyway. lines); these portions of the curves are of course physically
However, numerical results (Li and B a~ant, 1994) show that, meaningless since there exists no maximum deflection at all.
in order to obtain good accuracy for large (dimensionless) struc- The explanation is that these portions correspond to cases with
tural sizes, it seems important to assume the cohesive stress to negative cr, whereas our analytical expressions are valid only
vary linearly from node to node in the process zone, rather when the crack-opening displacement in the process zone is
than a piece-wise constant manner. The assumption of linear less than the crack-opening threshold wc (at which the stress is
variation of cohesive stress between the nodes leads to a tridiag- reduced to zero). With careful observation, one finds that, when
onal matrix connecting the nodal values of cohesive stresses to the condition of stability limit is satisfied, the critical size D is
the cohesive nodal forces (in detail, see Li and Ba~ant, 1994). actually the size at which the crack-mouth-opening displace-
Numerically, the differences in the maximum load values calcu- ment becomes equal to the threshold we. Above the critical size
lated by the eigenvalue analysis and by the load-deflection (i.e., on right portions of the curves in Fig. 3), the obtained
curves are usually in the fifth or sixth digit for linear (or nearly maximum deflections are exactly what one would obtain if the
linear) softening laws. load-deflection curve were solved by the conventional method,
that is, by solving the basic equations step by step for each
different cohesive crack lengths.
7 Numerical Solution of the Maximum Deflection
Numerical examples for the peak load solution using the
eigenvalue approach have been given in previous papers (e.g.,
t
Li and Hong, 1992; Li and Ba~ant, 1994). Therefore, we will
discuss only the numerical solution of maximum deflection, 10
which characterizes ductility of a structure. Although the maxi-
mum load solution and the maximum deflection solution are

c~0 = 0 . 2

A ? ,
0.5
i
1
D/Lo
,
5
p
10

Fig. 1 Beam under combined lateral load and axial load Fig. 3 Size effect curve for maximum deflection

9 6 8 / Vol. 6 2 , DECEMBER 1995 Transactions of the ASME

Downloaded From: http://pressurevesseltech.asmedigitalcollection.asme.org/ on 05/25/2015 Terms of Use: http://asme.org/terms

 As a check, we select, from the size effect curve, a maximum line deflection), can be determined through the eigenfunction
deflection value with its corresponding relative crack length obtained.
and its dimensionless size D. Then we use this dimensionless 3 Numerical examples of the solution of the maximum de-
size as the input and solve the process zone equation together flection as a function of the dimensionless beam depth demon-
with the crack-tip equation (~rtip = 1) for different crack lengths. strate that the maximum deflection solution ceases to be valid
In all the cases examined, the maximum deflection is found to if the structure dimension (e.g., beam depth) becomes smaller
be the same (within the numerical precision of the calculation) than a certain critical value. This critical value is characterized
and to occur at the same relative crack length. by the condition that the crack opening at the stability limit
reaches the threshold value at which the cohesive stress gets
reduced to zero.
8 Final Remarks and Conclusions
Acknowledgment
The cohesive crack model can be effectively analyzed in
Partial financial support for the present theory (in both Parts
terms of continuous influence functions. Under the assumption
I and II of the present study) was obtained under ONR grant
of a linear softening stress-displacement law, the criterion of
No. N00014-91-J-1109 to Northwestern University (monitored
stability limit, which has been analyzed by Ba~ant and Li
by Dr. Yapa Rajapakse). Additional financial support for some
(1995), becomes a linear eigenvalue problem when geometri-
ramifications to concrete fracture was obtained from the Center
cally similar structures are considered. The peak value of the
for Advanced Cement Based Materials at Northwestern Univer-
load parameter can be determined by solving the eigenvalue
sity.
problem. In this manner, the size effect of the cohesive crack
model becomes intimately related to the solutions of the eigen- References
value problem. There are some similarities between the eigen- Barenblatt, G. I., 1962, "The mathematical theory of equilibrium cracks in
value problem studied here and the eigenvalue problem for brittle fracture," Advances in Applied Mechanics, Vol. 7, pp. 55-129.
Ba~ant, Z. P., and Cedolin, L., 1991, Stability of Structures: Elastic, Inelastic,
the buckling load of a structure. Both eigenvalue problems are Fracture and Damage Theories, Oxford University Press, Oxford, U.K.
derived from the criterion of structural stability limit. Whereas, Ba~ant, Z. P., and Li, Y. N., 1995, "Stability of Cohesive Crack Model:
in the buckling problem, the eigenvalue is Euler's critical load, I. Energy Principles," ASME JOURNAL OF APPLIED MECHANICS, submitted for
in the cohesive crack model the eigenvalue is the structure size publication.
Hillerborg, A., Mod6er, M., and Petersson, P.-E., 1976, "Analysis of crack
for which the loading parameter is maximized at a given relative formation and crack growth in concrete by means of fracture mechanics and finite
cohesive crack length. The maximum load or load parameter elements," Cement and Concrete Research, Vol. 6, pp. 773-782.
can be calculated from the eigenfunctions. The following con- Li, Y. N., and Hong, A. P., 1992, ' 'Prediction of failure in notched infinite strips
clusions can be drawn: by the cohesive crack model," International Journal of Solids and Structures, Vol.
29, No. 23, pp. 2815-2828.
Li, Y. N., and Liang, R. Y., 1993, "The theory of boundary eigenvalue problem
1 When geometrically similar structures are considered, the in the cohesive crack model and its application," Journal of the Mechanics and
criterion of stability limit becomes an eigenvalue problem. The Physics of Solids, Vol. 41, No. 2, pp. 331-350.
size for which a given relative crack length corresponds to either Li, Y. N., and Ba~ant, Z. P., 1994, "Eigenvalue analysis of size effect in
cohesive crack model," International Journal of Fracture, Vol. 66, pp. 213-226.
the maximum load or the maximum displacement is the first Li, Y. N., MUller, M., and W~rner, J.-D., 1994, "Two new approaches in the
eigenvalue of a homogeneous Fredholm integral equation. The fictitious crack model and their applications," Journal of Engineering Mechanics
size effect curve can thus be calculated efficiently. ASCE, Vol. 120, No. 12, pp. 2604-2620.
2 If the softening stress-displacement law for the cohesive Petersson, P.-E., 1981, "Crack growth and development of fracture zones in
plain concrete and similar materials," Doctoral Dissertation, Lund Institute of
crack is linear, the eigenvalue problem becomes linearized and Technology, Sweden.
can be solved independently. The critical value of the loading Tricomi, F. G., 1957, Integral equations, Interscience Publishers, New York,
parameter (either the maximum load or the maximum load- (also Dover, 1985).

Journal of Applied Mechanics DECEMBER 1995, Vol. 62 / 969

Downloaded From: http://pressurevesseltech.asmedigitalcollection.asme.org/ on 05/25/2015 Terms of Use: http://asme.org/terms