Computers and Structures 106–107 (2012) 81–90

Contents lists available at SciVerse ScienceDirect

Computers and Structures

journal homepage: www.elsevier.com/locate/compstruc

Simulating bond failure in reinforced concrete by a plasticity model

Daiane de S. Brisotto, Eduardo Bittencourt ⇑, Virginia M. R. d’A. Bessa
Departamento de Engenharia Civil, Universidade Federal do Rio Grande do Sul, Av. Osvaldo Aranha, 99, 90035-190 Porto Alegre, RS, Brazil

a r t i c l e i n f o a b s t r a c t

Article history: In this work, an axisymmetric plasticity model is used to simulate the concrete-steel interface behavior. A
Received 25 November 2011 nonlocal correction is here introduced in order to capture the degradation of the bond due to splitting
Accepted 23 April 2012 cracks. Damage of the interface is also modeled as a function of the rib spacing, allowing application
Available online 19 May 2012
of the model to different bar diameters. The model is able to capture the transition from splitting to
pull-out failure and to yielding of bars with the same set of predefined interface parameters, showing
Keywords: the predictive character of the model. The development of macroscopic cracks is also correctly simulated.
Reinforced concrete
Ó 2012 Elsevier Ltd. All rights reserved.
Finite element
Fracture
Plasticity
Bond-slip

1. Introduction the concrete between them are still not explicitly considered and
are modeled by interface elements. Only in the high-resolution
In reinforced concrete structures, the bond between concrete scale ribs are explicitly considered [12], but such scale has not been
and steel bars is no longer perfect when concrete damage around fully explored yet due to complexities related to availability of data
bars takes place, changing the structural behavior. Two distinct and computational costs involved.
damage processes can be identified at the concrete-steel interface Although simpler, the two lower resolution scales depend heav-
[1]. The first occurs when concrete in the vicinity of the bars is not ily on bond-slip constitutive equations. These equations are, in
sufficiently confined by pressure or transverse reinforcement, general, based on parameters obtained by curve fitting with exper-
causing a splitting failure. Basically this process is triggered by iments. We can again divide constitutive equations in two groups.
cracks emerging from the bar, along its length. These cracks even- In the first, bond-slip relation is explicitly predefined [8,13–16]. In
tually propagate outward, reducing mechanical interlock due to rib the second, a yield or failure surface is predefined at the interface,
bearing. Splitting can be also caused by shrinking of the bar diam- being the bond-slip relation an outcome of the associated elastic–
eter due to steel yielding, in cases of long anchorage conditions. plastic solution [3,9,17,18]. In these cases, in principle, the transi-
The second type of damage occurs when a good level of confine- tion between the failure modes can be captured automatically be-
ment is provided and is caused by crushing of concrete between cause splitting stresses are fully coupled with bond stresses and,
ribs. The damage process is completed when the concrete between for this reason, these theories tend to be more predictive. Despite
ribs is sheared off. This failure process is called pull-out. Bond this apparent advantage, most of the applications follow the first
stresses in this case are much higher and failure occurs at much group. In this work, we decided to use a plasticity model in order
larger level of slip than in the case of splitting. For this reason to explore its weak and strong points in a variety of applications.
the failure by pull-out is considered less brittle than by splitting The yield surfaces used in the plasticity model proposed by
[2]. Lundgren and Gylltoft [9] are defined explicitly for the two failure
Considering finite element simulations, in general the bond be- modes. In the present work, we follow a similar model due to this
tween concrete and steel is modeled at three different scales [3]. In unique feature.
the low-resolution scale, steel is considered as unidimensional fi- Two main modifications in the Lundgren–Gylltoft theory [9] are
nite elements [4–7]. In the mid-resolution scale, both concrete proposed in the present work. First, concrete rupture is considered
and steel are modeled by volumetric finite elements [8–11]. Dowel by a discrete model. Such model can provide a more precise repre-
and Poisson’s effects can be automatically captured in this case, as sentation of macroscopic cracks than smeared schemes, introduc-
well as the splitting by steel yielding. Ribs or lugs of the bars and ing details such as the crack opening. For instance, it is shown
that the crack opening may act upon the bond behavior [19,20]
⇑ Corresponding author. Tel.: +55 51 3308 4268; fax: +55 51 3308 3999. and may be an important variable when bar corrosion is simulated
E-mail address: eduardo.bittencourt@ufrgs.br (E. Bittencourt). [21]. The relevance of the use of a discrete fracture methodology

0045-7949/$ - see front matter Ó 2012 Elsevier Ltd. All rights reserved.
http://dx.doi.org/10.1016/j.compstruc.2012.04.009
82 D.S. Brisotto et al. / Computers and Structures 106–107 (2012) 81–90

for reinforced concrete simulations is also acknowledged by the considered. The equation for the elastic matrix (Dij) is defined as
ACI [22]. The great disadvantage of discrete fracture methodologies follows
is the impossibility to represent longitudinal (splitting) cracks "   #
K 11 sgn uet K 12
when simplified 2D applications are made, such as axisymmetrical Dij ¼ Ec : ð4Þ
representations. In the present work, in order to overcome this dif- 0 K 22
ficulty, the effect of these cracks are incorporated in the interface The off-diagonal term in Eq. (4) expresses the interaction caused by
 
model. To accomplish that, bond stresses are considered depen- inclined ribs. The sgn uet factor indicates that slip produces sym-
dent also on the concrete conditions at the vicinity of the bar, metric negative values of tn, or interface compression. Cox and Herr-
turning the interface formulation nonlocal. Second, the Lund- mann [3] also consider this term. Ec is the elastic modulus of the
gren–Gylltoft model is extrapolated to be applied to bar diameters concrete and K11, K12 and K22 are elastic constants of the interface.
different from the original diameter used (16 mm). This is possible Physically Ec K22 represents the tangential traction, per unit of tan-
because a new damage parameter is introduced, based on rib gential deformation, a rib withstands elastically. For the specific
spacing. In addition, parameters used in the model are reviewed case of a 16 mm bar diameter, experiments reported by Lundgren
and, when available, compared with values used by other and Gylltoft [9] suggest that K22 ffi 0.006 mm1. This value is also
methodologies. coincident with the value used by Keuser and Mehlhorn [8] in a lin-
All applications presented here are restricted to loading applied ear-elastic analysis.
monotonically, axisymmetrical geometries and concretes with For uen negative, EcK11 can be seen as a penalty factor, being an
normal strength (compressive strength around 30 MPa). Section 2 arbitrary large positive value to avoid interface interpenetration.
describes the Lundgren–Gylltoft model and the changes proposed In this case, Lundgren [24] uses K11 equal to 0.15 mm1. For uen po-
in it by the present work. In Section 3 is briefly described a discrete sitive, or interface opening, K11 has a residual value, indicating a
fracture model [23] to represent concrete fracture. Applications to cohesion related to the structure strength at the vicinity of the
different cases are presented and discussed in Section 4. In the first bar. For the specific case of a 16 mm bar diameter, Lundgren [24]
two applications (Sections 4.1.1 and 4.1.2), short anchorage condi- considers that K11 drops linearly to nearly 0.005 mm1 until
tions are considered with different levels of confinement, while in uen ¼ 0:06 mm and then is kept constant for larger openings. Finally
the third case (Section 4.1.3) long anchorage is regarded. In the it is established that the coupling parameter K12 has a limit value
fourth application (Section 4.2) the development of primary cracks defined as K22/l, where l is the friction coefficient of the interface.
in the concrete is taken into account. Final remarks are considered In practice, the value K12 = 0.9K22/l is adopted [9,24]. In Cox and
in Section 5. Herrmann [3], on the other side, a much smaller relation between
K12 and K22 is suggested.
2. Elastic–plastic interface formulations Elastic limits are defined by functions F1 and F2 as follows:

F 1 ¼ jt t j þ lt n ¼ 0; ð5Þ
Bond behavior is governed by tractions ti at the interface be-
tween concrete and steel. The component tangential to steel bars F2 ¼ t 2t þ t 2n þ ctn ¼ 0: ð6Þ
(tt) corresponds to bond tractions, while the component normal to
Representation of these functions is shown in Fig. 1. Function F1
bars (tn) corresponds to splitting tractions. In the elastic range ti is
 e represents friction caused by mechanical interlock. Adhesion is con-
a function of elastic deformations, defined by the slipping ut and
e sidered negligible because we only consider ribbed bars in the pres-
opening un of the interface, according to the following relation:
ent study. Function F2 describes the situation where the resultant of
t i ¼ Dij uej ; ð1Þ the tractions acts as inclined compressive struts and reaches a limit
stress c. Initially c is equal to the compressive strength of the con-
where, crete, fc. As splitting failure impairs mechanical interlock, it can be
    associated with F1 = 0 and as pull-out failure is linked to crushing
tn uen
ti ¼ ; uei ¼ ; ð2Þ of concrete between ribs, it can be associated to F2 = 0.
tt uet
Observe that friction coefficient l, in this model, is limited to a
and Dij is the elastic matrix. As application to the finite element maximum value equal to 1. Apparently this limitation is related to
method is aimed here, elements are introduced at the interface. function F2, that can not present a decreasing bond strength for a
These elements have four nodes where, initially, two are coincident decreasing compression.
with the other two creating a zero-thickness element. In this case
Dij can be considered a cohesion matrix containing the stiffness of tt
linear springs attached to nodes. Positive values of uen indicates
opening of the interface.
Elastic limits are defined by one or more functions (F), which
can be represented by lines in the space of tractions ti. For loading
μ
along the elastic limits (F = 0, dF = 0), we have also plastic or dissi- F2 F1
pative deformations upi . Incrementally, total deformations can be 1
calculated as:
e p
elastic tn
dui ¼ dui þ dui : ð3Þ c domain
p
Once a flow rule is established to calculate dui , we have a set of
equations that can be solved by return mapping procedures.

2.1. Description of the Lundgren–Gylltoft model

For the sake of completeness, some basic equations of the elas-

tic–plastic formulation proposed by Lundgren and Gylltoft [9] are
reviewed here. Subsequent updates done by Lundgren [24] are also Fig. 1. Representation of elastic domain and functions F1 and F2.
 D.S. Brisotto et al. / Computers and Structures 106–107 (2012) 81–90 83

Loading along elastic limits leads to dissipative or plastic defor-

mation upi . An associated flow rule is assumed for loading along F2
1
0.97 (a)
0.91
(F2 = 0; dF2 = 0), which leads to plastic deformations as follows:

p oF 2
dui ¼ dk ; ð7Þ
ot i
dk is a plastic multiplier. For loading along F1 (F1 = 0;dF1 = 0), a non-

c/fc
associated flow rule is assumed because, in friction problems, dissi-
pative effects are mainly associated to slipping. Then,

p oG
dui ¼ dk ; ð8Þ 0.3
oti
where
G ¼ jt t j þ gt n ¼ 0: ð9Þ
0
0 0.08 0.15 0.25 1
g is a value smaller than 1. A fixed value equal to 0.04 is suggested, d
based on calibrations.
In the Lundgren–Gylltoft model, c and l are only functions of
the plastic deformation, showing a softening behavior with it. In 1 (b)
addition, these functions are pointwise defined. In the next section
we discuss alternative definitions for both functions.

2.2. Redefinition of the parameters c and l

In this work, parameters c and l are considered a function of an

0.55
internal variable d, defined as follows [3]: μ
 p 
u
d ¼ min t ; 1 ; ð10Þ 0.4
lk
where lk is the rib spacing. Such a variable indicates that the dam-
age process is completed once slip grows up to the length of lk, halt-
ing then the evolution of yield surfaces. Here we consider c(d) and
l(d) according to Fig. 2(a) and (b), respectively, which are approxi- 0
mately coincident with the curves proposed in the Lundgren–Gyllt- 0 0.067 0.33 1
oft model, when a 16 mm bar diameter is considered. (Fig. 2(a)
d
corresponds to the uniaxial compression curve for concrete, with Fig. 2. Variation of (a) c/fc and (b) l with damage parameter d.
a factor between plastic strain and damage parameter, while
Fig. 2(b) is obtained through ‘‘ring tests’’: pull-out tests in concrete
e ¼ 0:486cy  0:257/s ; ð11Þ
cylinders confined by steel tubes. The stresses in the steel are used
to evaluate l [9]). If not mentioned otherwise, lk is taken equal to where e is measured radially from the interface and cy is the con-
0.65/s, where /s is the steel bar diameter. It follows that the same crete cover (see Fig. 3). Hoop stresses inside the cylinder defined
plastic deformation of the interface leads to different values of c and by e are not considered. Only stresses between e and cy are relevant
l, depending on the bar diameters. The model, however, is not for bond strength. Then a characteristic length (‘ = cy  e) emerges
intended to be applied to every bar diameter. We suggest here in this case.
applications to diameters ranging from 12 to 20 mm. An average value for hoop stresses r  h in the region defined
The Lundgren–Gylltoft model also depends on the modeling of between e and cy can be defined as follows:
longitudinal cracks in order to simulate correctly the splitting
failure. These cracks are originated by hoop tensile stresses, rh, 1X N
r h ¼ rn ; ð12Þ
and emerge from bars [25]. When applications considered are axi- N n¼1 h
symmetrical, it is necessary to know a priori the number of longi-
tudinal cracks formed. Alternatively, in Dias-da-Costa et al. [26] an where rnh are the smoothed hoop stresses at finite element nodes, N
enriched axisymmetrical element is proposed in order to take into is the number of nodes inside ‘ in a line radially taken at one inter-
account such cracks. In Rabczuk and Belytschko [11] an expression face node (see Fig. 3). For each node of the interface, a value for r
 h is
function of the relation between concrete cover size and /s is used. attributed.
In Lowes et al. [6] a nonlocal modeling is used to incorporate the In order to take into account the effect of hoop stresses on bond
effect of the surrounding concrete-steel state in the interface. strength, we consider l to depend on a second damage factor h as
In the present work we pursue an alternative and original prop- follows:
osition, where hoop stresses at the vicinity of the bar, once greater lðd; hÞ ¼ lðdÞ=h; ð13Þ
than the tensile strength of the concrete, ft, are used as a weaken-
ing factor for friction. This process is nonlocal in the sense that where
bond strength at one particular point of the interface is, in fact, 
r h  ft
dependent on the concrete state in a region defined by some char- h ¼ exp ; ð14Þ
3f t
acteristic length. In order to define it, we consider that longitudinal
cracks do not harm bond capacity as long as they are kept inside and h  i is the Macauley bracket. Observe that Eq. (14) was
some distance from bars. Tepfers [1] approximately defines this calibrated for the cases studied in this work, which correspond to
distance (e) as: concretes with normal strength and toughness.
84 D.S. Brisotto et al. / Computers and Structures 106–107 (2012) 81–90

Un 
a

a m1
symmetry dn ¼ akn ð1  kn Þa1þ1 kn þ 1 ; ð18Þ
line rmax m m
cy   n1
Ut b b
dt ¼ bk ð1  kt Þb1 þ1 kt þ 1 : ð19Þ
smax t n n
rmax is the maximum normal traction at the cohesive zone and smax
e its tangential counterpart. rmax can be considered approximately
equal to concrete tension strength (ft). Finally, kn and kt define a
fraction of dn and dt, respectively, where peak tractions take place.
The cohesive tractions are considered between concrete finite
elements. Besides capturing the development of cracks, they keep
finite elements together in tension (for low levels of stresses)
and prevent interpenetration in compression.

σh
1
σh
2
σh σ3h σNh 4. Numerical experimentation

In this section we simulate different experiments with the cur-

rent model. In all cases studied in this work loading is monotonic
concrete
steel

crescent. Steel is modeled as elastic–plastic, with application of

von Mises yield criterion, associated flow rule and isotropic
hardening. Bulk concrete is considered elastic-linear being all
non-linearity considered in concrete-steel and concrete-concrete
interfaces, as described in Sections 2 and 3, respectively.
interface The following expression was chosen for K12 in all cases:
  K 22
Fig. 3. Representation of nodal hoop stresses rih used to calculate an average K 12 ¼ 0:5 ; ð20Þ
value ðr h Þ at the interface. l
which is an average of values suggested in the literature [9,24,3].
As a main advantage, this procedure eliminates 3D fracture The following expression is used to calculate K22 and is based on
considerations. the tangential force the ribs support elastically [29] (values in mm):
1:687
3. Discrete crack model for concrete K 22 ¼ : ð21Þ
/2s
While the procedure described in Section 2.2 can deal with lon- It is interesting to note that, in Cox and Herrmann [3], elastic con-
gitudinal (splitting) cracks, transversal (primary) cracks that occur stants are also considered inversely proportional to /s, but linear
in the concrete still need to be modeled. Several discrete crack relations are proposed. We comment later on K11 parameter.
models are available specifically for concrete ([27,28], etc). Here In addition, in some comparisons with experiments, the elastic
we decided to use the PPR model [23], due to its greater generality. modulus Ec and tensile strength ft of the concrete are not provided.
For the cohesive crack zone the following tractions Tn and Tt, nor- In these cases the following relations [13] are used (values in MPa):
mal and tangential to the crack plane respectively, take place:  1=3
"     fc þ 8
Un
a m Dn a m Dn m1 Ec ¼ 18275
10
; ð22Þ
Tn ¼  m 1 þ
dn m dn a dn  2=3
fc
 a1  m # ft ¼ 1:4 : ð23Þ
Dn m Dn 10
a 1  þ ; ð15Þ
dn a dn
The Newton–Raphson method is used as the solution procedure. A
tangent matrix obtained by numerical perturbation is used. Only
 n "  n  n1
1 b j Dt j n j Dt j axisymmetrical finite elements are considered.
Tt ¼ n 1 þ
dt n dt b dt Two different groups of reinforced ties are considered in the se-
 b1  n # quence. In the first (Section 4.1), bar is pulled from one side, nor-
j Dt j n Dt mally known as pull-out test. Different modes of failure are
b 1  þ ; ð16Þ
dt b dt captured depending on dimensions and material properties used.
In the second (Section 4.2), bar is pulled from both sides and pri-
where Un is the specific fracture energy for mode I, Dn is the open-
mary cracks appear. This configuration is also known as tension-
ing and Dt is the sliding of the crack surfaces. (In the Eqs. (15) and
pull test.
(16), Un and Ut are considered equal, where Ut is the specific frac-
Parameter K11 in this work is considered independent of bar
ture energy for mode II). For concrete, material constants a and b
diameter with good fitting with experiments. However it showed
are typically greater than 2. Such value gives a steeper descent of
a dependence on boundary conditions. In Section 4.1, a residual va-
the traction versus opening curve after peak and a much less steep
lue of 0.011 mm1 is used for K11 when uen is positive, while, in Sec-
curve near final openings. m, n are calculated according to the
tion 4.2, this value is 0.05 mm1. For uen negative, in all cases,
following expressions:
K11 = 0.15 mm1.
aða  1Þk2n bðb  1Þk2t
m¼ 2
; n¼ : ð17Þ 4.1. Pull-out simulations
1  akn 1  bk2t
dn is the maximum crack opening and dt is the maximum crack slid- Different materials, dimensions, confinements and anchorage
ing and are calculated as: conditions are here considered in order to induce three distinct
 D.S. Brisotto et al. / Computers and Structures 106–107 (2012) 81–90 85

failure behaviors: pull-out failure, splitting failure and yielding of Table 1

the steel bar. Representation of the boundary value problem used Geometrical properties (in mm), see Fig. 4, and concrete properties (in MPa).

in this section is depicted in Fig. 4. /c is the concrete diameter. In References le lc /s /c fc ft Ec

some experimental cases, used here as comparison, the surround- Lundgren [30] 50 100 16 70 35.6 2.7 30.3  103
ing concrete is not cylindrical, as the dashed lines indicate in Fig. 4. Balazs and Koch [31] 80 80 16 160a 25.5 2.0 29.4  103
In these cases the minimum cover cy is used to define /c, as seen in Magnusson [32] 40 40 16 300 27.5 2.2 30.0  103
the figure. lc is the concrete length and le is the embedded or an- Baena et al. [7] 60 200 12 200a 26.5 2.7b 27.6  103c

chored bar length. Sections a and b mark the limits of the embed- a
Quadratic cross section.
b
ded bar and section a corresponds to the closest section to loading. Estimated by the Eq. (23).
c
In the simulations, only half of the problem is modeled due to axi- Estimated by the Eq. (22).

symmetry. Interfacial elements are always uniformly distributed

throughout le. The bar is modeled by four-node quadrilateral finite
elements and the concrete by three-node triangles in a cross-diag- pected behavior of pull-out failure. Interface opening, on the other
onal pattern. side, is kept on a very tight range (maximum value smaller than
0.01 mm) due to confining conditions. Maximum opening occurs
at the peak loading.
4.1.1. Pull-out failure Observe that in all cases studied so far, concrete fracture energy
Experiments from different researchers are tested in this sec- is not considered. The boundary conditions considered are such that
tion as shown in Table 1. In all cases the Poisson’s ratio for con- tensile stresses in the concrete are not large enough to cause macro-
crete, mc, is equal to 0.2. Considering steel, elastic modulus Es is scopic cracks. Failure process is then basically related to crushing of
equal to 200 GPa and Poisson’s ratio ms is equal to 0.3. concrete around bars. Fig. 7 shows the path of tractions for a typical
In all cases concrete cover is sufficient to prevent splitting, ex- point of the interface in the tn/tt space, again for Lundgren’s case
cept in Lundgren [30]. In this case, however, a steel tube with a [30]. From the path followed by tractions, it is clear that initially
1 mm thickness is used to encase concrete, which avoids splitting. they follow surface F1. Pull-out force did not attain the peak yet
Table 2 shows the number of finite elements used in the simula- and slipping is very small. In this phase, friction l does not suffer
tions. In the Lundgren’s case [30] three different meshes are tested any damage (Fig. 2(b)). Tractions then achieve a threshold when
in order to investigate mesh dependency. they reach surface F2. This instant corresponds approximately to
Fig. 5 shows results of the pull-out force versus slipping for the the peak of the pull-out force (Fig. 5), which depends on the com-
experiments and numerical results obtained by the present algo- pressive strength of the concrete. Damage builds up and starts
rithm. Fig. 5(a) shows the cases with bar diameter equal to slowly to decrease parameters c and l. Observe that tractions re-
16 mm and Fig. 5(b) shows the case with bar diameter equal to main approximately at the region of maximum bond (the intersec-
12 mm. After peak loading, the force decreases with slipping until tion of both yield surfaces) throughout the whole damage process,
it is approximately equal to rib spacing. Then behavior becomes indicating that both concrete compressive strength and friction
purely frictional. A reasonably good agreement is achieved in all are involved in the failure process, as expected.
cases, despite some numerical oscillations after peak for some
cases. Also, results do not present mesh dependency, as seen in 4.1.2. Splitting failure
the figure. Confining pressure can change splitting failure to pull-out fail-
Once the cover is large enough to induce pull-out failure, fur- ure when cover is not large enough to prevent the propagation of
ther increase in it has no effect on bond capacity [33]. This is also splitting cracks. In practice, confining pressure can be the result
shown in Fig. 5(a), where a lateral confining pressure is added for of the action of the supports of the beam over the reinforcement,
one of the cases studied [32]. The pressure (5 MPa) has only a min- when the beam ends are supported at their underside.
or effect on pull-out force, as expected. This condition is also sim- In order to demonstrate the potential of the present algorithm
ulated in Lundgren [30] with similar results. Therefore, apart from to capture this effect, we simulate three experiments done by Mal-
material properties, pull-out force in these cases depends basically var [34]. The three tests are identified as P0, 2 and 5, following the
on bar diameter and embedded length, increasing with both. denominations used originally by the author. In the first, no confin-
Fig. 6 shows a typical distribution of slipping throughout the ing lateral pressure p was used and failure by splitting was ob-
length of the bar, for different levels of loading. Lundgren’s case served. In the other two, confining pressure was introduced (see
[30] is considered. It can be observed that, before peak loading Table 3) decreasing the effect of splitting cracks. In all cases
(45 kN), slipping is insignificant or smaller than 1 mm. Surpassing /s = 20 mm, /c = 76.2 mm, le = 61 mm and lc = 101.6 mm (see
peak, slip of the bar reaches values around 5 and 17 mm for loads Fig. 4). Concrete properties are Ec = 37.2 GPa and mc = 0.2. Other
30 and 3.3 kN, respectively. The pattern of slipping is constant, concrete properties are defined in Table 3. Steel properties are
which corresponds to a rigid body motion of the bar. This is an ex- Es = 205 GPa and ms = 0.3. In all cases, 520 volumetric finite ele-
ments and 13 interfacial finite elements are used.
The average bond traction versus slipping curves are presented
in Fig. 8. Experimental and numerical results obtained by the pres-
ent algorithm are shown. Average bond traction is calculated as:
φc P
tt ¼ : ð24Þ
p/s le
φs cy As seen in Fig. 8, bond tractions increase significantly with the
b a
P confining pressure. Case 5 has nearly twice as much peak bond
le stress as case P0. Also, behavior of the case P0 is much more brittle,
with bond tractions dropping sharply to very low values after peak.
lc In cases 2 and 5 drop of tractions is not as sharp and a frictional
bond is still present when slipping is approximately equal to rib
Fig. 4. Representation of the pull-out test. spacing.
86 D.S. Brisotto et al. / Computers and Structures 106–107 (2012) 81–90

Table 2 Fig. 10 presents the path of tractions for case P0 for a typical
Number of volumetric and interfacial finite elements (FEs) used in the simulations. point of the interface. It can be seen that tractions never reach
References Number of vol. FEs Number of interf. FEs pull-out surface F2 and are confined by friction controlled F1. Fric-
Lundgren [30] 665 12 tion l decreases rapidly due to splitting cracks (Eqs. (14) and (13))
2525 25 which explains the lower values of bond tractions in this case.
9752 50
Balazs and Koch [31] 1824 16 4.1.3. Yielding of the steel bar
Magnusson [32] 1524 10 If cover and anchored length of the bar are large enough, failure
Baena et al. [7] 1960 12
of the bar by yielding may occur prior to concrete failure. This sit-
uation was experimentally studied by Magnusson [35], among oth-
ers. Recently the problem was studied in details by Mazzarolo et al.
If we observe the effect of the confining pressure on hoop stres- [16]. Two cases investigated by Magnusson [35] are considered
ses rh in Fig. 9, an explanation for the behavior observed can be in- here. Geometrical dimensions and concrete properties are shown
ferred. In the case P0, Fig. 9(c), rh is larger than ft throughout the in Table 4. In both cases /c = 400 mm, /s = 16 mm and lc = 480 mm
whole cover, indicating a complete fracture of the cover by split- (see Fig. 4).
ting cracks. In the case 5, Fig. 9(a), rh larger than ft is kept inside Steel properties are: Es = 200 GPa, ms = 0.3, initial yielding stress
length e, which is equal to 8.5 mm in this case (Eq. (11)), so split- fy = 580 MPa and linear hardening modulus hs = 880 MPa. In both
ting cracks cease to affect bond. Failure is then basically by pull- cases 2352 volumetric finite elements are used. In cases 1 and 2,
out. Case 2 has an intermediate behavior. the number of interfacial finite elements is 36 and 22, respectively.

70 (a)
65
60 Lundgren [30]
Balazs and Koch [31]
55 Mgnusson [32]
50 + Present model (p=5MPa)
Present model
45
40
P (kN)

35
30
+
25
+ + ++ +
+ + +
20 + +
+ + +
15 + + +
10 +
+ +
+ +
+ +
5+ + +
+ + + + + + +
0 +
0 5 10 15
slip (mm)
70
(b)
65
60
55
50
Baena et al. [7]
45 Present model
40
P (kN)

35
30
25
20
15
10
5
0
0 5 10 15
slip (mm)

Fig. 5. Pull-out force on the bar (kN) versus bar slip (mm) for (a) cases where /s = 16 mm and (b) a case where /s = 12 mm. Experimental and corresponding numerical results
by the present model are shown. In the Lundgren’s case [30], results for three different meshes are shown (Table 2).
 D.S. Brisotto et al. / Computers and Structures 106–107 (2012) 81–90 87

25
Present model
Malvar [34] - test P0
Malvar [34] - test 2
-15 20 Malvar [34] - test 5

average bond traction (MPa)

P=44.7 kN (pre-peak)
slip (mm)

P=30.86 kN (post-peak)
15
-10 P=3.30 kN (post-peak)

10
-5

b a
interface position
0
Fig. 6. Slipping distribution considering different levels of loading, related to
0 2 4 6 8 10 12
sections a and b (see Fig. 4). Lundgren’s case [30] considered.
slip (mm)

Fig. 8. Average bond traction (MPa), according to Eq. (24), versus bar slip (mm).
Experimental and corresponding numerical results by the present model are shown.

-15
before localization, the maximum opening is nearly fifty times lar-
ger than the opening in the pull-out cases (Section 4.1). The result
-10 is a substantial reduction on tn, which by consequence weakens
bond. In Fig. 12(b) is shown the evolution of bond traction tt at sec-
tion a. We see that tt has an abrupt drop after opening reaches a
-5
critical value. Qualitatively this behavior is quite similar to average
tt (MPa)

bond traction obtained in case P0 (Fig. 8), where failure occurred

0 by splitting.
-30 -20 -10 0
tn (MPa)
5 4.2. Tension-pull simulation

Boundary conditions used in previous cases are slightly changed

10 according to Fig. 13. In this case bar is pulled by both ends. This
problem represents the conditions of a reinforced concrete be-
tween two primary (transversal) cracks in a tension zone. The par-
ticular example considered here was experimentally tested by
Fig. 7. Path of tractions for a typical point of the interface. Lundgren’s case [30]
considered (solid lines indicate initial yield surfaces and dashed lines indicate final Doerr [36] and then simulated by Keuser and Mehlorn [8], Prasad
surfaces). and Krishnamoorthy [10], among others. Only one-quarter of the
problem is modeled, due to axisymmetry and mirror symmetry.
Concrete properties are Ec = 35 GPa, mc = 0.2, ft = 2.7 MPa and
Table 3
fc = 37.2 MPa, while steel properties are Es = 200 GPa, ms = 0.33,
Concrete properties and confining pressure (in MPa) according to Malvar [34].
fy = 420 MPa and hs = 2050 MPa. The number of volumetric and
Test fc ft p interfacial finite elements is 974 and 64, respectively. The bar
P0 44.2 4.98 0 and concrete are modeled by four-node quadrilateral finite
2 40.2 4.92 10.3 elements.
5 40.2 4.92 31.0
Doerr [36] shows that, in this case, primary cracks occur at the
middle of the specimen, transversally to the bar. To capture these
cracks, cohesive concrete-concrete elements are inserted in this
As reported in Magnusson [35], in the case 1 final collapse oc- plane. The following fracture properties are considered for the
curs by localization of plastic deformations in the bar and in case cohesive elements (see Eqs. (15) and (16)): Un = 100 N/m,
2 collapse occurs by pull-out, after bar yielding. Both situations a = b = 3, kt = 0.01, kn = 0.005 and rmax = 2.7 MPa. As only mode I
are captured by the present formulation. Experimental and numer- opening is operative, mode II properties are irrelevant.
ical results of pull-out force versus slipping are given in Fig. 11. A Results of the experiments and numerical simulations are
fairly good agreement between experiments and numerical simu- shown in Fig. 14, for three levels of loading: P = 20, 40 and 70 kN.
lations is again obtained. The figure shows the distribution of tensile forces in the bar
Fig. 12(a) shows the variation of the interface opening through- (numerically these forces are obtained integrating stresses in the
out the length of the bar for different loading levels (case 1 consid- bar direction). For the lowest loading (20 kN), concrete remains un-
ered). Three different loading values are considered: a pre-peak cracked and forces are minimum in the bar at the center (x = 0).
value (P = 30.6 kN), a value at beginning of yielding (P = 118.2 kN) Bond is considered perfect in this region. A crack then opens in
and a value right before localization (P = 120 kN). A substantial the concrete at the symmetry line for the intermediate loading
interface opening can be seen near section a (see Fig. 4). The open- (40 kN). Steel bar forces increase in the region due to concrete frac-
ing increases substantially after yielding due to increased reduc- ture. Finally, for the highest loading (70 kN), yielding of the bar
tion of the cross section of the bar by Poisson’s effect. Right takes place at the loading section and at the crack.
88 D.S. Brisotto et al. / Computers and Structures 106–107 (2012) 81–90

(a) (b) (c)

e e

σh
5.0
2.5
-0.0
-2.5
-5.0
-7.5
-10.0
-12.5
-15.0
-17.5
-20.0
-22.6
-25.1
-27.6
-30.1

Fig. 9. Hoop stresses for cases (a) 5, (b) 2 and (c) P0. Values in MPa.

-30
70 kN). In the early stages of loading (P = 20 kN), bond changes
approximately linearly from zero to a maximum value near loading
-20 section. After the central crack is formed (P = 40 kN) a redistribu-
tion of tractions takes place, including a change in the sign of the
-10 bond tractions for 0 < x < 125 mm. Observe that this change corre-
sponds to an unloading of the tractions, even though the external
tt (MPa)

loading is a monotonically increasing function. Bond traction is

0 null for x  125 mm and assumes nearly the same distribution in
-40 -30 -20 -10 0
tn (MPa) the intervals 0 < x < 125 and 125 < x < 250, except with contrary
signs. When loading increases, bond tractions also increase, with
10
maximum values reaching about 15 MPa, for the final loading
(P = 70 kN). Observe also that a second generation of cracks was
20 expected to form at x  125, if the steel bar were still elastic. This
is a typical behavior expected in these cases [1].

30
5. Discussion and final remarks
Fig. 10. Path of tractions for case P0 (solid lines indicate initial yield surfaces and
dashed lines indicate final surfaces). The Lundgren–Gylltoft elastic–plastic formulation of the con-
crete-steel interface is used here to study axisymmetric reinforced
Table 4
concrete failure. The formulation was adapted to be used in a large
Test characteristics according to Magnusson [35]. range of bar diameters and also to incorporate 3D aspects of the con-
crete fracture (longitudinal cracks). Two new damage parameters, d
Test le (mm) fc (MPa) ft (MPa) Ec (GPa)
and h, are introduced in this regard. In particular, parameter h,
1 360 27.6 2.2 30.2 which introduces the damage caused by splitting cracks in the inter-
2 220 30.6 2.4 31.3
face formulation, rendered the formulation nonlocal. In addition,
discrete fracture mechanics is considered to take into account pri-
This is a critical test for the concrete-steel interface model. If its mary cracks, permitting a more refined picture of the failure process.
stiffness is too low, all stresses are carried by the steel and the gra- The model was able to capture some important aspects of the
dient of stresses observed in Fig. 14 tends to disappear. Concrete is concrete-steel behavior observed in experiments from different
not fully loaded and cracks do not occur. If concrete-steel interface sources of the literature, such as:
stiffness is too high, exaggerated gradients are observed in steel
stresses and cracks appear sooner than expected. This is the case  The transition from pull-out to splitting failure due to
when a perfect bond is considered for the concrete-steel interface. changes in confinement.
This conservative approach is the usual assumption made in gen-  The more brittle behavior associated to splitting. The bond
eral in structural analysis. strength reduction, the sudden decrease in bond after the
Fig. 15 shows distributions of bond tractions obtained numeri- peak and the poor residual bond, for relatively small values
cally for the three levels of loading considered in Fig. 14 (20, 40 and of slipping, are all captured by the model.
 D.S. Brisotto et al. / Computers and Structures 106–107 (2012) 81–90 89

(a)
140

120
y
x
100 P P
le= 500 mm

80
P (kN)

60 300 mm 300 mm

40
Fig. 13. Tension-pull test. /s = 16 mm; /c = 150 mm.
Magnusson [35]
20 Present model

present model
0 70 Doerr [36] - 20 kN
0 5 10 15 20
slip (mm) Doerr [36] - 40 kN
Doerr [36] - 70 kN
60
140
(b)

steel force (kN)

50
120
40

100
30

80 20
P (kN)

60 10

50 100 150 200 250

40 x (mm)
Magnusson [35]
Fig. 14. Force distribution along the bar (kN). x is measured from the symmetry line
20 Present model
(Fig. 13). Experimental and corresponding numerical results by the present model
are shown.

0
0 5 10 15 20
slip (mm)  The splitting associated to steel yielding, which is accompa-
Fig. 11. Pull-out force (kN) versus bar slip (mm): (a) test 1; (b) test 2. Experimental
nied by a considerable interface opening and also a sudden
and corresponding numerical results by the present model are shown. decrease in bond, locally.
 The changes in bond tractions that follow the occurrence of
 The transition in the failure process for different anchorage macroscopic cracks in the concrete. As a consequence, the
conditions: from concrete collapse (short anchorage) to methodology is able to localize the sections where succes-
steel collapse (long anchorage). sive generations of macroscopic cracks are formed.

(a) (b)
0.4 -15
interface opening (mm)

P=30.6 kN
bond traction (MPa)

P=118.2 kN
P=120 kN -10

0.2

-5

0 0
b a 0 -1 -2
interface position slip (mm)
Fig. 12. (a) Distribution of the interface opening for different load levels, (b) bond traction (MPa) versus slip (mm) at section a. Test 1 considered.
90 D.S. Brisotto et al. / Computers and Structures 106–107 (2012) 81–90

15 [7] Baena M, Torres L, Turon A, Barris C. Experimental study of bond behavior

between concrete and FRP bars using a pull-out test. Compos Part B Eng
20 kN 2009;40:784–97.
10 40 kN [8] Keuser M, Mehlhorn G. Finite element models for bond problems. J Struct Eng
70 kN 1987;113:2160–73.
bond traction (MPa)

[9] Lundgren K, Gylltoft K. A model for the bond between concrete and
5 reinforcement. Mag Concrete Res 2000;52:53–63.
[10] Prasad MVKV, Krishnamoorthy CS. Computational model for discrete crack
growth in plain and reinforced concrete. Comput Methods Appl Mech Eng
00 50 100 150 200 250 2002;191:2699–725.
[11] Rabczuk T, Belytchko T. Application of particle methods to static fracture of
reinforced concrete structures. Int J Fract 2006;137:19–49.
-5
[12] Ingraffea AR, Gerstle WH, Gergely P, Saouma V. Fracture mechanics of bond in
reiforced concrete. J Struct Eng ASCE 1984;110:871–90.
-10 [13] CEB-FIP, Bulletin No. 213/214, Model Code 1990, London; 1993.
[14] Jendele L, Cervenka J. Finite element modelling of reinforcement with bond.
Comput Struct 2006;84:17801791.
-15 x (mm) [15] Huang K. Modelling the bond between concrete and reinforcing steel in a fire.
Eng Struct 2010;32:3660–9.
Fig. 15. Numerically obtained bond traction distribution along the interface (MPa) [16] Mazzarolo E, Scotta R, Berto L, Saetta A. Long anchorage bondslip formulation
for (a) P = 20 kN, (b) P = 40 kN and (c) P = 70 kN. x is measured from the symmetry for modeling of r.c. elements and joints. Eng Struct 2012;34:330–41.
[17] Ibrahimbegovic A, Boulkertous A, Davenne L, Brancherie D. Modelling of
line (Fig. 13).
reinforced-concrete structures providing crack-spacing based on X-FEM ED-
FEM, and novel operator split solution procedure. Int J Numer Methods Eng
2010;83:452–81.
All these aspects of the failure process are captured essentially [18] Galvez JC, Benitez JM, Casati MJ, Tork BS, Cendon DA. Cohesive-frictional
model for bond and splitting action of prestressing wire. Int J Numer Anal
because the methodology correctly represents, at least qualita- Methods 2011;35:1257–77.
tively, the changes and interactions between the tractions involved [19] d’Avila VMR, Brisotto DS, Bittencourt E. Numerical simulation of cracking in
at the concrete-steel interface. This coupling between tractions is reinforced concrete members by an embedded model. Eng Comput
2008;25:739–63.
an intrinsic characteristic of plasticity models. In addition, none [20] Law DW, Tang D, Molyneaux TKC, Gravina R. Impact of crack width on bond:
of the examples studied presented mesh dependency. confined and unconfined rebar. Mater Struct 2011;44:1287–96.
On the other side, results of the algorithm are critically depen- [21] Vidal T, Castel A, Francois R. Analyzing crack width to predict corrosion in
reinforced concrete. Cement Concrete Res 2004;34:165–74.
dent on elastic parameters K11, K12 and K22. Although suggestions
[22] ACI Committee 446. Finite element analysis of fracture in concrete structures:
are presented to define them, e.g. Eqs. (20) and (21), some curve state-of-the-art report. ACI 446.3R-91. Detroit: American Concrete Institute;
fitting with experiments may be still necessary in order to cor- 1997.
[23] Park K, Paulino GH, Roesler JR. A unified potential-based cohesive model of
rectly attribute values for them. In particular, numerical parameter
mixed-mode fracture. J Mech Phys Solids 2009;57:891–908.
K11 showed to be specially critical for the simulations done. How- [24] Lundgren K. Modelling the effect of corrosion on bond in reinforced concrete.
ever, the experience we had in the present work indicates that the Mag Concrete Res 2002;54:165–73.
transition among different failure modes does not depend on [25] Tepfers R. Cracking of concrete cover along anchored deformed reinforcing
bars. Mag Concrete Res 1979;31:3–12.
changes in these parameters. [26] Dias-da-Costa D, Veludo J, Alfaiate J, Julio E. An element enriched formulation
for simulation of splitting failure. Eng Fract Mech 2011;78:301–16.
[27] Cendon DA, Galvez JC, Elices M, Planas J. Modeling the fracture of concrete
Acknowledgments under mixed loading. Int J Fract 2000;103:293310.
[28] Lens LN, Bittencourt E, d’Avila VMR. Constitutive models for cohesive zones in
mixed mode fracture of plain concrete. Eng Fract Mech 2009;76:2281–97.
We are pleased to acknowledge support from Brazilian Govern- [29] Lundgren K. Modelling of bond: theoretical model and analyses, Report 99:5.
ment through CNPq and CAPES Fellowships. Göteborg: Division of Concrete Structures, Chalmers University of Technology;
1999.
[30] Lundgren K. Pull-out tests of steel-encased specimens subjected to reversed
cyclic loading. Mater Struct 2000;33:450–6.
References [31] Balazs G, Koch R. Bond characteristics under reversed cyclic loading. Otto Graf
J 1995;6:47–62.
[1] Tepfers R. A theory of bond applied to overlapped tensile reinforcement splices [32] Magnusson J. Bond and anchorage of deformed bars in high-strength concrete.
for deformed bars. Ph.D. Thesis. Göteborg: Chalmers University of Technology; Licentiate Thesis. Publication 97:1, N.1113. Göteborg: Chalmers University of
1973. Technology; 1997.
[2] Bazant ZP, Sener S. Size effect in pullout tests. ACI Mater J 1988;85:347–51. [33] Untrauer RE, Henry RL. Influence of normal pressure on bond strength. ACI J
[3] Cox JV, Herrmann LR. Development of a plasticity bond model for steel 1965;62:577–86.
reinforced. Mech Cohes Frict Mater 1988;3:155–80. [34] Malvar LJ. Bond of reinforcement under controlled confinement. ACI Mater J
[4] Salari MR, Spacone E. Finite element formulations of one-dimensional 1992;89:593–601.
elements with bond-slip. Eng Struct 2001;23:815–26. [35] Magnusson J. Bond and anchorage of ribbed bars in high-strength concrete.
[5] Romdhane MRB, Ulm F-J. Computational mechanics of the steel-concrete Ph.D. Thesis. Göteborg: Division of Concrete Structures, Chalmers University of
interface. Int J Numer Anal Methods 2002;26:99–120. Technology; 2000.
[6] Lowes L, Moehle JP, Govindjee S. Concrete-steel bond model for use in finite [36] Doerr K. Ein Beitrag zur Berechnung von Stahlbetonscheiben unter besonderer
element modeling of reinforced concrete structures. ACI Struct J Berücksichtigung des Verbundverhaltens. PhD Thesis. Technischeltoch Schale,
2004;101:501–11. Darmstadt; 1980.