Engineering Structures 32 (2010) 69–79

Contents lists available at ScienceDirect

Engineering Structures
journal homepage: www.elsevier.com/locate/engstruct

Simplified monotonic moment–curvature relation considering fixed-end

rotation and axial force effect
Hyo-Gyoung Kwak ∗ , Sun-Pil Kim
Department of Civil and Environmental Engineering, Korea Advanced Institute of Science and Technology, 373-1 Guseong-dong, Yuseong-gu, Daejeon 305-701, Republic of Korea

article info abstract

Article history: A simple analytical procedure to analyze the reinforced concrete (RC) beams with a cracked section
Received 27 July 2007 is proposed on the basis of the simplified moment–curvature relations of RC sections. Unlike previous
Received in revised form analytical models which result in the overestimation of stiffness and underestimation of structural
19 August 2009
deformations induced from assuming perfect bond condition between steel and concrete, the proposed
Accepted 20 August 2009
Available online 22 September 2009
analytical procedure considers fixed-end rotation caused by anchorage. Furthermore, the proposed
analytical procedure, compared with previous numerical models, promotes the effectiveness of analysis
Keywords:
by reflecting several factors which can influence the nonlinearity of an RC structure into the simplified
Reinforced concrete moment–curvature relation. Finally, correlation studies between analytical and experimental results are
Anchorage slip conducted to establish the applicability of the proposed analytical procedure to the nonlinear analysis of
Fixed-end rotation RC structures.
Moment–curvature relation © 2009 Elsevier Ltd. All rights reserved.

1. Introduction numerical integrations and prediction of material states must be

conducted through entire imaginary layers, and the stiffness of
Since the structural behavior of reinforced concrete (RC) struc- structures is usually overestimated by ignoring the anchorage slip
tures represents complex features according to the variation of nu- of the main bar at the beam–column joint. The modified stiffness
merous influencing factors including concrete cracking, nonlinear approach, on the other hand, is based on an overall simplified mo-
properties of concrete, and yielding of reinforcing steel, an accurate ment–curvature relation reflecting the various stages of material
prediction of nonlinear responses in RC structures requires an ex- behavior [8–10]. In this approach, stiffness is determined as the
act consideration of these nonlinear effects. However, this is time slope of the moment–curvature relation which is constructed by a
consuming and makes the nonlinear analysis of an entire struc- section analysis, and a nonlinear analysis of RC beams can be con-
ture impossible. Nevertheless, in recent years, the increasing need ducted based on this moment–curvature relation. Therefore, this
for the assessment of the strength and serviceability of existing method, in comparison with the layered section approach, can re-
structures and newly designed critical structures has encouraged duce the calculation time and storage space in the nonlinear anal-
the application of nonlinear analysis and, in advance, the devel- ysis of RC beam. Nevertheless, this approach established, to date,
opment of advanced analytical methods capable of representing still has limitations for simulating the structural behavior of RC
the cracked behavior of RC structures under all possible loading structures because of the exclusion of the fixed-end rotation effect
conditions. caused by the anchorage slip in the beam–column joint.
Material nonlinear analyses of RC beams that comprise the pri- To address these limitations, an improved nonlinear approach
mary members of RC structures are generally conducted by using is introduced in this paper on the basis of the simplified mo-
one of the two basically different approaches: (1) the layered sec- ment–curvature relation of RC sections, which is uniquely defined
tion approach (see Fig. 1); and (2) the modified stiffness approach according to the dimensions of the concrete section, the material
(see Fig. 2). The former is based on the idealized stress–strain rela- properties of concrete and steel, and the axial force. However, ax-
tions for concrete and reinforcing steel in which a finite element ial force variation expected to be very small as compared with
is divided into imaginary concrete and steel layers in a section. axial force by dead load and super-imposed dead loads is not
This approach has been broadly used by many previous investiga- considered in this paper because only beams and columns sub-
tors [1–7]. However, the application of this approach to large struc- jected to constant axial force are focused on. In addition, to sim-
tures with many degrees of freedom has some limitations, because ulate the concentrative deformation at the beam–column joint,
equivalent stiffness (EIeq ) in the plastic hinge length is introduced
in the finite element formulation. Finally, validity of the proposed
∗ Corresponding author. Tel.: +82 42 869 3621; fax: +82 42 869 3610. algorithm is established by comparing the analytical predictions
E-mail address: khg@kaist.ac.kr (H.-G. Kwak). with experimental results and previous analytical studies.

0141-0296/$ – see front matter © 2009 Elsevier Ltd. All rights reserved.
doi:10.1016/j.engstruct.2009.08.017
70 H.-G. Kwak, S.-P. Kim / Engineering Structures 32 (2010) 69–79

εcm Cc
d’ di ε’s Cs
element i φ εci

t d

εs Ts

Fig. 1. Layered section approach.

M
EI′x 1
My
Eb2
EIx
EIx 1
1 EI ′x
1 φ
1 Eb1
1 EIx
EIx u
u1
My
Fig. 4. Bond stress–slip relation.

Fig. 2. Modified stiffness approach.

in improving the numerical stability of the solution. More details
2. Material models for the material models for concrete and steel can be found else-
where [12–15].
To accomplish section analysis to define moment–curvature re- In contrast to the concrete and steel whose material properties
lation, material models for concrete and steel must be defined are uniquely defined, the bond stress–slip relation depends on the
previously. Many mathematical models for concrete are currently relative deformation of concrete and steel. As well known, the rela-
used in the analysis of RC structures. Among these models, the tion between bond stress and slip depend on many factors includ-
monotonic envelope curve introduced by Kent and Park and later ing location, surface condition, and load step and so on. Therefore,
extended by Scott et al. [11] is adopted in this study for a layered it is practically impossible to establish a local bond stress–slip rela-
section approach for its simplicity and computational efficiency. tion. Concurrently, average bond stress–slip relation which is mea-
In this model, as shown in Fig. 3(a), the monotonic concrete sured over equal to five bar diameter is used. Moreover, the result
stress–strain relation in the compression region is described by is very sensitive to experimental error because the bond stress is
three regions. However, the stress–strain relation in the tensile derived from the change in steel stress, and the bond–slip relation
region is ignored in strength calculations in this study, because also depends on the position of the bars, the surface condition of
concrete has low tensile strength, generally less than 20% of the the bars, the loading stage, the boundary conditions and the an-
compressive strength, and it makes a negligibly small contribution chorage length of the bars. In spite of these difficulties, several ex-
to the strength and energy absorption capacity of an cracked RC perimental bond stress–slip relations have been proposed [16–19].
section. Reinforcing steel is modeled as a linear elastic, linear strain In this study, the bilinear bond stress–slip model in Fig. 4 is adopted
hardening material with yield stress σy , as shown in Fig. 3(b). The which is simplified based on bond stress–slip model [16] with bond
reasons for this approximation are: (1) computational convenience slip and stress at yielding state (u1 = 0.1 cm, τ1 = 14.85 MPa)
of the model; and (2) the behavior of RC members is greatly af- obtained from experimental studies. This model gives good ap-
fected by the yielding of reinforcing steel. It is, therefore, advis- proximation of the actual behavior for cases that do not exhibit
able to take advantage of the strain hardening behavior of steel significant bond slip and associated bond damage.

σc σ
f c′ Kf ′c CONFINED CONCRETE
ES2 1
UNCONFINED CONCRETE σy
Z2
ES1
Z1 – εu 1
ε
0 εu
0.2 Kf ′c

εc
ε0=0.002 εu – σy

(a) Concrete. (b) Steel.

Fig. 3. Stress–strain relation.

 H.-G. Kwak, S.-P. Kim / Engineering Structures 32 (2010) 69–79 71

M relations corresponding to various levels of axial force. Up to reach-

ing the balanced axial force Pb , the yield moment of a section in-
creases in proportion to the axial load P and the failure curvature
My considerably decreases by the presence of the axial load. In par-
B
ticular, the application of large seismic loads acting horizontally
requires an RC column to resist a load combination of a relatively
EI
Mcr A small axial load and large moment. Generally, a load combination
of (P , M ) will be located in a region upper-bounded by the axial
load of 0.5Pb at a P − M interaction diagram, where Pb is the bal-
φ anced axial force (see Fig. 6(b)).
φcr φy To implement the axial force effect into a moment–curvature
relation, a linear interpolation between two boundary points of
Fig. 5. Typical moment–curvature relation of an RC section. (0, MyB ) and (0.5Pb , 0.5Mb ) is used on the basis of the assumption
that the moment capacity is linearly proportional to the applied
3. Basic moment–curvature relation axial load P until P reaches to 0.5Pb (see Fig. 7). In reality, the
axial force acting on an RC column continually changes and this
To ensure ductile behavior in practice, steel contents less than
could be significant especially when a frame structure is subjected
the balanced design value are always used for flexural members.
to horizontal load, because lateral drift accompanying the floor
The typical moment–curvature relation for a lightly reinforced
concrete section with one top and one bottom layer of reinforcing rotation will occur. Nevertheless, this effect is not considered in
steel (under-reinforced concrete section) can be idealized to the this paper because only beams and columns subjected to constant
trilinear relation shown in Fig. 5. The first stage is to the cracking axial force are focused in which the variation is expected to very
(point A in Fig. 5), the second to yield of the tension steel (point B small as compared with axial force by dead load and super-impo-
in Fig. 5), and the third to the limit of useful strain in the concrete. sed dead loads.
The moment and curvature at first yield of the tensile steel
(point B in Fig. 5) should also be calculated using the defined 4.2. Anchorage-slip effect
stress–strain relations for concrete and steel. Based on normal
force equilibrium, a section analysis is carried out by first assuming Since bond stresses in RC structures arise from changes in the
that the tension steel reaches the yielding point. (1) During the steel stress along the length, the influence of bonding becomes
section analysis, for fixed curvature, assumed neutral axis should more pronounced in a cracked region. In a simplified analysis of
be changed until the difference between internal tensile force T RC structures, complete compatibility of strains between concrete
and compressive force C is less than the given tolerance. And then, and steel is usually assumed, which implies a perfect bond
if final steel strain is less than the yielding strain, this iteration is condition. However, this assumption is realistic only in regions
repeated again for an increase of curvature. Through successive where negligible stress transfer between the two components
iteration, the moment and curvature at steel yielding are finally takes place. In regions of high transfer stresses along the interface
determined. The connection from initial point to yield point gives
between reinforcing steel and surrounding concrete, such as near
the simplified basic moment–curvature relation defined in this
cracks, the bond stress is related to the relative displacement
paper. In advance, the moment–curvature relation to the post-
between reinforcing steel and concrete. Therefore, the bond-slip
yielding stage is approximated as a straight line which conserves
effect must also be considered to simulate the structural behavior
same energy until ultimate point.
more exactly. In this regard, many studies have been carried out to
4. Modification of basic moment–curvature relation simulate this effect based on finite element analysis [16,19,20] and
experimental results [17,18].
4.1. Axial force effect Two different elements, namely the bond-link element and the
bond-zone element, have been proposed to date for inclusion of
An axial force has an important effect on the moment–curvature the bond-slip effect in the finite element analysis of RC structures
relation of an RC section. Fig. 6 shows the moment–curvature [21,22]. However, the use of these elements requires a double

0.0 0.5 1.0 1.5 2.0 2.5 3.0

Curvature (x10-3rad/cm) Moment (ton-cm)
(a) Moment–curvature relation. (b) P-M diagram.

Fig. 6. Axial force effect on moment–curvature relation.

72 H.-G. Kwak, S.-P. Kim / Engineering Structures 32 (2010) 69–79

short-span beams may be greatly affected both by the shear and

also by the details of anchoring the beam reinforcement. In par-
ticular, slippage of the main bars from the anchorage zone (∆total
in Fig. 8(a)) accompanies the rotation of the beam fixed-end, θFE .
This cannot be simulated by any mechanical model, and the rigid
body deformation, which accounts for approximately 50% of the
total deformation [23], may increase as the deformation increases.
This phenomenon will be enlarged in the case of an under-rein-
forced concrete beam. Consequently, its exclusion may lead to an
overestimation of both the stiffness and the energy absorption
capacity of the structure. Accordingly, to simulate the structural
behavior more exactly, it is necessary to take into account the
fixed-end rotation caused by the anchorage slip.
The anchorage slip ∆total at a beam–column joint can be divided
into two slip components according to the difference in the
bond-slip behavior; ∆axial induced from pull-out of the anchored
reinforcement embedded through a column at beam–column joint
and/or footing at column-footing joint and ∆bending accompanied
Fig. 7. Modification of moment–curvature relation considering axial force effect.
by the enlargement of bending cracks at the end face of a beam.
Since ∆axial and ∆bending are basically caused by the axial force
and bending moment, respectively, the anchorage slip ∆total can
node to represent the relative slip between reinforcing steel and
be found by superposition of ∆axial and ∆bending after separate
concrete. In a beam structure defined by both end nodes along
calculation of each slip component.
the direction, it is impossible to use the double node at each end
node. To address this limitation in adopting the bond model, a nu- ∆total = ∆axial + ∆bending . (1)
merical algorithm that includes the bond-slip effect in the mo-
ment–curvature relation is proposed in this study. 4.3. Calculation of bond-slip component ∆axial
Unlike the critical region located in the vicinity of the beam
mid-span, as well as at the ends of long-span beams, the behav- A part of an RC member subjected to uniaxial tension is shown
ior of the critical region at the beam–column joint of relatively in Fig. 9. When the axial load Ph is applied, the far ends represent

Δtotal

ΘFE

Anchored h
d
Bars
c

δ FE

(a) Anchorage slip. (b) Fixed-end rotation.

Fig. 8. Rigid body deformation at the beam–column joint.

εx
εsx
εcx εs0
ε sl
Ph Ph

ld ld

(a) Tension member. (b) Strain distribution.

ε sx - ε cx

Δx

(c) Bond-slip distribution.

Fig. 9. Behavior of tension member.

 H.-G. Kwak, S.-P. Kim / Engineering Structures 32 (2010) 69–79 73

the straight portion in front of the hook is assumed to be ∆0 , then

x=1d, dΔx/dx=εsx-εcx the constant C2 will be ∆0 . In advance, as shown in Fig. 9(b), the
x=0, Δx=Δ0 strain difference at the column face that represents the maximum
value can be determined from d∆x /dx = εsx − εcx by substituting
Ph the anchorage length ld into x:
d∆
(x = ld ) = k1 C1 cosh(k1 ld ) + k1 ∆0 sinh(K1 ld ) = εsl − εcl . (4)
dx
hooked Since the concrete strain εc is negligibly small, compared with
anchorage the steel strain, the other constant C1 can be calculated as C1 =
of steel BEAM (εsl − k1 ∆0 sinh(k1 ld ))/(k1 cosh(k1 ld )). Accordingly, the general
solution is finally determined, and the steel strain at x = 0 can
COLUMN
be computed as:
d∆ εsl
Fig. 10. Embedded steel in RC column. εs0 = (x = 0) = − ∆0 k1 tanh(k1 ld ). (5)
dx cosh(k1 ld )
a fully cracked state with a steel strain of εsl , which means there The behavior of hooked bars in tension was studied experimen-
tally by Soroushian et al. [24]. They proposed the following relation
is no load carrying by concrete. In advance, the tensile force Ph
between the tensile force (Ph ) and the slip at ∆0 the starting point
is transferred from the steel bar to the concrete by bond stress,
of the hook (x = 0) on the basis of experimental observations:
and the value of the bond stress is zero at the inner end of the
0.2
transfer length ld . This means that there is no bond slip at the ∆0

central position. Moreover, it can be assumed that the strains in Ph = Es εs0 As = Phu , Phu = 271(0.05db − 0.25) (6)
2.54
steel and concrete are equal to each other at x = 0, and the strain
value corresponds to εs0 . where the force term Ph is in kN, and the slip ∆0 and bar diameter
From the strain distribution, the local slip ∆x can be defined as db are in mm.
the total difference in elongation between the reinforcement and Now the slip of anchored reinforcement, ∆axial , can be
calculated by substituting ∆0 determined from Eqs. (5) and (6) into
the concrete matrix measured over the length to a distance x from
the general solution. The slip when the reinforcement yields at the
the mid-span.
column face (x = ld ) can finally be represented by
Z x
∆x = (εsx − εcx )dx εsy ∆0
(2) ∆axial = tanh(k1 ld ) + . (7)
0 k1 cosh(k1 ld )
where εsx and εcx are the strain distributions of steel and concrete, In the case of interior beam–column joint in which anchorage
respectively (see Fig. 9). length ld is enough and/or beams at both sides are bended
On the basis of force equilibrium and the relation of Eq. (2) symmetrically with respect to interior beam–column joint, ∆axial
with the assumption that local bond stress is linearly dependent can be assumed as zero. Therefore, it is necessary to calculate ∆axial
on the local slip, the following well-known governing differential and Eq. (7) can be simplified as
equation for bond slip can be obtained. εsy
∆axial = tanh(k1 ld ). (8)
k1
d2 ∆x Eb Σ0 (1 + ηρ) d2 ∆x
− ∆x = − k21 ∆x =0 (3)
dx2 As Es dx2
4.4. Calculation of bond-slip component ∆bending
where η = Es /Ec , the steel ratio ρ = As /A − c, Σ0 is the perimeter
of the steel bar, Eb is the slope of the bond stress–slip relation, and Fig. 11 shows an infinitesimal beam element of length dx with
E and A mean elastic modulus and sectional area, respectively. axial force. If the linear bond stress–slip relation is assumed as
The general solution of Eq. (3) is given by ∆x = C1 sinh k1 x + shown in Fig. 4, the variations of the axial force components of steel
C2 cosh k1 x, where C1 and C2 are constants that have to be located at the tension (dPs ) and compression (dPs0 ) region can be
determined from the boundary conditions. Focusing on the hooked represented by [25]:
anchorage of reinforcing bars in Fig. 10, bond slip represents the dPs = fb · dx · Σ0 = Eb · ∆ · dx · Σ0 = NAs · Es · dεs (9)
maximum value at the point where the bar enters the column
at beam–column joint and/or footing at column-footing joint and dPs = N As · Es · dεs = N As · Es (c − d )dφx
0 0 0 0 0 0 0 0 0
(10)
slight bond slip may also occur at the starting point of the straight where N is the number of steels, fb is the bond stress at the steel
portion in front of the hook. If the bond slip at the starting point of interface and ()0 denotes the properties of steel in the compression

σc σc +d σc

P′s P′s + dP′s

Paxial Paxial

Ps Ps + dPs
fb

dx Ps Ps + dPs
(a) RC element. (b) Steel interface.

Fig. 11. Free body diagram for an RC element.

74 H.-G. Kwak, S.-P. Kim / Engineering Structures 32 (2010) 69–79

d′ curvature at x = l and can be calculated on the basis of the force

N′A′s equilibrium of Pcl + Psl0 = Psl + Paxial at the beam–column joint, in
c which Paxial is the applied axial load. Finally, the bond slip induced
from the bending behavior, ∆bending , is represented by:
y
ϕx (1 − γ (1 + λ))εsy
d
∆bending = C4 cosh(k2 l) = cosh(k2 l)
k2 sinh(k2 l)
(1 − γ (1 + λ))εsy
= (15)
NA s k2 tanh(k2 l)
where λ = Paxial /Psl .
Since Eq. (15), however, includes an unknown variable γ (see
Eq. (13)), the bond slip is not directly determined. This means that
the introduction of a formula related to the bond slip is required; a
Fig. 12. Section geometry. relation introduced by Gergely and Lutz [26], being popularly used
to determine the allowable crack width of cracked RC structures, is
part (see Fig. 12). Since the longitudinal strains are directly propor- adopted in this study:
tional to the distance from the neutral axis of zero strain, the vari- p
ation of curvature φx representing the gradient of the strain profile wmax = 1.08βc fs 3 dc A × 10−5 (mm) (16)
at the section is also the variation of strain. where βc is the ratio of distances from the tension face and from the
Until the reinforcing steel reaches its yielding stress, the con- steel centroid to the neutral axis, dc is the thickness of the concrete
crete strain at the extreme compression fiber is not large enough cover, and A is the concrete area surrounding one bar.
to assume that the elastic modulus of concrete Ec is constant across Because it can be assumed that crack width w is equivalent to
the concrete section. By ignoring the tensile force in the concrete, two times the bond slip at the cracked location (0.5w = ∆bending ),
which makes a negligibly small contribution after cracking, the the unknown variable γ can be determined from Eqs. (15) and
variation of compression force of concrete in terms of curvature (16) with assumed βc , which has a recommended value of 1.2
φx can be determined as: in a beam without axial force. Then, the neutral axis depth c is
Z Z c calculated from the relation of r = (1 −α)/{α 2 /2ηρ +β(α− d0 /d)}
1
dPc = Ec dεc dAc = Ec · dφx · y · bdy = Ec · b · dφx · c 2 (11) in Eq. (13). In advance, βc is recalculated and again compared
0 2 with the assumed value, and these iterations are repeated until
where b is the width of section, and c is the distance from the the difference between βc assumed and calculated is less than
extreme compression fiber to the neutral axis. Namely, the vari- the given tolerance. Generally, the neutral axis depth c decreases
ation of curvature along the length can be expressed by dφx /dx = with an increase of the bending moment acting on a section.
Eb · Σ0 · ∆x /(Ec · bc 2 /2 + A0s Es0 (c − d0 )) on the basis of the force Nevertheless, a constant value of c determined on the basis of the
equilibrium condition, dPc + dPs0 = dPs (see Fig. 11). The variation yielding of reinforcing steel at the beam–column joint (x = l) is
of concrete strain at the steel location can be simplified in terms of assumed as a representative constant value along the entire span
curvature φx as dεc /dx = (d − c ) · dφx /dx. length for the computational convenience, because the neutral
When the bond slip ∆x at the steel–concrete interface is axis depth maintains an almost constant value from the initial
defined by the relative displacement between reinforcing steel and cracking up to the yielding of reinforcing steel and the behavior
concrete (∆ = us − ux ), the first order and second order differential of the beam–column joint mainly affect the global behavior of the
equations of bond slip lead to: structures.

d∆x dus duc 5. Calculation of equivalent stiffness

= − = εs − εc (12)
dx dx dx
  To account for the fixed-end rotation in this study, the reduced
d ∆x
2
d εs d εc 4Eb 1−α stiffness, EI, in the moment–curvature relation for the elements
 
= − = 1−   ∆x located at the ends of the beam within the range of the plastic
dx2 dx dx Es ds  α2
+β α− d0 
2ηρ d hinge length Lp is used. Among the various empirical expressions
= k21 (1 − γ )∆x (13) have been proposed, in the case of no axial load, the relatively
simple equation of Lp = 0.25d + 0.075z, proposed by Sawyer [27],
N 0 A0s Es0
where η = EEs , ρ = NA bd
s
,β = α = dc . NAs Es
, can effectively be used, where d and z are the effective depth of
c
Accordingly, the following governing differential equation for a section and the distance from the critical section to the point
the bond slip, which has the same form for an axial member (see of contraflexure, respectively. However, since the plastic hinge
Eq. (3)), is obtained. length increases in proportion to the applied axial load, it may be
difficult to estimate the plastic hinge length of an axially loaded
d2 ∆x member with only this simple equation. Accordingly, the plastic
− k22 ∆x = 0. (14) hinge length of Lp = xh, proposed by Bayrak and Sheikh [28], is
dx2
introduced as an upper limit value, where h is the section depth
The general solution to Eq. (14) is given by ∆x = C3 sinh k2 x + and x is an experimental parameter ranging from 0.9 to 1.0, and
C4 cosh k2 x, in which C3 and C4 are constants that have to be an inequality condition of Lp ≤ d is designated in this paper for
determined from the boundary conditions, and k22 means K12 (1−γ ). simplicity of analysis.
For a beam structure, it can be assumed that the steel and concrete If a beam with rotational stiffness kθ at both ends is subjected to
strains at the end of a structure (x = 0) are zero, because the a horizontal force P, as shown in Fig. 13(a), the corresponding hori-
bending moment by the applied lateral loads is zero. From the zontal drift ∆l in which shear deformations are being neglected for
boundary condition of d∆/dx(x = 0) = εs0 − εc0 ∼ = 0, a constant simplifying formulation can be obtained by
C3 is consequently determined as zero. The constant C4 can also
be determined from the boundary condition of, d∆/dx(x = l) = PL3 PL2
C4 k2 sinh(k2 l) = εsl − εcl = εsl − (d − c )φl , where φl denotes the
∆1 = + (17)
12EI 2kθ
 H.-G. Kwak, S.-P. Kim / Engineering Structures 32 (2010) 69–79 75

P P
EIeq Lp
Kθ

EI L EI L2

Kθ
EIeq Lp

(a) Simplified model with kθ . (b) Simplified model with EIeq .

Fig. 13. Consideration of the equivalent stiffness.

M 6. Solution algorithm

My For the analysis of RC structures, Timoshenko beam theory was

used in this study [29]. This theory is well established and widely
used in the analysis of beams, and details for the formulation of a
EI
beam element can be found elsewhere [29]. (3) In a typical Tim-
oshenko beam, it is usual to assume that cross sections normal
EIeq to the beam axis remain normal before deformation, but this is
not necessarily after deformation. In addition, the effects of shear
deformations are not taken into consideration in simulating non-
linear behavior since the normal bending stresses reach a maxi-
mum at the extreme fibers, where the transverse shear stresses are
φ at their lowest value, and reach a minimum at mid-depth of the
0 φy
beam, where the transverse shear stresses are highest. Thus, the
interaction between transverse shear stresses and normal bending
Fig. 14. Modification of monotonic envelope curve.
stresses is relatively small and can be ignored. This means that the
flexural rigidity EI is replaced by that corresponding to the curva-
where the first term is the contribution by the bending deforma- ture calculated from the nodal displacements by φ = (θi − θj )/L,
tion of the beam and the second term by the end rotational stiffness whereas the shear rigidity of a beam element GA is assumed to be
kθ . This rotational stiffness can be determined on the basis of ∆total unchanged, where θi and θj denote the rotational deformations at
calculated by ∆total = ∆axial + ∆bending . both end nodes, and L is the element length. Every nonlinear anal-
ysis algorithm consists of four basic steps: formation of the current
∆total My
θfe = , kθ = (18) stiffness matrix, solution of the equilibrium equations for the dis-
d−c θfe,y placement increments, state determination of all elements in the
where My is the moment and θfe,y is rigid body rotation at the model, and a convergence check. Since the global stiffness matrix
beam–column joint when reinforcing steels yield, respectively. of the structure depends on the displacement increments, the so-
When the same force acts on a beam with reduced stiffness EIeq lution of equilibrium equations is typically accomplished with an
at both ends, as shown in Fig. 13(b), the horizontal deflection ∆2 iterative method through the convergence check. The nonlinear so-
can also be calculated by the moment area method. lution scheme selected in this study uses a tangent stiffness matrix
at the beginning of each load step in combination with a constant
P {EIeq L32 + 2EIlp (4l2p + 6lp L2 + 3L32 )} stiffness matrix during the subsequent correction phase; that is,
∆2 = . (19) the incremental-iterative method.
12EIeq EI
All the remaining algorithms, from the construction of an ele-
From the equality condition of ∆1 ≡ ∆2 , the equivalent sti- ment stiffness matrix to the iteration at each load step, correspond
ffness EIeq can be determined by with those used in a classical nonlinear analysis of RC structures.
More details can be found elsewhere [16,30,18,20].
1 1 1
= + (20)
EIeq β · kθ · L EI 7. Numerical application

where β = α(1 − 2α + 4/3α 2 ), α = Lp /L. In order to establish the validity and applicability of the pro-
The same derivation procedure for a cantilevered beam is posed moment–curvature relation, correlation studies between
applied and the equivalent stiffness EIeq obtained in this case has analytical result (4) including layered section approach by authors
the same form as Eq. (20) except the parameter β has the form of and experimental studies are conducted. Among the many ex-
β = α(1 − α + 1/3α 2 ). The bending stiffness of elements located perimental results available in the literature, five RC specimens
within the plastic hinge length Lp from both end faces of a member are investigated and discussed, as these specimens represent typi-
will be represented by EIeq instead of EI which is still used at cal structural behaviors according to various effects such as steel
the other region. This means that the slope of moment–curvature ratio, boundary condition, and application of axial force. These
relation is finally modified (see Fig. 14). specimens are BEAMR6 and BEAMR4 experimented by Ma [31],
76 H.-G. Kwak, S.-P. Kim / Engineering Structures 32 (2010) 69–79

P
3.49

4#6

#2 40.64
1.91

4#6
14@10+18.9 = 158.9
1.91
22.86

(a) Details of BEAMR6.

P
3.49

4#6

#2 40.64
1.91

3#5
158.75
1.91
22.86

(b) Details of BEAMR4.

P
5.08

2#6

#2 30.48
3.81

2#6
5.08 87.6
15.24

(c) Details of BEAMS1.

Fig. 15. Details of beam members (unit: cm).

Table 1
Material and sectional properties used in the application.
SPECIMEN Ec (MPa) Es (MPa) fc0 (MPa) fy (MPa) ρ (Ast /bd) ρ 0 (Ac /bd) P (kN)

BEAMR6 24,804 200,569 31 451 0.014 0.014 0

BEAMR4 22,785 200,569 30 451 0.014 0.007 0
BEAMS1 27,770 200,569 34 496 0.012 0.012 0
COLUMN1 23,975 200,569 26 496 0.012 0.012 178
COLUMN2 30,447 200,569 42 496 0.010 0.010 45

specimen 00.147 (BEAMS1) and 40.048 (COLUMN1) experimented In these specimens, the plastic deformation is concentrated at the
by Wight and Sozen [32], and specimen 1 (COLUMN2) experi- end of a beam with narrow width, accompanying fixed-end rota-
mented by Low and Moehle [33]. The material and sectional prop- tion that occurs in addition to elastic rotation at the cracking stage.
erties of each specimen are summarized in Table 1. To simulate more exact structural behavior with the beam element
The geometry and cross-sectional dimensions of the three beam formulated on the basis of the average deformation in an element,
members are presented in Fig. 15. The first two specimens are con- a separate consideration of this region is required in the finite el-
nected to columns with hook in the anchorage of main reinforcing ement modeling. This is necessary because the ultimate capacity
bars, and bending mechanism is largely responsible for the struc- may be overestimated if the plastic hinge length is not precisely
tural behavior. On the other hand, the last specimen is connected to taken into consideration. Since the calculated plastic hinge length
RC column without hook and the structural response of this spec- Lp is determined as 20 cm for BEAMR6 and BEAMR4, and 15 cm for
imen appears to be affected more by shear force than the other BEAMS1, the specimens are modeled along the entire span with an
beam specimens because of its relatively small span to depth ratio. element of L = 10 cm except a right end element.
 H.-G. Kwak, S.-P. Kim / Engineering Structures 32 (2010) 69–79 77

150 120

100 80

Load(KN)
Load(KN)

50 40

Experiment Experiment
Layered Section Approach Section approach
Layered section Approach
This Study Study
This study

0 0
0 1 2 3 4 5 0 1 2 3 4 5
Displacement(cm) Displacement(cm)

(a) BEAMR6. (b) BEAMR4.

100
Load(KN)

50

Experiment
Layered section approach
This study

0
0 1 2 3 4
Displacement(cm)

(c) BEAMS1.

Fig. 16. Load–deflection relation for beam members.

Fig. 16 compares the load–deflection relations obtained by because the introduced numerical model does not reflect the stiff-
a layered section approach and the proposed model with the ness degradation caused by shear cracks accompanied with bend-
experimental results of BEAMR6, BEAMR4 and BEAMS1. The laye- ing cracks. Nevertheless, the proposed model can effectively be
red section approach gives very satisfactory predictions of the used in simulating the nonlinear response of RC beams.
elastic behavior before initial cracking and of the value of yielding The next two specimens, COLUMN1 and COLUMN2, are selected
moment itself, but the displacements corresponding to the yielding to demonstrate additional effects by the axial load, and details of
moment are underestimated. This means that the layered section these specimens are presented in Fig. 17.
approach, which adopts a perfect bond assumption, is limited in Since λ in Eq. (15) and βc in Eq. (16) represent the different
terms of describing the cracking behavior of RC beams that are values as the axial load acts (λ = 0.631, βc = 1.4 for COLUMN1,
accompanied by a fixed-end rotation concentrated at the end of λ = 0.484, βc = 1.24 for COLUMN2), ∆bending calculated in Eq. (15)
the beam. will be changed, and RC columns show different slip behavior from
In contrast, the introduced numerical model, which considers that of RC beams. In addition, the specimens are modeled along the
the fixed-end rotation effect according to the aforementioned entire span with an element of l = 5 cm on the basis of the plastic
modification procedure, provides good agreement with experi- hinge length, excluding a right end element.
mental results through the entire loading steps. The elastic stiff- The responses represented in Fig. 18 compare the load–defle-
nesses before the initial cracking are underestimated because ction relations obtained by the layered section approach and the
the fixed-end rotation is taken into account by changing the av- proposed model with the experimental results. The results of
erage bending stiffness EI without any further consideration of the present study display very satisfactory agreement with the
uncracked or cracked section states (see Fig. 14). However, these measured data. Meanwhile, the numerical results not considering
differences at the elastic loading steps do not represent a remark- the axial force effect represent an underestimation of the stiffness
able influence on the entire structural response ranging from ini- and the ultimate strength of RC columns. These differences will
tial cracking to large deformations after yielding of steel. In spite be enlarged as the magnitude of the axial load increases. Finally,
of relatively accurate consideration of the fixed-end rotation ef- the proposed analytical method can also be effectively used to
fect in the introduced numerical model, a slight difference be- calculate the nonlinear behavior of RC columns subject to axial load
tween the analytical response and the experimental data still exist, as well as RC beams without axial force.
78 H.-G. Kwak, S.-P. Kim / Engineering Structures 32 (2010) 69–79

P
5.08

2#6

#2 30.48 P
3.81
2#6
5.08 87.6
15.24

(a) COLUMN1.

P
16.51

1.27

1.27 P
1.27

51.5
#2 #6
(b) COLUMN2.

Fig. 17. Details of column members (unit: cm).

120 30

80 20
Load(KN)
Load(KN)

40 10
Experiment Experiment
Layered section approach Layered section approach
This study This study
with axial force effect with axial force effect
without axial force effect without axial force effect

0 0
0 1 2 3 4 0 1 2 3
Displacement(cm) Displacement(cm)

(a) COLUMN1. (b) COLUMN2.

Fig. 18. Load–deflection relation for column members.

8. Conclusions model can be effectively used to simulate structural response of

RC columns.
This paper concentrates on the introduction of nonlinear analy-
Acknowledgment
sis of RC beams and columns based on simplified moment–curva-
ture relations. Unlike most mathematical or mechanical models
This research was supported by a grant the Smart Infra-
found in the literature, the proposed model has taken into account
Structure Technology Center funded by the Korea Science and
fixed-end rotation caused by anchorage slip at the fixed-end of Engineering Foundation.
a beam–column joint. Also, to consider axial force effect, the ba-
sic moment–curvature relation of the section which is uniquely
References
defined according to the dimensions of the concrete section and
the material properties of concrete and steel and axial force is [1] Kwak HG, Kim JK. Implementation of bond-slip effect in analyses of RC frames
modified. under cyclic loads using layered section method. Eng Struct 2006;28(12):
1715–27.
Correlation studies between analytical results and experimen- [2] Kwak HG, Kim JK. p − ∆ effect of slender RC columns under seismic load. Eng
tal values for the representative RC beams and columns have Struct 2007;29(11):3121–33.
[3] Park R, Kent DC, Sampson RA. Reinforced concrete members with cyclic
yielded the following conclusions: (1) to accurately predict the loading. J Struct Div, ASCE 1972;98(ST7):1341–60.
structural behavior of an RC beam-to-column sub-assemblage [4] Son JK, Fam A. Finite element modeling of hollow and concrete-filled fiber
where the nonlinear response is concentrated, a modification of composite tubes in flexure: Model development, verification and investigation
of tube parameters. Eng Struct 2008;30(10):2656–66.
the moment–curvature relation is strongly required to consider [5] Song HW, You DW, Byun KJ, Maekawa K. Finite element failure analysis of
fixed-end rotation caused by anchorage slip; (2) the proposed reinforced concrete T-girder bridges. Eng Struct 2002;24(2):151–62.
 H.-G. Kwak, S.-P. Kim / Engineering Structures 32 (2010) 69–79 79

[6] Taucer T, Spacone E, Filippou FC. A fiber beam–column element for seismic [19] Monti G, Filippou FC, Spacone E. Analysis of hysteretic behavior of anchored
response analysis of reinforced concrete structures, earthquake engrg. reinforcing bars. ACI Struct J 1997;94(2):248–61.
Research Center Report no. EERC 91-17. Berkeley (CA): Univ. of California; [20] Monti G, Filippou FC, Spacone E. Finite element for anchored bars under cyclic
1991. load reversals. J Struct Eng 1997;123(5):614–23.
[7] Zhu Z, Ahmad I, Mirmiran A. Fiber element modeling for seismic performance [21] de Groot AK, Kusters GMA, Monnier T. Numerical modeling of bond-slip
of bridge columns made of concrete-filled FRP tubes. Eng Struct 2006;28(14): behavior, heron, concrete mechanics, 1981; 26(1B): 6–38.
2023–35. [22] Ngo D, Scordelis AC. Finite element analysis of reinforced concrete beams. ACI
[8] Clough RW, Johnston SB. Effect of stiffness degradation on earthquake ductility
1967;64(3):152–63.
requirements. In: Proceedings of japan earthquake engineering symposium;
[23] Saatcioglu M, Ozcebe G. Response of reinforced concrete columns to simulated
1966.
seismic loading. ACI Struct J 1989;86(1):3–12.
[9] Roufaiel MSL, Meyer C. Analytical modeling of hysteretic behavior of
reinforced concrete frame. J Struct Eng 1987;113(3):429–44. [24] Soroushian P, Obaseki K, Nagi M, Rojas M. Pullout behaviour of hooked bars in
[10] Takeda T, Sozen MA, Nielsen NN. Reinforced concrete response to simulated exterior beam–column connections. ACI Struct J 1988;85(3):269–76.
earthquake. J Struct Div, ASCE 1970;96(12):2257–573. [25] Kwak HG, Kim SP. Nonlinear analysis of RC beams based on moment–
[11] Scott BD, Park R, Priestley MJN. Stress-strain behavior of concrete confi- curvature relation. Comput & Structures 2002;80(6):615–28.
ned by overlapping hoops at low and high strain rates. ACI 1982;79(1): [26] Gergely P, Lutz LA. Maximum crack width in reinforced concrete flexural
13–27. members. ACI special publication SP-20, Detroit: ACI; 1973. p. 87–117.
[12] Kwak HG, Filippou FC. Finite element analysis of reinforced concrete structures [27] Sawyer HA. Design of concrete frames for two failure states. In: Proceedings of
under monotonic loads. Report no. UCB/SEMM-90/14. Structural engineering the international symposium on the flexural mechanics of reinforced concrete,
mechanics and materials. Berkeley: University of California; 1990. ASCE-ACI, Miami, November, 1964. p. 405–431.
[13] Kwak HG, Kim SP. Nonlinear analysis of RC beams subject to cyclic loadings. [28] Bayrak O, Sheikh SA. High-strength concrete columns under simulated
J Struct Eng, ASCE 2001;127(12):1436–44. earthquake loading. ACI Struct J 1997;94(6):708–22.
[14] Kwak HG, Kim SP. Cyclic moment–curvature relation of RC beam. Mag Concr [29] Owen DRJ, Hinton E. Finite elements in plasticity. Pineridge Press Limited;
Res 2002;54(6):435–47. 1980.
[15] Kwak HG, Kim SP. Monotonic moment–curvature relation of RC beam. Mag
[30] Chen WF. Plasticity in reinforced concrete. McGraw-Hill Book Company; 1982.
Concr Res 2002;54(6):423–34.
[31] Ma SM, Bertero VV, Popov EP. Experimental and analytical studies on the
[16] Ayoub A, Filippou FC. Mixed Formulation of bond-slip problems under cyclic
hysteretic behavior of reinforced concrete rectangular and t-beam, earthquake
loads. J Struct Eng 1999;125(6):661–71.
[17] Eligehausen R, Popov EP, Bertero VV. Local bond stress-slip relationships of engrg. Research Center Report No. EERC 76-2. Berkeley (CA): Univ. Of
deformed bars under generalized excitations. Report No. UCB/EERC 83-23. California; 1976.
Berkeley: Earthquake Engineering Research Center, University of California; [32] Wight JK, Sozen MA. Strength decay of RC columns under shear reversals.
1983. J Struct Div, ASCE 1975;101(ST-5):1053–65.
[18] Hayashi S, Kokusho S. Bond behavior in the neighborhood of the crack. In: [33] Low SS, Moehle JP. Experimental study of reinforced concrete columns
Proceedings of the US–Japan joint seminar on finite element analysis of subjected to multi-axial cyclic loading, earthquake engrg. Research Center
reinforced concrete. 1985. p. 364–373. Report No. EERC 87-14, Berkeley (CA): Univ. of California; 1987.